Rearrange alphdiff calculation in main/cdngen.F90 and cdn/m_perp.f90

cdn/m_perp.f90

@@ -6,7 +6,7 @@
 
 MODULE m_m_perp
 CONTAINS 
-  SUBROUTINE m_perp(atoms,itype,noco,vr0, chmom,qa21,alphdiff)
+  SUBROUTINE m_perp(atoms,itype,iRepAtom,noco,vr0, chmom,qa21)
     !***********************************************************************
     ! calculates the perpendicular part of the local moment.
     ! if l_relax is true the angle of the output local moment is calculated
@@ -16,8 +16,8 @@ CONTAINS
     ! Philipp Kurz 2000-02-09
     !***********************************************************************
 
+    USE m_constants
     USE m_intgr, ONLY : intgr3
-    USE m_constants, ONLY : fpi_const
     USE m_polangle
     USE m_rotdenmat
     USE m_types
@@ -26,11 +26,10 @@ CONTAINS
     TYPE(t_atoms),INTENT(IN)     :: atoms
 
     !     .. Scalar Arguments ..
-    INTEGER, INTENT (IN) :: itype
+    INTEGER, INTENT (IN) :: itype, iRepAtom
     !     ..
     !     .. Array Arguments ..
     REAL, INTENT    (IN) :: chmom(:,:)!(atoms%ntype,dimension%jspd)
-    REAL, INTENT    (IN) :: alphdiff(atoms%ntype) 
     REAL, INTENT    (IN) :: vr0(:,:,:)!(atoms%jmtd,atoms%ntype,jspd)
     COMPLEX, INTENT (IN) :: qa21(atoms%ntype)
     !     ..
@@ -38,12 +37,20 @@ CONTAINS
     INTEGER iri
     REAL b_xavh,scale,b_con_outx,b_con_outy,mx,my,mz,&
          &     alphh,betah,mz_tmp,mx_mix,my_mix,mz_mix
-    REAL    rho11,rho22
+    REAL    rho11,rho22, alphdiff
     COMPLEX rho21
     !     ..
     !     .. Local Arrays ..
     REAL b_xc_h(atoms%jmtd),b_xav(atoms%ntype)
 
+    ! angles in nocoinp file are (alph-alphdiff)
+    IF (noco%l_ss) THEN
+       alphdiff = 2.0*pi_const*(noco%qss(1)*atoms%taual(1,iRepAtom) + &
+                                noco%qss(2)*atoms%taual(2,iRepAtom) + &
+                                noco%qss(3)*atoms%taual(3,iRepAtom) )
+    ELSE
+       alphdiff = 0.0
+    END IF
 
     !---> calculated the comp. of the local moment vector
     mx = 2*REAL(qa21(itype))
@@ -71,10 +78,10 @@ CONTAINS
        my = 2*AIMAG(rho21)
        mz = rho11 - rho22
        CALL pol_angle(mx,my,mz,betah,alphh)
-       WRITE  (6,8027) noco%beta(itype),noco%alph(itype)-alphdiff(itype)
-       WRITE (16,8027) noco%beta(itype),noco%alph(itype)-alphdiff(itype)
-       WRITE  (6,8028) betah,alphh-alphdiff(itype)
-       WRITE (16,8028) betah,alphh-alphdiff(itype)
+       WRITE  (6,8027) noco%beta(itype),noco%alph(itype)-alphdiff
+       WRITE (16,8027) noco%beta(itype),noco%alph(itype)-alphdiff
+       WRITE  (6,8028) betah,alphh-alphdiff
+       WRITE (16,8028) betah,alphh-alphdiff
 8027   FORMAT(2x,'-->',10x,' input noco%beta=',f9.5, '  input noco%alpha=',f9.5)
 8028   FORMAT(2x,'-->',10x,'output noco%beta=',f9.5, ' output noco%alpha=',f9.5)
 
@@ -101,8 +108,8 @@ CONTAINS
        ENDIF
        ! calculate angles alpha and beta in global frame
        CALL pol_angle(mx_mix,my_mix,mz_mix,betah,alphh)
