c
c
c
      function rint(r)
C
C     H2O            INTERPOLATIONS POTENTIAL
C
      implicit real*8 (a-h,o-z)
      real rlgdre*8
C
      dimension r(3),gamma(5),a(5,5),b(5,6),p(5),x(6)
C
      data d0,alpha0,r0/4.621d0,2.294d0,0.971d0/
      data gamma/45.0d0,75.0d0,104.52d0,135.0d0,180.0d0/
      data b/-0.69629d0,-0.48884d0,-0.45792d0,-0.44409d0,-0.41761d0,
     *       1.40886d0,1.39153d0,1.39691d0,1.46819d0,1.54303d0,
     *       0.12042d0,0.14350d0,0.17865d0,0.19215d0,0.19151d0,
     *       0.19193d0,0.15280d0,0.13400d0,0.09861d0,0.05140d0,
     *       0.62717d0,1.38163d0,1.43055d0,1.36016d0,1.22214d0,
     *       0.31021d0,0.12674d0,0.18001d0,0.38321d0,0.51142d0/
C
      logical nfirst
      data nfirst/.true./
C
      if (nfirst) then
      nfirst=.false.
C
C     BESTIMMUNG DER LEGENDRE ENTWICKLUNGS KOEFF.
C
      do 1 i=1,5
      do 1 j=1,5
 1    a(i,j)=rlgdre(j-1,cos(1.7453292d-2*gamma(i)))
C
      call mat56(a,5,b,6,det)
C
      end if
C
      rmax=dmax1(r(1),r(3))
      rmaxa0=rmax/0.529177249d0 - 3.5d0
      rmin=dmin1(r(1),r(3))
      cosgam=-(r(2)**2-r(1)**2-r(3)**2)/(2*r(1)*r(3))
C
      do 2 i=1,5
 2    p(i)=rlgdre(i-1,cosgam)
C
      do 3 j=1,6
      x(j)=0.0d0
      do 3 i=1,5
      x(j)=x(j)+b(i,j)*p(i)
 3    continue
C
      d=d0+x(1)*exp(-x(2)*rmaxa0-x(3)*rmaxa0**2)
      alpha=alpha0+x(4)*exp(-x(5)*rmaxa0-x(6)*rmaxa0**2)
C
C
      rint=d*(exp(-alpha*(rmin-r0))-1.0d0)**2-d
C
      return
      end