!-------------------------------------------------------------------
      ! function  : vli3_at
      !
      ! package   : Li3 Potential
      !
      ! Language  : Fortran 77
      !
      ! author    : F. Colavecchia (flavioc@lanl.gov)
      !           : G. A. Parker   (parker@phyast.nhn.ou.edu)
      !
      ! date      : 05/12/02
      ! revision  :
      !
      ! purpose   : compute Li3 Axilrod Teller term , eq. (14) of paper JCP 118 5484 (2003).
      !
      ! input     :    r    -> (dim=3) interatomic distances in bohr.
      !
      ! output    :    vli3_at -> li3 Axilrod Teller potential in a.u.  
      !
      !-------------------------------------------------------------------
      real*8 function vli3_at(r)
      implicit none
      !
      ! All distances in atomic units
      !
      real*8 r(3),cos1,cos2,cos3
      real*8 C,v,rp,coef,rp3,angfac
      parameter (coef=10000d0)
      
      vli3_at = 1.0d0/coef
      C = 170100d0
      cos1 = 0d0
      cos2 = 0d0
      cos3 = 0d0
      rp = r(1)*r(2)*r(3)
      rp3=rp**3
      if(rp.lt.1d-10) then
          if(abs(r(1)).lt.1d-10) cos1 = 1d0
          if(abs(r(2)).lt.1d-10) cos2 = 1d0
          if(abs(r(3)).lt.1d-10) cos3 = 1d0            
      else 
          cos1 = (r(2)**2+r(3)**2-r(1)**2)/(2*r(2)*r(3))
          cos2 = (r(1)**2+r(3)**2-r(2)**2)/(2*r(1)*r(3))
          cos3 = (r(2)**2+r(1)**2-r(3)**2)/(2*r(2)*r(1))
      end if

      
      angfac = (1d0+3d0*cos1*cos2*cos3)
      vli3_at=C/(coef*C+rp3)*angfac
      
      return
      end