From fb8ed3241ce7740a22059f24070803eebed7328a Mon Sep 17 00:00:00 2001 From: giannozz Date: Tue, 13 Dec 2005 16:38:50 +0000 Subject: [PATCH] - it is possible to specify j=l-1/2 and j=l+1/2 states for relativistic calculations (rel=2) using "config" to pass the electronic configuratio - more checks on what is read in input - some cleanup git-svn-id: http://qeforge.qe-forge.org/svn/q-e/trunk/espresso@2584 c92efa57-630b-4861-b058-cf58834340f0 --- atomic/el_config.f90 | 2 ++ atomic/ld1_readin.f90 | 61 +++++++++++++++++++++++++++++-------------- atomic/ld1_setup.f90 | 21 +-------------- atomic/ld1inc.f90 | 1 - atomic/run_test.f90 | 6 +++++ 5 files changed, 50 insertions(+), 41 deletions(-) diff --git a/atomic/el_config.f90 b/atomic/el_config.f90 index 3883530fb..0f8e0fb26 100644 --- a/atomic/el_config.f90 +++ b/atomic/el_config.f90 @@ -146,8 +146,10 @@ do n=nwfc+1,nwf if (oc(n) > 2*(2*ll(n)+1)) & & call errore('el_config','wrong occupancy:'//occup,n) enddo +! ! pseudopotentials: N quantum number .ne. wavefunction label ! we find the lowest N for each L and reassign N as it +! if (.not.all_elec) then do l=0,3 n0=1000 diff --git a/atomic/ld1_readin.f90 b/atomic/ld1_readin.f90 index 49f027d32..687a8606e 100644 --- a/atomic/ld1_readin.f90 +++ b/atomic/ld1_readin.f90 @@ -169,11 +169,11 @@ subroutine ld1_readin nspin = 1 else if(lsd == 1) then nspin = 2 + if (rel == 2) call errore('ld1_readin', & + & 'local spin density and spin-orbit not allowed',1) else call errore('ld1_readin','lsd not correct',1) endif - if (rel == 2 .and. lsd == 1) call errore('ld1_readin', & - & 'local spin density and spin-orbit not allowed',1) if (config == ' ') then call read_config (rel, lsd, nwf, el, nn, ll, oc, isw, jj) @@ -182,12 +182,11 @@ subroutine ld1_readin ! ! check same labels corresponding to different spin or j value ! - do n=1,nwf + jj(1:nwf)=0.d0 + do n=1,nwf do i=n+1,nwf if (el(i) == el(n)) then - if ( lsd==0 ) then - call errore('ld1_readin',el(i)//' appears twice',i) - else if (rel == 2) then + if (rel == 2) then if (ll(n) > 0) then jj(n) = ll(n) + (isw(n)-1.5) jj(i) = ll(i) + (isw(i)-1.5) @@ -198,12 +197,23 @@ subroutine ld1_readin else call errore('ld1_readin',el(i)//' appears twice',i) end if + else if ( lsd==0 ) then + call errore('ld1_readin',el(i)//' appears twice',i) end if end if end do end do + if (rel == 2) isw(1:nwf)=1 end if ! + ! In the spin polarized or relativistic case adjust the occupations + ! + if (lsd == 1) then + call occ_spin(nwf,nwfx,el,nn,ll,oc,isw) + else if (rel == 2) then + call occ_spinorb(nwf,nwfx,el,nn,ll,jj,oc,isw) + endif + ! ! generate the radial grid - note that if iswitch = 2 the radial grid ! is not generated but read from the pseudopotential file ! @@ -212,19 +222,12 @@ subroutine ld1_readin rhoc=0.0_dp endif ! - ! In the spin polarized case adjust the occupations - ! - if (lsd == 1) call occ_spin(nwf,nwfx,el,nn,ll,oc,isw) - ! - ! If vdW calculation, read the input parameters - ! if (iswitch == 1) then ! ! no more data needed for AE calculations ! return - ! - + ! else if (iswitch == 3) then ! ! reading input for PP generation @@ -243,13 +246,16 @@ subroutine ld1_readin 500 call errore('ld1_readin','reading inputp',abs(ios)) if (lloc < 0 .and. rcloc <=0.0_dp) & - call errore('ld1_readin','rcloc is negative',1) + call errore('ld1_readin','rcloc must be positive',1) if (pseudotype < 1.or.pseudotype > 3) & call errore('ld1_readin','specify correct pseudotype',1) ! call read_psconfig (rel, lsd, nwfs, els, nns, lls, ocs, & isws, jjs, enls, rcut, rcutus ) ! + if (rel==2) call occ_spinorbps & + (nwfs,nwfsx,els,nns,lls,jjs,ocs,rcut,rcutus,enls,isws) + ! lmax = maxval(lls(1:nwfs)) ! zdum = zed @@ -303,10 +309,6 @@ subroutine ld1_readin jjts=0.0_dp jjtsc=0.0_dp - do n=1,nwf - oc_old(n)=oc(n) - enddo - nconf=1 configts=' ' @@ -369,6 +371,11 @@ subroutine ld1_readin if (lsd == 1) then call occ_spin(nwftsc(nc),nwfsx,eltsc(1,nc),nntsc(1,nc),lltsc(1,nc),& octsc(1,nc), iswtsc(1,nc)) + else if (rel == 2) then + call occ_spinorb(nwftsc(nc),nwfsx,eltsc(1,nc), & + & nntsc(1,nc),lltsc(1,nc),jjtsc(1,nc),octsc(1,nc),iswtsc(1,nc)) + else + jjtsc=0.0_dp endif end do ! @@ -542,6 +549,7 @@ subroutine read_config(rel, lsd, nwf, el, nn, ll, oc, isw, jj) ! do n=1,nwf if (rel < 2) then + jj(n) = 0.0_dp if (lsd == 0) then read(5,*,err=20,end=20,iostat=ios) & el(n), nn(n), ll(n), oc(n) @@ -559,12 +567,25 @@ subroutine read_config(rel, lsd, nwf, el, nn, ll, oc, isw, jj) el(n), nn(n), ll(n), oc(n), jj(n) isw(n)=1 if ((abs(ll(n)+0.5_dp-jj(n)) > 1.e-3_dp) .and. & - (abs(ll(n)-0.5_dp-jj(n)) > 1.e-3_dp) .and. abs(jj(n)) > 1.e-3_dp) & + (abs(ll(n)-0.5_dp-jj(n)) > 1.e-3_dp) .and. abs(jj(n)) > 1.e-3_dp) & call errore('read_config','jj wrong',n) if (oc(n) > (2.0_dp*jj(n)+1.0_dp) .and. abs(jj(n)) > 1e-3_dp) & call errore('read_config','occupations wrong',n) 22 call errore('read_config','reading orbital (rel)',abs(ios)) endif + ! + ! Check: no two same wavefunctions + ! + do ncheck=1,n-1 + if ( el(ncheck) == el(n) .and. isw(ncheck) == isw(n) .and. & + jj(ncheck) == jj(n) ) then + call errore('read_config', & + 'same wavefunction '//el(n)//' appears twice',n) + endif + enddo + ! + ! More sanity checks + ! write(label,'(a2)') el(n) read (label,'(i1)') ncheck if (ncheck /= nn(n) .or. & diff --git a/atomic/ld1_setup.f90 b/atomic/ld1_setup.f90 index 17786731c..19a180bae 100644 --- a/atomic/ld1_setup.f90 +++ b/atomic/ld1_setup.f90 @@ -18,7 +18,7 @@ subroutine ld1_setup use funct, only : get_iexch, dft_is_meta !, set_dft_from_name implicit none - integer n, n1, nc, nwf0, m, nwftot, ios + integer n, n1, nc, m, nwftot, ios logical ok, hf, oep ! ! transform dft in a series of codes for the exchange and @@ -37,25 +37,10 @@ subroutine ld1_setup ! call set_sl3(sl3,lmx) ! - ! split the spin-orbit states - ! - if (rel == 2) then - call occ_spinorb(nwf,nwfx,el,nn,ll,jj,oc,isw) - do n=1,nwf - oc_old(n)=oc(n) - enddo - else - jj=0.d0 - endif - ! ! make the correspondence all-electron pseudopotential - ! and split the spin-orbit states ! if (iswitch >= 2) then do nc=1, nconf - if (rel.eq.2) & - call occ_spinorb(nwftsc(nc),nwfsx,eltsc(1,nc), & - & nntsc(1,nc),lltsc(1,nc),jjtsc(1,nc),octsc(1,nc),iswtsc(1,nc)) do n=1,nwftsc(nc) nstoaec(n,nc)=0 do n1=1,nwf @@ -86,13 +71,9 @@ subroutine ld1_setup ! if (iswitch == 3) then isws=1 - nwf0=nwfs - if (rel==2) & - call occ_spinorbps(nwfs,nwfsx,els,nns,lls,jjs,ocs,rcut,rcutus,enls,isws) do n=1,nwf core_state(n)=.true. enddo - do n=1,nwfs nstoae(n)=0 do n1=1,nwf diff --git a/atomic/ld1inc.f90 b/atomic/ld1inc.f90 index 93efc4a3d..fe0bcafe7 100644 --- a/atomic/ld1inc.f90 +++ b/atomic/ld1inc.f90 @@ -26,7 +26,6 @@ module ld1inc real(DP) :: & jj(nwfx), & ! the total angular momentum oc(nwfx), & ! the occupations of the all-electron atom - oc_old(nwfx), & ! saves the occupations of the all-electron zed, & ! the ionic charge enne, & ! the number of electrons sl3(0:lmx2,0:lmx2,0:lmx2) diff --git a/atomic/run_test.f90 b/atomic/run_test.f90 index a4201ad67..2df6e0614 100644 --- a/atomic/run_test.f90 +++ b/atomic/run_test.f90 @@ -16,6 +16,8 @@ subroutine run_test use ld1inc implicit none + real(DP) :: oc_old(nwfx) + integer & n, & ! counter on wavefunctions n1,& ! counter on mesh points @@ -26,6 +28,10 @@ subroutine run_test character(len=1) :: nch real(DP) :: dum + do n=1,nwf + oc_old(n)=oc(n) + enddo + file_tests = trim(prefix)//'.test' open(unit=13, file=file_tests, iostat=ios, err=1111, status='unknown') 1111 call errore('ld1_setup','opening file_tests',abs(ios))