Even more simplifications in main/cdngen.F90

cdn/cdntot.f90

@@ -4,11 +4,8 @@
 !     vacuum, and mt regions      c.l.fu
 !     ********************************************************
       CONTAINS
-      SUBROUTINE cdntot(&
-     &                  stars,atoms,sym,&
-     &                  vacuum,input,cell,oneD,&
-     &                  qpw,rho,rht,l_printData,&
-     &                  qtot,qistot)
+      SUBROUTINE cdntot(stars,atoms,sym,vacuum,input,cell,oneD,&
+                        den,l_printData,qtot,qistot)
 
       USE m_intgr, ONLY : intgr3
       USE m_constants
@@ -19,27 +16,20 @@
       USE m_convol
       USE m_xmlOutput
       IMPLICIT NONE
-!     ..
+
 !     .. Scalar Arguments ..
-      TYPE(t_stars),INTENT(IN) :: stars
-      TYPE(t_atoms),INTENT(IN) :: atoms
-      TYPE(t_sym),INTENT(IN)   :: sym
-      TYPE(t_vacuum),INTENT(IN):: vacuum
-      TYPE(t_input),INTENT(IN) :: input
-      TYPE(t_oneD),INTENT(IN)  :: oneD
-      TYPE(t_cell),INTENT(IN)  :: cell
-      LOGICAL,INTENT(IN)       :: l_printData
-      REAL,    INTENT (OUT):: qtot,qistot
-!     ..
-!     .. Array Arguments ..
-      COMPLEX, INTENT (IN) :: qpw(stars%ng3,input%jspins)
-      REAL,    INTENT (IN) :: rho(:,0:,:,:) !(atoms%jmtd,0:sphhar%nlhd,atoms%ntype,input%jspins)
-      REAL,    INTENT (IN) :: rht(vacuum%nmzd,2,input%jspins)
-!-odim
-!+odim
-!     ..
+      TYPE(t_stars),INTENT(IN)  :: stars
+      TYPE(t_atoms),INTENT(IN)  :: atoms
+      TYPE(t_sym),INTENT(IN)    :: sym
+      TYPE(t_vacuum),INTENT(IN) :: vacuum
+      TYPE(t_input),INTENT(IN)  :: input
+      TYPE(t_oneD),INTENT(IN)   :: oneD
+      TYPE(t_cell),INTENT(IN)   :: cell
+      TYPE(t_potden),INTENT(IN) :: den
+      LOGICAL,INTENT(IN)        :: l_printData
+      REAL,INTENT(OUT)          :: qtot,qistot
+
 !     .. Local Scalars ..
-    ! COMPLEX x
       COMPLEX x(stars%ng3)
       REAL q,qis,w,mtCharge
       INTEGER i,ivac,j,jspin,n,nz
@@ -66,7 +56,7 @@
 !     -----mt charge
          CALL timestart("MT")
          DO 10 n = 1,atoms%ntype
-            CALL intgr3(rho(:,0,n,jspin),atoms%rmsh(:,n),atoms%dx(n),atoms%jri(n),w)
+            CALL intgr3(den%mt(:,0,n,jspin),atoms%rmsh(:,n),atoms%dx(n),atoms%jri(n),w)
             qmt(n) = w*sfp_const
             q = q + atoms%neq(n)*qmt(n)
    10    CONTINUE
@@ -77,9 +67,9 @@
                DO nz = 1,vacuum%nmz
                   IF (oneD%odi%d1) THEN
                      rht1(nz,ivac,jspin) = (cell%z1+(nz-1)*vacuum%delz)*&
-     &                    rht(nz,ivac,jspin)
+     &                    den%vacz(nz,ivac,jspin)
                   ELSE
-                     rht1(nz,ivac,jspin) =  rht(nz,ivac,jspin)
+                     rht1(nz,ivac,jspin) =  den%vacz(nz,ivac,jspin)
                   END IF
                END DO
                CALL qsf(vacuum%delz,rht1(1,ivac,jspin),q2,vacuum%nmz,0)
