Dear Gavin, Great idea to create asapot.f. Many thanks. A modified subroutine also extends clm: asa_rho_pot.f where CVOUTx statements are accompanied by RHOx statements in the following way, allocate(RHOx(jrx(JATOM),ncom+3),CVOUTx(jrx(JATOM),ncom+3)) RHOx=0.d0 RHOx(jj,1)=RHO(1) CVOUTx=(0.d0,0.d0)
RHOx(jj,LM1)=RHOx(jj,LM1) + RHOK(J)*BES(L)*YKA(LMMULT1) CVOUTx(jj,LM1)=CVOUTx(jj,LM1) + POTK(J)*BES(L)*YKA(LMMULT1) clm(jj,LM1,JATOM)=pi*(4.d0*pi)*RHOx(jj,LM1)*IMAG**L V(jj,LM1,JATOM)=(4.d0*pi)*CVOUTx(jj,LM1)*IMAG**L*VOUTF deallocate(CVOUTx,RHOx) Value and gradient of rho4pir2=clm(jj,1,JATOM) look perfectly continuous. I have not tested RHOx(jj,LM1) with LM>1. Another problem is to choose rext(JATOM) in an intrinsic way. The following two conditions work impeccably together with Mattheiss's superposition method and LEED applications: Q(JATOM) = sum_R0^R(jj) clm(jj,1,JATOM) grad = V(jj,1,JATOM)-V(jj-1,1,JATOM) if(Q(JATOM)>=Z(JATOM).or.grad<0.d0)then rext(JATOM)=R(jj) exit endif Since the ASA sphere donates electrons to the interstice, it is probable that Q(JATOM) is <Z(JATOM). It is unlikely that V bends over into another atom. Best regards, John On Tue, 2012-10-09 at 13:33 -0600, Gavin Abo wrote: > Thanks Prof. Blaha, after your explanation, I was able to reproduce > John's plot with the attached subroutine called from lapw0. > _______________________________________________ > Wien mailing list > Wien at zeus.theochem.tuwien.ac.at > http://zeus.theochem.tuwien.ac.at/mailman/listinfo/wien