Line data Source code
1 : !--------------------------------------------------------------------------------
2 : ! Copyright (c) 2018 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 : ! This routine allows to rotate the cdn in a way that the direction of magnetization aligns with the direction of the spin quantization axis.
7 : ! This routine also allows to reverse the rotation by using the angles stored in atoms (nococonv%alph,nococonv%beta) which are generated by the
8 : ! routine magnMomFromDen.
9 : !
10 : ! Robin Hilgers, Nov '19 Adaption to new nococonv type in Feb '20, Added RelaxMixing parameter + Allow Relaxation of alpha and Beta individually Jul '20'
11 : MODULE m_Relaxspinaxismagn
12 :
13 : USE m_magnMomFromDen
14 : USE m_types
15 : USE m_types_fleurinput
16 : USE m_flipcdn
17 : USE m_constants
18 : USE m_polangle
19 :
20 : IMPLICIT NONE
21 :
22 : CONTAINS
23 :
24 : !Rotates cdn to global frame at initialization before the scf loop.
25 0 : SUBROUTINE initRelax(noco,nococonv,atoms,input,vacuum,sphhar,stars,sym ,cell,den)
26 : TYPE(t_input), INTENT(IN) :: input
27 : TYPE(t_atoms), INTENT(IN) :: atoms
28 : TYPE(t_noco), INTENT(IN) :: noco
29 : TYPE(t_nococonv), INTENT(INOUT) :: nococonv
30 : TYPE(t_stars), INTENT(IN) :: stars
31 : TYPE(t_vacuum), INTENT(IN) :: vacuum
32 : TYPE(t_sphhar), INTENT(IN) :: sphhar
33 : TYPE(t_sym), INTENT(IN) :: sym
34 :
35 : TYPE(t_cell), INTENT(IN) :: cell
36 : TYPE(t_potden), INTENT(INOUT) :: den
37 :
38 0 : REAL ::zeros(atoms%ntype)
39 0 : zeros(:)=0.0
40 0 : ALLOCATE(nococonv%alphPrev(atoms%ntype),nococonv%betaPrev(atoms%ntype))
41 0 : nococonv%alphPrev=noco%alph_inp
42 0 : nococonv%betaPrev=noco%beta_inp
43 :
44 0 : CALL flipcdn(atoms,input,vacuum,sphhar,stars,sym,noco ,cell,zeros,noco%beta_inp,den)
45 0 : CALL flipcdn(atoms,input,vacuum,sphhar,stars,sym,noco ,cell,noco%alph_inp,zeros,den)
46 0 : nococonv%alph=zeros
47 0 : nococonv%beta=zeros
48 :
49 0 : END SUBROUTINE initRelax
50 :
51 664 : SUBROUTINE precond_noco(it,vacuum,sphhar,stars,sym ,cell,noco,nococonv,input,atoms,inden,outden,fsm)
52 : use m_types_mixvector
53 : INTEGER,INTENT(IN) :: it
54 : TYPE(t_input), INTENT(IN) :: input
55 : TYPE(t_atoms), INTENT(IN) :: atoms
56 : TYPE(t_noco), INTENT(IN) :: noco
57 : TYPE(t_nococonv), INTENT(IN) :: nococonv
58 : TYPE(t_stars), INTENT(IN) :: stars
59 : TYPE(t_vacuum), INTENT(IN) :: vacuum
60 : TYPE(t_sphhar), INTENT(IN) :: sphhar
61 : TYPE(t_sym), INTENT(IN) :: sym
62 :
63 : TYPE(t_cell), INTENT(IN) :: cell
64 : TYPE(t_potden), INTENT(IN) :: inden
65 : TYPE(t_potden), INTENT(INOUT) :: outden
66 : TYPE(t_mixvector), INTENT(INOUT) :: fsm
67 :
68 1814 : if (.not.(noco%l_noco.and.any(noco%l_alignMT))) return
69 :
70 6 : select case (noco%mag_mixing_scheme)
71 : case(1)
72 0 : if (it>1) return
73 0 : call precond_noco_anglerotate(vacuum,sphhar,stars,sym ,cell,noco,nococonv,input,atoms,inden,outden,fsm)
74 : case(2)
75 0 : call precond_noco_anglerotate(vacuum,sphhar,stars,sym ,cell,noco,nococonv,input,atoms,inden,outden,fsm)
76 : case(3)
77 6 : call precond_noco_densitymatrix(vacuum,sphhar,stars,sym ,cell,noco,nococonv,input,atoms,inden,outden,fsm)
78 : end select
79 : END subroutine precond_noco
80 :
81 : !Preconditioner to control relaxation of the direction of the magnetic moment
82 0 : SUBROUTINE precond_noco_anglerotate(vacuum,sphhar,stars,sym ,cell,noco,nococonv,input,atoms,inden,outden,fsm)
83 : use m_types_mixvector
84 : TYPE(t_input), INTENT(IN) :: input
85 : TYPE(t_atoms), INTENT(IN) :: atoms
86 : TYPE(t_noco), INTENT(IN) :: noco
87 : TYPE(t_nococonv), INTENT(IN) :: nococonv
88 : TYPE(t_stars), INTENT(IN) :: stars
89 : TYPE(t_vacuum), INTENT(IN) :: vacuum
90 : TYPE(t_sphhar), INTENT(IN) :: sphhar
91 : TYPE(t_sym), INTENT(IN) :: sym
92 :
93 : TYPE(t_cell), INTENT(IN) :: cell
94 : TYPE(t_potden), INTENT(IN) :: inden
95 : TYPE(t_potden), INTENT(INOUT) :: outden
96 : TYPE(t_mixvector), INTENT(OUT) :: fsm
97 :
98 0 : real,dimension(atoms%ntype) :: dphi,dtheta,zeros
99 0 : TYPE(t_potden) :: delta_den,outden_rot
100 : integer :: n
101 0 : zeros(:) = 0.0
102 : !Put outden in local frame of inden
103 0 : CALL flipcdn(atoms,input,vacuum,sphhar,stars,sym,noco ,cell,-nococonv%alphPrev,zeros,outden)
104 0 : CALL flipcdn(atoms,input,vacuum,sphhar,stars,sym,noco ,cell,zeros,-nococonv%betaPrev,outden)
105 : !rotation angle
106 0 : CALL gimmeAngles(input,atoms,noco,vacuum,sphhar,stars,outden,dPhi,dtheta)
107 : !dphi = dphi *(noco%mix_RelaxWeightOffD-1.0)
108 0 : dtheta = dtheta *(noco%mix_RelaxWeightOffD-1.0)
109 :
110 : !if (any(abs(dphi)>2.0).or.any(abs(dtheta)>2.0)) THEN
111 : ! print *,"No precond"
112 : ! dphi=0.0
113 : !dtheta=0.0
114 : !endif
115 :
116 : !Scale Off-diagonal parts
117 : !DO n=1,atoms%ntype
118 : ! outden%mt(:,0:,n,3)=outden%mt(:,0:,n,3)*noco%mix_RelaxWeightOffD(n)
119 : !outden%mt(:,0:,n,4)=outden%mt(:,0:,n,4)*noco%mix_RelaxWeightOffD(n)
120 : !ENDDO
121 :
122 : !CALL cureTooSmallAngles(atoms,dphi,dtheta)
123 0 : CALL flipcdn(atoms,input,vacuum,sphhar,stars,sym,noco ,cell,zeros,dtheta,outden)
124 : !CALL flipcdn(atoms,input,vacuum,sphhar,stars,sym,noco ,cell,dphi,zeros,outden)
125 : !Rotate back in global frame
126 0 : CALL flipcdn(atoms,input,vacuum,sphhar,stars,sym,noco ,cell,zeros,nococonv%betaPrev,outden)
127 0 : CALL flipcdn(atoms,input,vacuum,sphhar,stars,sym,noco ,cell,nococonv%alphPrev,zeros,outden)
128 :
129 0 : CALL gimmeAngles(input,atoms,noco,vacuum,sphhar,stars,outden,dPhi,dtheta)
130 :
131 0 : call delta_den%subPotDen(outden,inden)
132 : !CALL flipcdn(atoms,input,vacuum,sphhar,stars,sym,noco ,cell,-dphi,zeros,delta_den)
133 : !CALL flipcdn(atoms,input,vacuum,sphhar,stars,sym,noco ,cell,zeros,-dtheta,delta_den)
134 :
135 0 : call fsm%alloc()
136 0 : call fsm%from_density(delta_den,vacuum%nmzxyd)
137 :
138 0 : END SUBROUTINE precond_noco_anglerotate
139 :
140 : !Preconditioner to control relaxation of the direction of the magnetic moment
141 0 : SUBROUTINE precond_noco_densitymatrix(vacuum,sphhar,stars,sym ,cell,noco,nococonv,input,atoms,inden,outden,fsm)
142 : use m_types_mixvector
143 : TYPE(t_input), INTENT(IN) :: input
144 : TYPE(t_atoms), INTENT(IN) :: atoms
145 : TYPE(t_noco), INTENT(IN) :: noco
146 : TYPE(t_nococonv), INTENT(IN) :: nococonv
147 : TYPE(t_stars), INTENT(IN) :: stars
148 : TYPE(t_vacuum), INTENT(IN) :: vacuum
149 : TYPE(t_sphhar), INTENT(IN) :: sphhar
150 : TYPE(t_sym), INTENT(IN) :: sym
151 :
152 : TYPE(t_cell), INTENT(IN) :: cell
153 : TYPE(t_potden), INTENT(IN) :: inden,outden
154 : TYPE(t_mixvector), INTENT(OUT) :: fsm
155 :
156 0 : TYPE(t_potden) :: delta_den
157 : integer :: n
158 0 : real,dimension(atoms%ntype) :: zeros,theta,phi
159 0 : zeros(:) = 0.0
160 0 : call delta_den%subPotDen(outden,inden)
161 : !Put in local frame
162 0 : CALL gimmeAngles(input,atoms,noco,vacuum,sphhar,stars,inden,Phi,theta)
163 0 : zeros(:) = 0.0
164 0 : CALL flipcdn(atoms,input,vacuum,sphhar,stars,sym,noco ,cell,-phi,zeros,delta_den)
165 0 : CALL flipcdn(atoms,input,vacuum,sphhar,stars,sym,noco ,cell,zeros,-theta,delta_den)
166 : !Scale Off-diagonal parts
167 0 : DO n=1,atoms%ntype
168 0 : delta_den%mt(:,0:,n,3)=delta_den%mt(:,0:,n,3)*noco%mix_RelaxWeightOffD(n)
169 0 : delta_den%mt(:,0:,n,4)=delta_den%mt(:,0:,n,4)*noco%mix_RelaxWeightOffD(n)
170 : ENDDO
171 : !Put back in global frame
172 0 : zeros(:) = 0.0
173 0 : CALL flipcdn(atoms,input,vacuum,sphhar,stars,sym,noco ,cell,zeros,theta,delta_den)
174 0 : CALL flipcdn(atoms,input,vacuum,sphhar,stars,sym,noco ,cell,phi,zeros,delta_den)
175 :
176 0 : call fsm%alloc()
177 0 : call fsm%from_density(delta_den,vacuum%nmzxyd)
178 :
179 :
180 0 : END SUBROUTINE precond_noco_densitymatrix
181 :
182 : !Purges to small angles below 10^-4 rad to 0. => Stabilizes convergence.
183 0 : SUBROUTINE cureTooSmallAngles(atoms,angleA,angleB)
184 : TYPE(t_atoms),INTENT(IN) :: atoms
185 : REAL ,INTENT(INOUT) :: angleA(:)
186 : REAL ,INTENT(INOUT), OPTIONAL :: angleB(:)
187 :
188 : REAL :: eps
189 : INTEGER :: i
190 0 : eps=0.0001
191 0 : DO i=1, atoms%ntype
192 0 : IF (abs(angleA(i)).LE.eps) angleA(i)=0.0
193 0 : IF(PRESENT(angleB).AND.abs(angleB(i)).LE.eps) angleB(i)=0.0
194 : END DO
195 :
196 0 : END SUBROUTINE cureTooSmallAngles
197 :
198 : !Calculates angles from magnetization and assigns correct sign to be used in the rotation (flipcdn) routine properly.
199 4 : SUBROUTINE Gimmeangles(Input,Atoms,Noco,Vacuum,Sphhar,Stars,Den,Phitemp,Thetatemp)
200 : use m_types_nococonv
201 : TYPE(t_input) ,INTENT(IN) :: input
202 : TYPE(t_atoms) ,INTENT(IN) :: atoms
203 : TYPE(t_noco) ,INTENT(IN) :: noco
204 : TYPE(t_stars) ,INTENT(IN) :: stars
205 : TYPE(t_vacuum),INTENT(IN) :: vacuum
206 : TYPE(t_sphhar),INTENT(IN) :: sphhar
207 : TYPE(t_potden),INTENT(IN) :: den
208 : REAL ,INTENT(OUT) :: phiTemp(atoms%ntype),thetaTemp(atoms%ntype)
209 :
210 4 : REAL :: moments(3,atoms%ntype)
211 :
212 16 : type(t_nococonv):: nococonv
213 4 : call nococonv%avg_moments(den,atoms,moments,thetatemp,phitemp)
214 : !!CALL magnMomFromDen(input,atoms,noco,den,moments,thetaTemp,phiTemp)
215 12 : phiTemp(:)=(-1)*phiTemp(:)
216 12 : thetaTemp(:)=(-1)*thetaTemp(:)
217 :
218 4 : END SUBROUTINE gimmeAngles
219 :
220 : !Rotates from global frame into current local frame
221 666 : SUBROUTINE toLocalSpinFrame(fmpi,vacuum,sphhar,stars&
222 : ,sym ,cell,noco,nococonv,input,atoms,l_adjust,den,l_update_nococonv)
223 : use m_constants
224 : TYPE(t_mpi),INTENT(IN) :: fmpi
225 : TYPE(t_input), INTENT(IN) :: input
226 : TYPE(t_atoms), INTENT(IN) :: atoms
227 : TYPE(t_noco), INTENT(IN) :: noco
228 : TYPE(t_nococonv), INTENT(INOUT) :: nococonv
229 : TYPE(t_stars),INTENT(IN) :: stars
230 : TYPE(t_vacuum),INTENT(IN) :: vacuum
231 : TYPE(t_sphhar),INTENT(IN) :: sphhar
232 : TYPE(t_sym),INTENT(IN) :: sym
233 :
234 : TYPE(t_cell),INTENT(IN) :: cell
235 : LOGICAL,INTENT(IN) :: l_adjust
236 : TYPE(t_potden),INTENT(INOUT) :: den
237 : LOGICAL,OPTIONAL,INTENT(IN) :: l_update_nococonv
238 :
239 666 : REAL :: zeros(atoms%ntype),alph_old(atoms%ntype),dalph
240 : integer :: n
241 1816 : if (.not.any(noco%l_alignMT)) RETURN
242 8 : if (fmpi%irank==0) THEN
243 12 : zeros(:) = 0.0
244 :
245 12 : alph_old=nococonv%alphprev
246 4 : if (l_adjust) then
247 : !if (.not.allocated(nococonv%alphPrev)) allocate(nococonv%alphprev(atoms%ntype),nococonv%betaprev(atoms%ntype))
248 4 : call Gimmeangles(input,atoms,noco,vacuum,sphhar,stars,den,nococonv%alphPrev,nococonv%betaPrev)
249 4 : write(oUnit,*) "Current magnetization angles"
250 12 : DO n=1,atoms%ntype
251 12 : write(oUnit,*) n,nococonv%alphPrev(n),nococonv%betaprev(n)
252 : ENDDO
253 :
254 : endif
255 : !Now try to minimize difference to previous angles
256 : !DO n=1,atoms%ntype
257 : ! dalph=abs(alph_old(n)-nococonv%alphPrev(n))
258 : ! if (abs(nococonv%alph(n)-nococonv%alphprev(n)-Pi_const)<dalph) THEN
259 : ! nococonv%alphprev(n)=nococonv%alphprev(n)+pi_const
260 : ! nococonv%betaprev(n)=-1*nococonv%betaprev(n)
261 : ! elseif (abs(nococonv%alph(n)-nococonv%alphprev(n)+Pi_const)<dalph) THEN
262 : ! nococonv%alphprev(n)=nococonv%alphprev(n)-pi_const
263 : ! nococonv%betaprev(n)=-1*nococonv%betaprev(n)
264 : ! endif
265 : !enddo
266 :
267 :
268 24 : CAlL flipcdn(atoms,input,vacuum,sphhar,stars,sym,noco ,cell,merge(nococonv%alphprev,zeros,noco%l_alignMT),merge(nococonv%betaprev,zeros,noco%l_alignMT),Den,toGlobal=.false.)
269 :
270 : !CALL flipcdn(atoms,input,vacuum,sphhar,stars,sym,noco ,cell,merge(-nococonv%alphPrev,zeros,noco%l_alignMT),zeros,den)
271 : !CALL flipcdn(atoms,input,vacuum,sphhar,stars,sym,noco ,cell,zeros,merge(-nococonv%betaPreV,zeros,noco%l_alignMT),den)
272 4 : if (present(l_update_nococonv)) then
273 4 : if (l_update_nococonv) THEN
274 16 : nococonv%alph=merge(nococonv%alphPrev,nococonv%alph,noco%l_alignMT)
275 16 : nococonv%beta=merge(nococonv%betaPrev,nococonv%beta,noco%l_alignMT)
276 12 : nococonv%alphPrev=0.0
277 12 : nococonv%betaPrev=0.0
278 : ENDIF
279 : ENDIF
280 : endif
281 8 : call den%distribute(fmpi%mpi_comm)
282 8 : call nococonv%mpi_bc(fmpi%mpi_comm)
283 :
284 :
285 : END SUBROUTINE
286 :
287 : !Rotates into the global frame so mixing can be performed without any restrictions. => Compatible with anderson mixing scheme.
288 1338 : SUBROUTINE toGlobalSpinFrame(noco,nococonv,vacuum,sphhar,stars&
289 : ,sym ,cell,input,atoms, den,fmpi,l_update_nococonv)
290 :
291 : TYPE(t_mpi),INTENT(IN),OPTIONAL :: fmpi
292 : TYPE(t_input), INTENT(IN) :: input
293 : TYPE(t_atoms), INTENT(IN) :: atoms
294 : TYPE(t_noco), INTENT(IN) :: noco
295 : TYPE(t_nococonv), INTENT(INOUT) :: nococonv
296 : TYPE(t_stars),INTENT(IN) :: stars
297 : TYPE(t_vacuum),INTENT(IN) :: vacuum
298 : TYPE(t_sphhar),INTENT(IN) :: sphhar
299 : TYPE(t_sym),INTENT(IN) :: sym
300 :
301 : TYPE(t_cell),INTENT(IN) :: cell
302 : TYPE(t_potden), INTENT(INOUT) :: Den
303 : LOGICAL,OPTIONAL,INTENT(IN) :: l_update_nococonv
304 :
305 1338 : REAL :: zeros(atoms%ntype)
306 :
307 : LOGICAL l_irank0
308 3654 : if (.not.any(noco%l_alignMT)) RETURN
309 :
310 14 : l_irank0=.true.
311 14 : if (present(fmpi)) l_irank0=fmpi%irank==0
312 :
313 14 : if (l_irank0) then
314 24 : zeros(:)=0.0
315 48 : CAlL flipcdn(atoms,input,vacuum,sphhar,stars,sym,noco ,cell,merge(nococonv%alph,zeros,noco%l_alignMT),merge(nococonv%beta,zeros,noco%l_alignMT),Den,toGlobal=.true.)
316 : !CAlL flipcdn(atoms,input,vacuum,sphhar,stars,sym,noco ,cell,zeros,merge(nococonv%beta,zeros,noco%l_alignMT),Den)
317 : !CALL flipcdn(atoms,input,vacuum,sphhar,stars,sym,noco ,cell,merge(nococonv%alph,zeros,noco%l_alignMT),zeros,Den)
318 : ! Nococonv is zero now since rotation has been reverted.
319 8 : if (present(l_update_nococonv)) THEN
320 3 : if (l_update_nococonv) THEN
321 12 : nococonv%alphPrev=merge(nococonv%alph,nococonv%alphPrev,noco%l_alignMT)
322 12 : nococonv%betaPrev=merge(nococonv%beta,nococonv%betaPrev,noco%l_alignMT)
323 12 : nococonv%alph=merge(zeros,nococonv%alph,noco%l_alignMT)
324 12 : nococonv%beta=merge(zeros,nococonv%beta,noco%l_alignMT)
325 : ENDIF
326 : ENDIF
327 : ENDIF
328 14 : if (present(fmpi)) then
329 12 : call den%distribute(fmpi%mpi_comm)
330 12 : call nococonv%mpi_bc(fmpi%mpi_comm)
331 : endif
332 :
333 : END SUBROUTINE
334 :
335 : END MODULE m_RelaxSpinAxisMagn
336 :
337 : #ifdef CPP_NEVER
338 : !Relaxation routine for only relaxing alpha. How it works: See two routines below.
339 : SUBROUTINE alphaRelax(vacuum,sphhar,stars&
340 : ,sym ,cell,noco,nococonv,input,atoms,den)
341 : TYPE(t_input), INTENT(IN) :: input
342 : TYPE(t_atoms), INTENT(IN) :: atoms
343 : TYPE(t_noco), INTENT(IN) :: noco
344 : TYPE(t_nococonv),INTENT(INOUT) :: nococonv
345 : TYPE(t_stars), INTENT(IN) :: stars
346 : TYPE(t_vacuum), INTENT(IN) :: vacuum
347 : TYPE(t_sphhar), INTENT(IN) :: sphhar
348 : TYPE(t_sym), INTENT(IN) :: sym
349 :
350 : TYPE(t_cell), INTENT(IN) :: cell
351 : TYPE(t_potden), INTENT(INOUT) :: den
352 :
353 : REAL :: moments(3,atoms%ntype)
354 : REAL :: diffT(atoms%ntype),diffP(atoms%ntype), zeros(atoms%ntype)
355 : zeros(:)=0.0
356 :
357 : CALL gimmeAngles(input,atoms,noco,vacuum,sphhar,stars,den,diffP,diffT)
358 : CALL cureTooSmallAngles(atoms,diffT,diffP)
359 : diffP=diffP-nococonv%alphPrev
360 : diffP=diffP*noco%mix_RelaxWeightOffD
361 : CALL flipcdn(atoms,input,vacuum,sphhar,stars,sym,noco ,cell,-diffP-nococonv%alphPrev,zeros,den)
362 : CALL flipcdn(atoms,input,vacuum,sphhar,stars,sym,noco ,cell,zeros,-nococonv%betaRlx,den)
363 : nococonv%beta=nococonv%betaRlx
364 : nococonv%betaPrev=nococonv%betaRlx
365 : nococonv%alph=nococonv%alphPrev+diffP
366 : nococonv%alphPrev=nococonv%alph
367 : CALL cureTooSmallAngles(atoms,nococonv%alph)
368 :
369 : END SUBROUTINE alphaRelax
370 :
371 : !Relaxation routine for only relaxing beta. How it works: See one routine below.
372 : SUBROUTINE betaRelax(vacuum,sphhar,stars&
373 : ,sym ,cell,noco,nococonv,input,atoms,den)
374 : TYPE(t_input), INTENT(IN) :: input
375 : TYPE(t_atoms), INTENT(IN) :: atoms
376 : TYPE(t_noco), INTENT(IN) :: noco
377 : TYPE(t_nococonv),INTENT(INOUT) :: nococonv
378 : TYPE(t_stars), INTENT(IN) :: stars
379 : TYPE(t_vacuum), INTENT(IN) :: vacuum
380 : TYPE(t_sphhar), INTENT(IN) :: sphhar
381 : TYPE(t_sym), INTENT(IN) :: sym
382 :
383 : TYPE(t_cell), INTENT(IN) :: cell
384 : TYPE(t_potden), INTENT(INOUT) :: den
385 :
386 : LOGICAL :: nonZeroAngles
387 : REAL :: moments(3,atoms%ntype)
388 : REAL :: diffT(atoms%ntype),diffP(atoms%ntype), zeros(atoms%ntype)
389 : zeros(:)=0.0
390 :
391 : CALL gimmeAngles(input,atoms,noco,vacuum,sphhar,stars,den,diffP,diffT)
392 : CALL cureTooSmallAngles(atoms,diffT,diffP)
393 : diffT=diffT-nococonv%betaPrev
394 : diffT=diffT*noco%mix_RelaxWeightOffD
395 : CALL flipcdn(atoms,input,vacuum,sphhar,stars,sym,noco ,cell,-nococonv%alphRlx,zeros,den)
396 : CALL flipcdn(atoms,input,vacuum,sphhar,stars,sym,noco ,cell,zeros,-diffT-nococonv%betaPrev,den)
397 : nococonv%beta=nococonv%betaPrev+diffT
398 : nococonv%betaPrev=nococonv%beta
399 : nococonv%alph=nococonv%alphRlx
400 : CALL cureTooSmallAngles(atoms,nococonv%beta)
401 :
402 : END SUBROUTINE betaRelax
403 :
404 : !Relaxation routine for both angles at the same time. Calculates the angles by which the magnetization direction changed
405 : !based on the angle nococonv%*Prev which store the angles from the previous iteration. The resulting angular difference is then weighted (mix_RelaxWeightOffD)
406 : !and the cdn will be manipulated so that the magnetization direction is rotated towards the new magnetization direction by AngularDifference*Weight.
407 : !If both angles are relaxed: weight=1 direction of magnetization is || to SQA after relaxation.
408 : SUBROUTINE bothRelax(vacuum,sphhar,stars&
409 : ,sym ,cell,noco,nococonv,input,atoms,den)
410 : TYPE(t_input), INTENT(IN) :: input
411 : TYPE(t_atoms), INTENT(IN) :: atoms
412 : TYPE(t_noco), INTENT(IN) :: noco
413 : TYPE(t_nococonv),INTENT(INOUT) :: nococonv
414 : TYPE(t_stars), INTENT(IN) :: stars
415 : TYPE(t_vacuum), INTENT(IN) :: vacuum
416 : TYPE(t_sphhar), INTENT(IN) :: sphhar
417 : TYPE(t_sym), INTENT(IN) :: sym
418 :
419 : TYPE(t_cell), INTENT(IN) :: cell
420 : TYPE(t_potden), INTENT(INOUT) :: den
421 :
422 : LOGICAL :: nonZeroAngles
423 : REAL :: moments(3,atoms%ntype)
424 : REAL :: diffT(atoms%ntype),diffP(atoms%ntype), zeros(atoms%ntype)
425 : zeros(:)=0.0
426 :
427 : CALL gimmeAngles(input,atoms,noco,vacuum,sphhar,stars,den,diffP,diffT)
428 : CALL cureTooSmallAngles(atoms,diffT,diffP)
429 : diffT=diffT-nococonv%betaPrev
430 : diffT=diffT*noco%mix_RelaxWeightOffD
431 : diffP=diffP-nococonv%alphPrev
432 : diffP=diffP*noco%mix_RelaxWeightOffD
433 : CALL flipcdn(atoms,input,vacuum,sphhar,stars,sym,noco ,cell,-diffP-nococonv%alphPrev,zeros,den)
434 : CALL flipcdn(atoms,input,vacuum,sphhar,stars,sym,noco ,cell,zeros,-diffT-nococonv%betaPrev,den)
435 : nococonv%beta=nococonv%betaPrev+diffT
436 : nococonv%betaPrev=nococonv%beta
437 : nococonv%alph=nococonv%alphPrev+diffP
438 : nococonv%alphPrev=nococonv%alph
439 : CALL cureTooSmallAngles(atoms,nococonv%beta,nococonv%alph)
440 :
441 : END SUBROUTINE bothRelax
442 : #endif
|