[Wien] setrmt_lapw #4
John Rundgren
jru at kth.se
Tue Oct 2 16:09:31 CEST 2012
Dear Professor Blaha,
Attached are vcoul curves in Ryd, LM=0,0, with rmt=1.92,1.74 from
setrmt_lapw and rmt=2.0,1.6 from the TiO2 example in UG.
vcoul(r<RMT) is practically invariant with respect to rmt, whereas
vcoul(r>rmt) is sensitive to rmt. The latter behavior is remarkable.
Is something overlooked among the modifications? For instance, must clm
be extended first and vcoul afterwards?
Best regards,
John
On Tue, 2012-10-02 at 11:54 +0200, John Rundgren wrote:
> Dear Professor Blaha,
>
> Many thanks for the final remarks on the radial extension of vcoul. A
> similar extension of clmsum is also needed (no spin).
>
> In lapw0.F a statement concerning clm is found, line 563,
> READ(8,2021) (CLM(J,LM1,JATOM), J=1,JRI(JATOM))
> leading down to lines 1109:1141 where summations occur just as in the
> Fourier extension of vcoul, namely,
> DO J=2,NKK ...
> DO LM1=1,LLMM ...
> CQ(lm1)=CQ(lm1)+RHOKT*YKA(LM1)*BES(LL+1)
>
> I should be most obliged for a piece of code in which CQ is transformed
> into
> new_clm(j,lm1,jatom), j=1,jrx(jatom), jrx(jatom)>jri(jatom)
>
> Best regards,
> John Rundgren
>
>
> On Sun, 2012-09-30 at 20:35 +0200, Peter Blaha wrote:
> > I can see that now the potential is at least continuous, but the slopes
> > still do not match.
> >
> > I don't have the equations with me at home, but as far as I can see:
> >
> > I hope, you have redimensioned the arrays r, jrx properly.
> >
> > > elseif(r(irad)>rmt(jatom))then
> > > SURFIN=CVOUT(irad,LM1)*(-4.D0*PI)*(RMT(JATOM)/R(IRAD))**L* &
> > > IMAG**L*VOUTF
> > > endif
> >
> > I don't think this is ok. The modification of cvout by r(i)/r(RMT) is valid only
> > inside the spheres, where the coulomb potential depends on on both, the
> > charge inside the sphere and the (scaled) potential at RMT.
> >
> > Outside the sphere, we only need the bessel expansion of the plane waves, i.e.
> >
> > SURFIN=CVOUT(irad,LM1)*(-4.D0*PI)*IMAG**L*VOUTF
> > ----------------------------------------------------------------------------------
> > PS: We want the coefficients V_LM(r) of V(r) = SUM_LM (V_LM(r) Y_LM(r-hat)) (1)
> >
> > In the interstitial the potential is given as:
> > V(r) = 1/vol sum_K (c_K exp (iKr) ) (2)
> >
> >
> > We expand the exp (iKr) = 4 pi sum_LM (i**L j_L(Kr) Y_LM(K) * Y_LM(r) (3)
> >
> > we set equ (1) = equ (2); put (3) into (2); and multiply this equality by Y_LM(r)
> > and integrate over the angles.
> > Then from (1) only one LM remains, i.e. only V_LM(r) remains (thats the coefficient you are looking for) and
> >
> > V_LM(r) = 4 pi/vol sum_K [ c_K i**L j_L(Kr) Y_LM(K) ]
> >
> > This is (besides the constants) CVOUT(jj,LM1)=CVOUT(jj,LM1) + POTK(J)*BES(L)*YKA(LMMULT1)
> >
> > Depending how you later on "use" the calculated quantities, it may need an additional
> > r**2, or a factor of 1/sqrt(4pi), which is equal to Y_00 ! or similar....
> >
> > and this is probably added in the SURFIN line.
> >
> > Best regards
> > Peter
> >
> >
> > Am 30.09.2012 16:03, schrieb John Rundgren:
> > > Dear Professor Blaha,
> > >
> > > I should be most grateful if you would criticize the following attempt
> > > to extend vcoul beyond RMT. Below are my modifications following your
> > > advice of 09/14/2012. TiO2 in UG is used for testing.
> > >
> > > 1) External radii > RMT are defined:
> > > rext=(/2.6d0,2.2d0/) !for TiO2.
> > >
> > > do JATOM=NSTART(myid),NSTOP(myid)
> > > do J=1,jri(jatom)+100
> > > R(J)=R0(JATOM)*EXP(dble(J-1)*DX(JATOM))
> > > if(r(j)>rext(jatom))then
> > > jrx(jatom)=j-1
> > > exit
> > > endif
> > > enddo
> > > enddo
> > > !cvout depends on jrx(jatom) below.
> > > allocate(cvout(maxval(jrx),ncom+3))
> > >
> > > 2) The part of code where exp(iKr) is expanded in Bessel functions.
> > > CVOUT=(0.d0,0.d0)
> > >
> > > do jj=jri(jatom),jrx(jatom)
> > > r(jj)=r0(jatom)*exp(dx(jatom)*dble(jj-1))
> > >
> > > DO J=2,NKK
> > > ARG=r(jj)*ABSK(J) !RMT(JATOM) replaced by r(jj).
> > > CALL SPHBES(lmax2+1,ARG,BES)
> > > CALL YKAV(J,JATOM,YKA,LMMTMX,LMMULT,LLMM)
> > > LMMULT1=0
> > >
> > > DO LM1=1,LLMM
> > > L=IABS( LM(1,LM1,JATOM) )
> > > M=LM(2,LM1,JATOM)
> > > LMMULT1=LMMULT1+1
> > > IF (LMMULT1.NE.1) THEN
> > > IF ((L.EQ.IABS(LMMULT(1,LMMULT1-1,JATOM))).AND. &
> > > (M.EQ.LMMULT(2,LMMULT1-1,JATOM))) THEN
> > > LMMULT1=LMMULT1-1
> > > endif
> > > endif
> > > CVOUT(jj,LM1)=CVOUT(jj,LM1) + POTK(J)*BES(L)*YKA(LMMULT1)
> > > ENDDO !LM1
> > >
> > > IF(SWITCH.EQ.'COUL') call
> > > rhopw(rhok(j),yka,ABSK(J),rhopw0r,jrj,rhopw0)
> > >
> > > ENDDO !J
> > > enddo !jj
> > > jrj=jrx(jatom)
> > >
> > > !cvout = constant in interval [R0,RMT] as in original code.
> > > do lm1=1,llmm
> > > do jj=1,jri(jatom)-1
> > > cvout(jj,lm1)=cvout(jri(jatom),lm1)
> > > enddo
> > > enddo
> > >
> > > 3) Statement
> > > CVOUT(lm1)=CVOUT(lm1)/IMAG1*SQRT2
> > > is replaced by
> > > CVOUT(1:jrx(jatom),lm1)=CVOUT(1:jrx(jatom),lm1)/IMAG1*SQRT2
> > >
> > > 4) Statement SURFIN=... is modified. Is **L due to a Legendre expansion?
> > > if(r(irad)<=rmt(jatom))then
> > > SURFIN=CVOUT(irad,LM1)*(-4.D0*PI)*(R(IRAD)/RMT(JATOM))**L* &
> > > IMAG**L*VOUTF
> > > elseif(r(irad)>rmt(jatom))then
> > > SURFIN=CVOUT(irad,LM1)*(-4.D0*PI)*(RMT(JATOM)/R(IRAD))**L* &
> > > IMAG**L*VOUTF
> > > endif
> > >
> > > The modifications are meant to be valid for all LM's. Please, see the
> > > Attachment depicting TiO2_R.vcoul in Rydberg units, LM=0,0. The test
> > > indicates that "run SCF" is not perturbed by the modifications made.
> > >
> > > I wonder what your decision will be about the modifications.
> > > Best regards,
> > > John Rundgren
> > >
> > >
> > >
> > > _______________________________________________
> > > Wien mailing list
> > > Wien at zeus.theochem.tuwien.ac.at
> > > http://zeus.theochem.tuwien.ac.at/mailman/listinfo/wien
> > >
> >
>
>
> _______________________________________________
> Wien mailing list
> Wien at zeus.theochem.tuwien.ac.at
> http://zeus.theochem.tuwien.ac.at/mailman/listinfo/wien
-------------- next part --------------
A non-text attachment was scrubbed...
Name: VC-setrmt-UG.ps
Type: application/postscript
Size: 18149 bytes
Desc: not available
URL: <http://zeus.theochem.tuwien.ac.at/pipermail/wien/attachments/20121002/401dadda/attachment.ps>
More information about the Wien
mailing list