@@ -93,7 +83,7 @@
          END IF
 !     -----is region
          IF (.not.judft_was_Argument("-oldfix")) THEN
-            CALL convol(stars,x,qpw(:,jspin),stars%ufft)
+            CALL convol(stars,x,den%pw(:,jspin),stars%ufft)
             qis = x(1)*cell%omtil
          ELSE
           qis = 0.
@@ -104,14 +94,14 @@
 !     >                 symor,tau,mrot,rmt,sk3,bmat,ig2,ig,
 !     >                 kv3(1,j),
 !     <                 x)
-!            qis = qis + qpw(j,jspin)*x*nstr(j)
+!            qis = qis + den%pw(j,jspin)*x*nstr(j)
 !   30    CONTINUE
          CALL pwint_all(&
      &                 stars,atoms,sym,oneD,&
      &                 cell,&
      &                 x)
          DO j = 1,stars%ng3
-             qis = qis + qpw(j,jspin)*x(j)*stars%nstr(j)
+             qis = qis + den%pw(j,jspin)*x(j)*stars%nstr(j)
          ENDDO
          endif
          qistot = qistot + qis

cdn/cdnval.F90

@@ -980,6 +980,16 @@ CONTAINS
           !-for
        END DO ! end of loop ispin = jsp_start,jsp_end
        CALL closeXMLElement('mtCharges')
+
+       IF(vacuum%nvac.EQ.1) THEN
+          den%vacz(:,2,:) = den%vacz(:,1,:)
+          IF (sym%invs) THEN
+             den%vacxy(:,:,2,:) = CONJG(den%vacxy(:,:,1,:))
+          ELSE
+             den%vacxy(:,:,2,:) = den%vacxy(:,:,1,:)
+          END IF
+       END IF
+
     END IF ! end of (mpi%irank==0)
     !+t3e
     !Note: no deallocation anymore, we rely on Fortran08 :-)
@@ -988,5 +998,6 @@ CONTAINS
        IF ((banddos%dos.OR.banddos%vacdos).AND.(banddos%ndir/=-2))  CALL juDFT_end("DOS OK",mpi%irank)
        IF (vacuum%nstm.EQ.3)  CALL juDFT_end("VACWAVE OK",mpi%irank)
     END IF
+
   END SUBROUTINE cdnval
 END MODULE m_cdnval

global/qfix.f90

@@ -51,8 +51,7 @@ CONTAINS
     ! qfix==0 means no qfix was given in inp.xml. 
     ! In this case do nothing except when forced to fix!
     
-    CALL cdntot(stars,atoms,sym, vacuum,input,cell,oneD,&
-                den%pw,den%mt,den%vacz,.TRUE., qtot,qis)
+    CALL cdntot(stars,atoms,sym,vacuum,input,cell,oneD,den,.TRUE.,qtot,qis)
 
     !The total nucleii charge
     zc=SUM(atoms%neq(:)*atoms%zatom(:))
@@ -91,9 +90,9 @@ CONTAINS
     END IF
 
     IF (ABS(fix-1.0)<1.E-6) RETURN !no second calculation of cdntot as nothing was fixed
+
     CALL openXMLElementNoAttributes('fixedCharges')
-    CALL cdntot(stars,atoms,sym, vacuum,input,cell,oneD,&
-                den%pw,den%mt,den%vacz,l_printData, qtot,qis)
+    CALL cdntot(stars,atoms,sym,vacuum,input,cell,oneD,den,l_printData,qtot,qis)
     CALL closeXMLElement('fixedCharges')
 
     IF (fix>1.1) CALL juDFT_WARN("You lost too much charge")

main/cdngen.F90

