anchor.diff:  
     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)
  
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:       enddo
        intr=ir-1
        if (intr.lt.0.or.intr.ge.NRps(ik))stop 'bad intr at prep_ps'
        xx = r - radnl(intr,ik) 
        do l=0,Mlps(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)
  