From a93a5b2a7a315ae8d9d06b6974a4658602ce2106 Mon Sep 17 00:00:00 2001 From: Gregor Michalicek Date: Wed, 28 Mar 2018 18:57:12 +0200 Subject: [PATCH] Some more simplifications in main/cdngen.F90 --- cdn_mt/magMoms.f90 | 21 ++++++++++-- global/qfix.f90 | 37 +++++++++++---------- main/cdngen.F90 | 79 ++++++++++++++------------------------------ main/fleur.F90 | 2 +- main/mix.F90 | 2 +- optional/pldngen.f90 | 2 +- optional/stden.f90 | 2 +- 7 files changed, 67 insertions(+), 78 deletions(-) diff --git a/cdn_mt/magMoms.f90 b/cdn_mt/magMoms.f90 index c1165843..da2e83ac 100644 --- a/cdn_mt/magMoms.f90 +++ b/cdn_mt/magMoms.f90 @@ -8,7 +8,7 @@ MODULE m_magMoms CONTAINS -SUBROUTINE magMoms(dimension,input,atoms,noco,vTot,chmom,qa21) +SUBROUTINE magMoms(dimension,input,atoms,noco,vTot,stdn,svdn,chmom,qa21) USE m_types USE m_xmlOutput @@ -24,11 +24,28 @@ SUBROUTINE magMoms(dimension,input,atoms,noco,vTot,chmom,qa21) REAL, INTENT(INOUT) :: chmom(atoms%ntype,dimension%jspd) COMPLEX, INTENT(IN) :: qa21(atoms%ntype) + REAL, INTENT(IN) :: stdn(atoms%ntype,dimension%jspd) + REAL, INTENT(IN) :: svdn(atoms%ntype,dimension%jspd) INTEGER :: iType, j, iRepAtom - REAL :: smom + REAL :: sval,stot,scor,smom CHARACTER(LEN=20) :: attributes(4) + WRITE (6,FMT=8000) + WRITE (16,FMT=8000) + DO iType = 1,atoms%ntype + sval = svdn(iType,1) - svdn(iType,input%jspins) + stot = stdn(iType,1) - stdn(iType,input%jspins) + scor = stot - sval + WRITE (6,FMT=8010) iType,stot,sval,scor,svdn(iType,1),stdn(iType,1) + WRITE (16,FMT=8010) iType,stot,sval,scor,svdn(iType,1),stdn(iType,1) + END DO + + 8000 FORMAT (/,/,10x,'spin density at the nucleus:',/,10x,'type',t25,& + 'input%total',t42,'valence',t65,'core',t90,& + 'majority valence and input%total density',/) + 8010 FORMAT (i13,2x,3e20.8,5x,2e20.8) + WRITE (6,FMT=8020) WRITE (16,FMT=8020) diff --git a/global/qfix.f90 b/global/qfix.f90 index d02b6d09..5674c070 100644 --- a/global/qfix.f90 +++ b/global/qfix.f90 @@ -11,8 +11,8 @@ MODULE m_qfix ! qfix file no longer supported! CONTAINS - SUBROUTINE qfix( stars,atoms,sym,vacuum,& - sphhar,input,cell,oneD,den,l_printData,force_fix,fix) + SUBROUTINE qfix(stars,atoms,sym,vacuum,sphhar,input,cell,oneD,& + den,l_noco,l_printData,force_fix,fix) USE m_types USE m_cdntot @@ -29,7 +29,7 @@ CONTAINS TYPE(t_oneD),INTENT(IN) :: oneD TYPE(t_cell),INTENT(IN) :: cell TYPE(t_potden),INTENT(INOUT) :: den - LOGICAL,INTENT(IN) :: l_printData,force_fix + LOGICAL,INTENT(IN) :: l_noco,l_printData,force_fix REAL, INTENT (OUT) :: fix ! .. Local Scalars .. @@ -51,14 +51,13 @@ 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%pw,den%mt,den%vacz,.TRUE., qtot,qis) !The total nucleii charge zc=SUM(atoms%neq(:)*atoms%zatom(:)) zc = zc + 2*input%efield%sigma - IF (fixtotal) THEN !-roa fix = zc/qtot @@ -76,23 +75,27 @@ CONTAINS den%vacxy(:vacuum%nmzxy,:oneD%odi%nq2-1,:vacuum%nvac,:) END IF WRITE (6,FMT=8000) zc,fix - 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 closeXMLElement('fixedCharges') - !+roa ELSE fix = (zc - qtot) / qis + 1. den%pw(:stars%ng3,:) = fix*den%pw(:stars%ng3,:) WRITE (6,FMT=8001) zc,fix - 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 closeXMLElement('fixedCharges') ENDIF + IF (l_noco) THEN + !fix also the off-diagonal part of the density matrix + den%pw(:stars%ng3,3) = fix*den%pw(:stars%ng3,3) + IF (input%film.AND.fixtotal) THEN + den%vacz(:,:,3:4) = fix*den%vacz(:,:,3:4) + den%vacxy(:,:,:,3) = fix*den%vacxy(:,:,:,3) + END IF + 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 closeXMLElement('fixedCharges') + IF (fix>1.1) CALL juDFT_WARN("You lost too much charge") IF (fix<.9) CALL juDFT_WARN("You gained too much charge") diff --git a/main/cdngen.F90 b/main/cdngen.F90 index 7bd94b34..c7ea772c 100644 --- a/main/cdngen.F90 +++ b/main/cdngen.F90 @@ -70,7 +70,7 @@ SUBROUTINE cdngen(eig_id,mpi,input,banddos,sliceplot,vacuum,& TYPE(t_noco) :: noco_new !Local Scalars - REAL fix,qtot,scor,stot,sval,dummy + REAL fix,qtot,dummy INTEGER ivac,j,jspin,jspmax,k,iType LOGICAL l_enpara @@ -82,7 +82,7 @@ SUBROUTINE cdngen(eig_id,mpi,input,banddos,sliceplot,vacuum,& !pk non-collinear (start) INTEGER igq2_fft(0:stars%kq1_fft*stars%kq2_fft-1) - COMPLEX,ALLOCATABLE :: qa21(:), cdomvz(:,:) + COMPLEX,ALLOCATABLE :: qa21(:) !pk non-collinear (end) IF (mpi%irank.EQ.0) THEN @@ -133,55 +133,6 @@ SUBROUTINE cdngen(eig_id,mpi,input,banddos,sliceplot,vacuum,& stars,cell,sphhar,atoms,vTot,outDen,stdn,svdn) IF (mpi%irank.EQ.0) THEN - !block 2 unnecessary for slicing: begin - IF (.NOT.sliceplot%slice) THEN - CALL openXMLElementNoAttributes('allElectronCharges') - CALL qfix(stars,atoms,sym,vacuum, sphhar,input,cell,oneD,outDen,.TRUE.,.true.,fix) - CALL closeXMLElement('allElectronCharges') - !pk non-collinear (start) - IF (noco%l_noco) THEN - !fix also the off-diagonal part of the density matrix - outDen%pw(:stars%ng3,3) = fix*outDen%pw(:stars%ng3,3) - IF (input%film) THEN - outDen%vacz(:,:,3:4) = fix*outDen%vacz(:,:,3:4) - outDen%vacxy(:,:,:,3) = fix*outDen%vacxy(:,:,:,3) - END IF - END IF - !pk non-collinear (end) - - !spin densities at the nucleus - !and magnetic moment in the spheres - IF (input%jspins.EQ.2) THEN - WRITE (6,FMT=8000) - WRITE (16,FMT=8000) - DO iType = 1,atoms%ntype - sval = svdn(iType,1) - svdn(iType,input%jspins) - stot = stdn(iType,1) - stdn(iType,input%jspins) - scor = stot - sval - WRITE (6,FMT=8010) iType,stot,sval,scor,svdn(iType,1),stdn(iType,1) - WRITE (16,FMT=8010) iType,stot,sval,scor,svdn(iType,1),stdn(iType,1) - END DO - - noco_new = noco - - CALL magMoms(dimension,input,atoms,noco_new,vTot,chmom,qa21) - - !Generate and save the new nocoinp file if the directions of the local - !moments are relaxed or a constraint B-field is calculated. - IF (ANY(noco%l_relax(:atoms%ntype)).OR.noco%l_constr) THEN - CALL genNewNocoInp(input,atoms,jij,noco,noco_new) - END IF - - IF (noco%l_soc) CALL orbMagMoms(dimension,atoms,noco,clmom) - END IF - !block 2 unnecessary for slicing: end - END IF ! .NOT.sliceplot%slice - - 8000 FORMAT (/,/,10x,'spin density at the nucleus:',/,10x,'type',t25,& - 'input%total',t42,'valence',t65,'core',t90,& - 'majority valence and input%total density',/) - 8010 FORMAT (i13,2x,3e20.8,5x,2e20.8) - IF(vacuum%nvac.EQ.1) THEN outDen%vacz(:,2,:) = outDen%vacz(:,1,:) IF (sym%invs) THEN @@ -195,11 +146,31 @@ SUBROUTINE cdngen(eig_id,mpi,input,banddos,sliceplot,vacuum,& 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 - - ENDIF ! mpi%irank.EQ.0 + 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 closeXMLElement('allElectronCharges') + + IF (input%jspins.EQ.2) THEN + noco_new = noco + + !Calculate and write out spin densities at the nucleus and magnetic moments in the spheres + CALL magMoms(dimension,input,atoms,noco_new,vTot,stdn,svdn,chmom,qa21) + + !Generate and save the new nocoinp file if the directions of the local + !moments are relaxed or a constraint B-field is calculated. + IF (ANY(noco%l_relax(:atoms%ntype)).OR.noco%l_constr) THEN + CALL genNewNocoInp(input,atoms,jij,noco,noco_new) + END IF + + IF (noco%l_soc) CALL orbMagMoms(dimension,atoms,noco,clmom) + END IF + END IF ! mpi%irank.EQ.0 + #ifdef CPP_MPI CALL mpi_bc_potden(mpi,stars,sphhar,atoms,input,vacuum,oneD,noco,outDen) #endif @@ -207,8 +178,6 @@ SUBROUTINE cdngen(eig_id,mpi,input,banddos,sliceplot,vacuum,& DEALLOCATE (qvac,qvlay,qa21) DEALLOCATE (igq_fft) - RETURN - END SUBROUTINE cdngen END MODULE m_cdngen diff --git a/main/fleur.F90 b/main/fleur.F90 index d51fe780..54ad8961 100644 --- a/main/fleur.F90 +++ b/main/fleur.F90 @@ -184,7 +184,7 @@ CONTAINS CALL readDensity(stars,vacuum,atoms,cell,sphhar,input,sym,oneD,archiveType,CDN_INPUT_DEN_const,& 0,fermiEnergyTemp,l_qfix,inDen) CALL timestart("Qfix") - CALL qfix(stars,atoms,sym,vacuum, sphhar,input,cell,oneD,inDen,.FALSE.,.false.,fix) + CALL qfix(stars,atoms,sym,vacuum, sphhar,input,cell,oneD,inDen,noco%l_noco,.FALSE.,.false.,fix) CALL timestop("Qfix") CALL writeDensity(stars,vacuum,atoms,cell,sphhar,input,sym,oneD,archiveType,CDN_INPUT_DEN_const,& 0,-1.0,0.0,.FALSE.,inDen) diff --git a/main/mix.F90 b/main/mix.F90 index ae67ba93..96ac20d5 100644 --- a/main/mix.F90 +++ b/main/mix.F90 @@ -243,7 +243,7 @@ SUBROUTINE mix(stars,atoms,sphhar,vacuum,input,sym,cell,noco,oneD,& CALL closeXMLElement('densityConvergence') !fix charge of the new density - CALL qfix(stars,atoms,sym,vacuum, sphhar,input,cell,oneD,inDen,.FALSE.,.false., fix) + CALL qfix(stars,atoms,sym,vacuum, sphhar,input,cell,oneD,inDen,noco%l_noco,.FALSE.,.false., fix) IF(atoms%n_u.NE.n_u_keep) THEN inDen%mmpMat = n_mmpTemp diff --git a/optional/pldngen.f90 b/optional/pldngen.f90 index 003a6460..e23c6159 100644 --- a/optional/pldngen.f90 +++ b/optional/pldngen.f90 @@ -153,7 +153,7 @@ SUBROUTINE pldngen(sym,stars,atoms,sphhar,vacuum,& den%vacz(:,:,4) = AIMAG(cdomvz(:,:)) den%vacxy(:,:,:,3) = cdomvxy END IF - CALL qfix(stars,atoms,sym,vacuum,sphhar,input,cell,oneD,den,.FALSE.,.true.,fix) + CALL qfix(stars,atoms,sym,vacuum,sphhar,input,cell,oneD,den,noco%l_noco,.FALSE.,.true.,fix) rho(:,0:,1:,:input%jspins) = den%mt(:,0:,1:,:input%jspins) qpw(1:,:input%jspins) = den%pw(1:,:input%jspins) rht(1:,1:,:input%jspins) = den%vacz(1:,1:,:input%jspins) diff --git a/optional/stden.f90 b/optional/stden.f90 index d8b90eb4..b7106f91 100644 --- a/optional/stden.f90 +++ b/optional/stden.f90 @@ -204,7 +204,7 @@ SUBROUTINE stden(mpi,sphhar,stars,atoms,sym,DIMENSION,vacuum,& IF (mpi%irank == 0) THEN ! Check the normalization of total density - CALL qfix(stars,atoms,sym,vacuum,sphhar,input,cell,oneD,den,.FALSE.,.true.,fix) + CALL qfix(stars,atoms,sym,vacuum,sphhar,input,cell,oneD,den,.FALSE.,.FALSE.,.true.,fix) z=SUM(atoms%neq(:)*atoms%zatom(:)) IF (ABS(fix*z-z)>0.5) THEN CALL judft_warn("Starting density not charge neutral",hint= & -- GitLab