anchor.diff:        enddo
         intr=ir-1
         if (intr.lt.0.or.intr.ge.NRps(ik))stop 'bad intr at prep_ps'
  -      ratio1=(r-radnl(intr,ik))/(radnl(intr+1,ik)-radnl(intr,ik))
  -      ratio2=1-ratio1
  +      xx = r - radnl(intr,ik) 
         do l=0,Mlps(ik)
  -        uVr(l)=ratio1*udVtbl(intr+1,l,ik)+ratio2*udVtbl(intr,l,ik)
  -        duVr(l)=ratio1*dudVtbl(intr+1,l,ik)+ratio2*dudVtbl(intr,l,ik)
  +         uVr(l) = udVtbl_a(intr,l)*xx**3 + udVtbl_b(intr,l)*xx**2 &
  +                 +udVtbl_c(intr,l)*xx    + udVtbl_d(intr,l)
  +         duVr(l)=dudVtbl_a(intr,l)*xx**3 +dudVtbl_b(intr,l)*xx**2 &
  +                +dudVtbl_c(intr,l)*xx    +dudVtbl_d(intr,l)         
         enddo
         lm=0
         do l=0,Mlps(ik)
  
positive.code:         iuV(lma_tbl(lm,a))=inorm(l,ik)
        enddo
      enddo
    enddo
  
    return
  End Subroutine prep_ps_periodic
  !--------10--------20--------30--------40--------50--------60--------70--------80--------90--------100-------110-------120--------130
  
negative.code: 
    do a=1,NI
      ik=Kion(a)
  
      allocate(xn(0:NRps(ik)-1),yn(0:NRps(ik)-1),an(0:NRps(ik)-2) &
           ,bn(0:NRps(ik)-2),cn(0:NRps(ik)-2),dn(0:NRps(ik)-2))
      
      xn(0:NRps(ik)-1) = radnl(1:NRps(ik),ik)
      do l=0,Mlps(ik)
         yn(0:NRps(ik)-1) = udVtbl(1:NRps(ik),l,ik)
         call spline(NRps(ik),xn,yn,an,bn,cn,dn)
         udVtbl_a(1:NRps(ik)-1,l) = an(0:NRps(ik)-2)
         udVtbl_b(1:NRps(ik)-1,l) = bn(0:NRps(ik)-2)
         udVtbl_c(1:NRps(ik)-1,l) = cn(0:NRps(ik)-2)
         udVtbl_d(1:NRps(ik)-1,l) = dn(0:NRps(ik)-2)
  
         yn(0:NRps(ik)-1) = dudVtbl(1:NRps(ik),l,ik)
         call spline(NRps(ik),xn,yn,an,bn,cn,dn)
         dudVtbl_a(1:NRps(ik)-1,l) = an(0:NRps(ik)-2)
         dudVtbl_b(1:NRps(ik)-1,l) = bn(0:NRps(ik)-2)
         dudVtbl_c(1:NRps(ik)-1,l) = cn(0:NRps(ik)-2)
         dudVtbl_d(1:NRps(ik)-1,l) = dn(0:NRps(ik)-2)        
      end do
      
      do j=1,Mps(a)
        x=Lx(Jxyz(j,a))*Hx-(Rion(1,a)+Jxx(j,a)*aLx)
        y=Ly(Jxyz(j,a))*Hy-(Rion(2,a)+Jyy(j,a)*aLy)
  