@@ -71,11 +71,11 @@ SUBROUTINE cdngen(eig_id,mpi,input,banddos,sliceplot,vacuum,&
 
    !Local Scalars
    REAL fix,qtot,dummy
-   INTEGER ivac,j,jspin,jspmax,k,iType
+   INTEGER jspin,jspmax
    LOGICAL l_enpara
 
    !Local Arrays
-   REAL stdn(atoms%ntype,dimension%jspd),svdn(atoms%ntype,dimension%jspd),alpha_l(atoms%ntype)
+   REAL stdn(atoms%ntype,dimension%jspd),svdn(atoms%ntype,dimension%jspd)
    REAL chmom(atoms%ntype,dimension%jspd),clmom(3,atoms%ntype,dimension%jspd)
    INTEGER,ALLOCATABLE :: igq_fft(:)
    REAL   ,ALLOCATABLE :: qvac(:,:,:,:),qvlay(:,:,:,:,:)
@@ -117,42 +117,32 @@ SUBROUTINE cdngen(eig_id,mpi,input,banddos,sliceplot,vacuum,&
    DO jspin = 1,jspmax
       CALL timestart("cdngen: cdnval")
       CALL cdnval(eig_id,&
-                  mpi,kpts,jspin,sliceplot,noco, input,banddos,cell,atoms,enpara,stars, vacuum,dimension,&
+                  mpi,kpts,jspin,sliceplot,noco,input,banddos,cell,atoms,enpara,stars,vacuum,dimension,&
                   sphhar,sym,obsolete,igq_fft,vTot,oneD,coreSpecInput,&
-                  outDen,results,qvac,qvlay,qa21, chmom,clmom)
+                  outDen,results,qvac,qvlay,qa21,chmom,clmom)
       CALL timestop("cdngen: cdnval")
    END DO
 
    IF (mpi%irank.EQ.0) THEN
       IF (l_enpara) CLOSE (40)
-      CALL cdntot(stars,atoms,sym, vacuum,input,cell,oneD, outDen%pw,outDen%mt,outDen%vacz,.TRUE., qtot,dummy)
+      CALL cdntot(stars,atoms,sym, vacuum,input,cell,oneD,outDen,.TRUE.,qtot,dummy)
       CALL closeXMLElement('valenceDensity')
    END IF ! mpi%irank = 0
 
    CALL cdncore(results,mpi,dimension,oneD,sliceplot,input,vacuum,noco,sym,&
                 stars,cell,sphhar,atoms,vTot,outDen,stdn,svdn)
 
-   IF (mpi%irank.EQ.0) THEN
-      IF(vacuum%nvac.EQ.1) THEN
-         outDen%vacz(:,2,:) = outDen%vacz(:,1,:)
-         IF (sym%invs) THEN
-            outDen%vacxy(:,:,2,:) = CONJG(outDen%vacxy(:,:,1,:))
-         ELSE
-            outDen%vacxy(:,:,2,:) = outDen%vacxy(:,:,1,:)
-         END IF
-      END IF
-
-      IF (sliceplot%slice) THEN
+   IF (sliceplot%slice) THEN
+      IF (mpi%irank.EQ.0) THEN
          CALL writeDensity(stars,vacuum,atoms,cell,sphhar,input,sym,oneD,archiveType,CDN_INPUT_DEN_const,&
                            1,-1.0,0.0,.FALSE.,outDen,'cdn_slice')
       END IF
+      CALL juDFT_end("slice OK",mpi%irank)
    END IF
 
-   IF (sliceplot%slice) CALL juDFT_end("slice OK",mpi%irank)
-
    IF (mpi%irank.EQ.0) THEN
       CALL openXMLElementNoAttributes('allElectronCharges')
-      CALL qfix(stars,atoms,sym,vacuum, sphhar,input,cell,oneD,outDen,noco%l_noco,.TRUE.,.true.,fix)
+      CALL qfix(stars,atoms,sym,vacuum,sphhar,input,cell,oneD,outDen,noco%l_noco,.TRUE.,.true.,fix)
       CALL closeXMLElement('allElectronCharges')
 
       IF (input%jspins.EQ.2) THEN