Line data Source code
1 : !--------------------------------------------------------------------------------
2 : ! Copyright (c) 2016 Peter Grünberg Institut, Forschungszentrum Jülich, Germany
3 : ! This file is part of FLEUR and available as free software under the conditions
4 : ! of the MIT license as expressed in the LICENSE file in more detail.
5 : !--------------------------------------------------------------------------------
6 :
7 : MODULE m_get_mt_perturbation
8 : USE m_juDFT
9 : USE m_polangle
10 : USE m_types
11 : USE m_constants
12 : USE m_mt_tofrom_grid
13 :
14 : IMPLICIT NONE
15 :
16 : CONTAINS
17 :
18 0 : SUBROUTINE get_mt_local_perturbation(atoms,sphhar,sym,noco,den,den1,den1im)
19 : TYPE(t_atoms), INTENT(IN) :: atoms
20 : TYPE(t_sphhar), INTENT(IN) :: sphhar
21 : TYPE(t_sym), INTENT(IN) :: sym
22 : TYPE(t_noco), INTENT(IN) :: noco
23 : TYPE(t_potden), INTENT(INOUT) :: den, den1, den1im
24 :
25 0 : TYPE(t_gradients) :: grad
26 :
27 : INTEGER :: n, nsp, imesh, i
28 : REAL :: rho_11, rho_22, m
29 : REAL :: rhotot, rho_up, rho_down, theta, phi
30 : COMPLEX :: n1, mx1, my1, mz1, m1, t1, p1, rho1_up, rho1_down
31 : REAL :: eps=1E-10
32 : REAL, ALLOCATABLE :: ch(:,:), chre(:,:), chim(:,:)
33 :
34 0 : nsp=atoms%nsp()
35 0 : ALLOCATE(ch(nsp*atoms%jmtd,4), den1%theta_mt(nsp*atoms%jmtd,atoms%ntype), &
36 0 : den1%phi_mt(nsp*atoms%jmtd,atoms%ntype))
37 0 : ALLOCATE(chre(nsp*atoms%jmtd,4), chim(nsp*atoms%jmtd,4))
38 0 : ALLOCATE(den1im%theta_mt(nsp*atoms%jmtd,atoms%ntype), den1im%phi_mt(nsp*atoms%jmtd,atoms%ntype))
39 :
40 0 : CALL init_mt_grid(4,atoms,sphhar,.FALSE.,sym)
41 :
42 0 : DO n=1,atoms%ntype
43 0 : CALL mt_to_grid(.FALSE., 2, atoms, sym, sphhar, .FALSE., den%mt(:,0:,n,:), n, noco, grad, ch)
44 0 : CALL mt_to_grid(.FALSE., 4, atoms, sym, sphhar, .FALSE., den1%mt(:,0:,n,:), n, noco, grad, chre)
45 0 : CALL mt_to_grid(.FALSE., 4, atoms, sym, sphhar, .FALSE., den1im%mt(:,0:,n,:), n, noco, grad, chim)
46 :
47 0 : DO imesh = 1, nsp*atoms%jri(n)
48 0 : rho_11 = ch(imesh,1)
49 0 : rho_22 = ch(imesh,2)
50 :
51 0 : m = rho_11 - rho_22
52 :
53 : ! Calculate perturbed total and magnetization density
54 0 : n1 = chre(imesh,1) + ImagUnit * chim(imesh,1)
55 0 : mx1 = chre(imesh,2) + ImagUnit * chim(imesh,2)
56 0 : my1 = chre(imesh,3) + ImagUnit * chim(imesh,3)
57 0 : mz1 = chre(imesh,4) + ImagUnit * chim(imesh,4)
58 :
59 0 : theta = den%theta_mt(imesh,n)
60 0 : phi = den%phi_mt(imesh,n)
61 :
62 : ! Calculate the perturbed absolute value of the magnetization density
63 0 : m1 = cmplx(0.0, 0.0)
64 0 : m1 = m1 + mx1 * SIN(theta) * COS(phi)
65 0 : m1 = m1 + my1 * SIN(theta) * SIN(phi)
66 0 : m1 = m1 + mz1 * COS(theta)
67 :
68 : ! Calculate the perturbed angles
69 0 : t1 = cmplx(0.0, 0.0)
70 0 : t1 = t1 + mx1 * COS(theta) * COS(phi) / m
71 0 : t1 = t1 + my1 * COS(theta) * SIN(phi) / m
72 0 : t1 = t1 - mz1 * SIN(theta) / m
73 :
74 0 : p1 = cmplx(0.0, 0.0)
75 0 : p1 = p1 - mx1 * SIN(theta) * SIN(phi) / m
76 0 : p1 = p1 + my1 * SIN(theta) * COS(phi) / m
77 :
78 0 : rho1_up = (n1 + m1)/2
79 0 : rho1_down = (n1 - m1)/2
80 :
81 0 : chre(imesh,1) = REAL(rho1_up )
82 0 : chre(imesh,2) = REAL(rho1_down)
83 0 : chim(imesh,1) = AIMAG(rho1_up )
84 0 : chim(imesh,2) = AIMAG(rho1_down)
85 0 : den1%theta_mt(imesh,n) = REAL(t1)
86 0 : den1%phi_mt(imesh,n) = REAL(p1)
87 0 : den1im%theta_mt(imesh,n) = AIMAG(t1)
88 0 : den1im%phi_mt(imesh,n) = AIMAG(p1)
89 : END DO
90 0 : den1%mt(:,0:,n,:)=0.0
91 0 : den1im%mt(:,0:,n,:)=0.0
92 0 : CALL mt_from_grid(atoms,sym,sphhar,n,2,chre,den1%mt(:,0:,n,:))
93 0 : CALL mt_from_grid(atoms,sym,sphhar,n,2,chim,den1im%mt(:,0:,n,:))
94 0 : DO i=1,atoms%jri(n)
95 0 : den1%mt(i,:,n,:)=den1%mt(i,:,n,:)*atoms%rmsh(i,n)**2
96 0 : den1im%mt(i,:,n,:)=den1im%mt(i,:,n,:)*atoms%rmsh(i,n)**2
97 : END DO
98 : END DO
99 :
100 0 : CALL finish_mt_grid()
101 :
102 0 : END SUBROUTINE get_mt_local_perturbation
103 :
104 0 : SUBROUTINE get_mt_global_perturbation(atoms,sphhar,sym,den,den1,den1im,noco,vtot,vtot1,vtot1im)
105 : TYPE(t_atoms), INTENT(IN) :: atoms
106 : TYPE(t_sphhar), INTENT(IN) :: sphhar
107 : TYPE(t_sym), INTENT(IN) :: sym
108 : TYPE(t_potden), INTENT(IN) :: den, den1, den1im
109 : TYPE(t_noco), INTENT(IN) :: noco
110 : TYPE(t_potden), INTENT(IN) :: vtot
111 : TYPE(t_potden), INTENT(INOUT) :: vtot1,vtot1im
112 :
113 0 : TYPE(t_gradients) :: grad
114 0 : TYPE(t_potden) :: vtot_loc
115 :
116 : INTEGER :: n, nsp, imesh, i
117 : REAL :: theta, phi, v11, v22
118 : COMPLEX :: v1up, v1down, v1eff, b1eff
119 0 : REAL, ALLOCATABLE :: ch(:,:), chtmp(:,:), chretmp(:,:), chimtmp(:,:)
120 :
121 : COMPLEX :: a11, a22, a21, a12, av11, av22, av21, av12
122 : COMPLEX :: t1, p1, v21, v12, v1mat11, v1mat22, v1mat21, v1mat12
123 :
124 0 : vtot_loc = vtot
125 0 : nsp=atoms%nsp()
126 0 : ALLOCATE(chtmp(nsp*atoms%jmtd,4))
127 0 : ALLOCATE(chretmp(nsp*atoms%jmtd,4))
128 0 : ALLOCATE(chimtmp(nsp*atoms%jmtd,4))
129 : !TODO: Make sure the indices for rho1 are 1,2,3,4 == n1,mx1,my1,mz1
130 0 : CALL init_mt_grid(4,atoms,sphhar,.FALSE.,sym)
131 0 : DO n=1,atoms%ntype
132 :
133 0 : DO i=1,atoms%jri(n)
134 0 : vtot_loc%mt(i,:,n,:)=vtot%mt(i,:,n,:)*atoms%rmsh(i,n)**2
135 0 : vtot1%mt(i,:,n,:)=vtot1%mt(i,:,n,:)*atoms%rmsh(i,n)**2
136 0 : vtot1im%mt(i,:,n,:)=vtot1im%mt(i,:,n,:)*atoms%rmsh(i,n)**2
137 : END DO
138 :
139 0 : CALL mt_to_grid(.FALSE.,4,atoms,sym,sphhar,.FALSE.,vtot_loc%mt(:,0:,n,:),n,noco,grad,chtmp(:,1:2))
140 0 : CALL mt_to_grid(.FALSE.,2,atoms,sym,sphhar,.FALSE.,vtot1%mt(:,0:,n,:),n,noco,grad,chretmp(:,1:2))
141 0 : CALL mt_to_grid(.FALSE.,2,atoms,sym,sphhar,.FALSE.,vtot1im%mt(:,0:,n,:),n,noco,grad,chimtmp(:,1:2))
142 :
143 0 : DO imesh = 1, nsp*atoms%jri(n)
144 0 : v11 = chtmp(imesh, 1)
145 0 : v22 = chtmp(imesh, 2)
146 0 : v21 = chtmp(imesh, 3) + ImagUnit * chtmp(imesh, 4)
147 0 : v12 = chtmp(imesh, 3) - ImagUnit * chtmp(imesh, 4)
148 :
149 0 : v1up = chretmp(imesh,1) + ImagUnit * chimtmp(imesh,1)
150 0 : v1down = chretmp(imesh,2) + ImagUnit * chimtmp(imesh,2)
151 :
152 0 : theta = den%theta_mt(imesh,n)
153 0 : phi = den%phi_mt(imesh,n)
154 :
155 0 : v1eff = (v1up + v1down)/2.0
156 0 : b1eff = (v1up - v1down)/2.0
157 :
158 0 : v1mat11 = v1eff + b1eff * COS(theta) !11
159 0 : v1mat22 = v1eff - b1eff * COS(theta) !22
160 0 : v1mat21 = b1eff * SIN(theta)*EXP( Imagunit*phi) !21
161 0 : v1mat12 = b1eff * SIN(theta)*EXP(-Imagunit*phi) !12
162 :
163 0 : t1 = den1%theta_mt(imesh,n) + ImagUnit * den1im%theta_mt(imesh,n)
164 0 : p1 = den1%phi_mt(imesh,n) + ImagUnit * den1im%phi_mt(imesh,n)
165 :
166 0 : a11 = -ImagUnit * p1 / 2.0
167 0 : a22 = ImagUnit * p1 / 2.0
168 0 : a21 = EXP( ImagUnit*phi) * t1 / 2.0
169 0 : a12 = -EXP(-ImagUnit*phi) * t1 / 2.0
170 :
171 0 : av11 = a11 * v11 + a12 * v21 !11
172 0 : av22 = a21 * v12 + a22 * v22 !22
173 0 : av21 = a21 * v11 + a22 * v21 !21
174 0 : av12 = a11 * v12 + a12 * v22 !12
175 :
176 0 : v1mat11 = v1mat11 + av11 + CONJG(av11)
177 0 : v1mat22 = v1mat22 + av22 + CONJG(av22)
178 0 : v1mat21 = v1mat21 + av21 + CONJG(av12)
179 0 : v1mat12 = v1mat12 + av12 + CONJG(av21)
180 :
181 0 : chretmp(imesh, 1) = REAL(v1mat11)
182 0 : chretmp(imesh, 2) = REAL(v1mat22)
183 0 : chretmp(imesh, 3) = REAL(v1mat21)
184 0 : chretmp(imesh, 4) = REAL(v1mat12)
185 :
186 0 : chimtmp(imesh, 1) = AIMAG(v1mat11)
187 0 : chimtmp(imesh, 2) = AIMAG(v1mat22)
188 0 : chimtmp(imesh, 3) = AIMAG(v1mat21)
189 0 : chimtmp(imesh, 4) = AIMAG(v1mat12)
190 :
191 : END DO
192 :
193 0 : vtot1%mt(:,0:,n,:)=0.0
194 0 : vtot1im%mt(:,0:,n,:)=0.0
195 :
196 0 : CALL mt_from_grid(atoms,sym,sphhar,n,4,chretmp,vtot1%mt(:,0:,n,:))
197 0 : CALL mt_from_grid(atoms,sym,sphhar,n,4,chimtmp,vtot1im%mt(:,0:,n,:))
198 :
199 : END DO
200 :
201 0 : CALL finish_mt_grid()
202 :
203 0 : END SUBROUTINE get_mt_global_perturbation
204 :
205 : END MODULE m_get_mt_perturbation
|