-       WRITE  (6,8029) betah,alphh-alphdiff(itype)
-       WRITE (16,8029) betah,alphh-alphdiff(itype)
+       WRITE  (6,8029) betah,alphh-alphdiff
+       WRITE (16,8029) betah,alphh-alphdiff
 8029   FORMAT(2x,'-->',10x,' new noco%beta  =',f9.5, '  new noco%alpha  =',f9.5)
        noco%alph(itype) = alphh
        noco%beta(itype) = betah

main/cdngen.F90

@@ -91,7 +91,7 @@ SUBROUTINE cdngen(eig_id,mpi,input,banddos,sliceplot,vacuum,&
    REAL   ,ALLOCATABLE :: qvac(:,:,:,:),qvlay(:,:,:,:,:)
 
    !pk non-collinear (start)
-   REAL    rhoint,momint,alphdiff(atoms%ntype)
+   REAL    rhoint,momint,alphdiff
    INTEGER igq2_fft(0:stars%kq1_fft*stars%kq2_fft-1)
    COMPLEX,ALLOCATABLE :: qa21(:), cdomvz(:,:)
    !pk non-collinear (end)
@@ -326,24 +326,17 @@ SUBROUTINE cdngen(eig_id,mpi,input,banddos,sliceplot,vacuum,&
                WRITE (16,FMT=8010) n,stot,sval,scor,svdn(n,1),stdn(n,1)
             END DO
 
+            noco_new = noco
+
             CALL magMoms(dimension,input,atoms,chmom)
 
-            noco_new = noco
             IF (noco%l_mperp) THEN
-               ! angles in nocoinp file are (alph-alphdiff)
                iatom = 1
                DO n = 1,atoms%ntype
-                  IF (noco%l_ss) THEN
-                     alphdiff(n)= 2.*pi_const*(noco%qss(1)*atoms%taual(1,iatom) + &
-                                               noco%qss(2)*atoms%taual(2,iatom) + &
-                                               noco%qss(3)*atoms%taual(3,iatom) )
-                  ELSE
-                     alphdiff(n)= 0.
-                  END IF
                   !calculate the perpendicular part of the local moment
                   !and relax the angle of the local moment or calculate
                   !the constraint B-field.
-                  CALL m_perp(atoms,n,noco_new,vTot%mt(:,0,:,:),chmom,qa21,alphdiff)
+                  CALL m_perp(atoms,n,iatom,noco_new,vTot%mt(:,0,:,:),chmom,qa21)
                   iatom= iatom + atoms%neq(n)
                END DO
             END IF
@@ -351,7 +344,6 @@ SUBROUTINE cdngen(eig_id,mpi,input,banddos,sliceplot,vacuum,&
             !save the new nocoinp file if the dierctions of the local
             !moments are relaxed or a constraint B-field is calculated.
             l_relax_any = .false.
-            iatom = 1
             DO itype = 1,atoms%ntype
                l_relax_any = l_relax_any.OR.noco%l_relax(itype)
             END DO
@@ -359,9 +351,13 @@ SUBROUTINE cdngen(eig_id,mpi,input,banddos,sliceplot,vacuum,&
                IF (.not. noco%l_mperp) THEN
                   CALL juDFT_error ("(l_relax_any.OR.noco).AND.(.NOT. )" ,calledby ="cdngen")
                END IF
+               iatom = 1
                DO itype = 1, atoms%ntype
                   IF (noco%l_ss) THEN
-                     noco_new%alph(itype) = noco%alph(itype) - alphdiff(itype)
+                     alphdiff = 2.0*pi_const*(noco%qss(1)*atoms%taual(1,iatom) + &
+                                              noco%qss(2)*atoms%taual(2,iatom) + &
+                                              noco%qss(3)*atoms%taual(3,iatom) )
+                     noco_new%alph(itype) = noco%alph(itype) - alphdiff
                      DO WHILE (noco_new%alph(n) > +pi_const)
                         noco_new%alph(n)= noco_new%alph(n) - 2.*pi_const
                      END DO
@@ -371,6 +367,7 @@ SUBROUTINE cdngen(eig_id,mpi,input,banddos,sliceplot,vacuum,&
                   ELSE
                      noco_new%alph(itype) = noco%alph(itype)
                   END IF
+                  iatom= iatom + atoms%neq(n)
                END DO
 
                OPEN (24,file='nocoinp',form='formatted', status='old')