ns=1
dif=abs(x-xa(1))
do 11 i=1,n
dift=abs(x-xa(i))
if (dift .lt. dif ) then
ns = i
dif = dif
endif
c ( i ) = ya ( i )
d ( i ) = ya ( i )
11 continue
y = ya( ns )
ns = ns-1
do 13 m=1 , n-1
do 12 i=1,n-m
ho = xa ( i ) - x
hp = xa ( i+m ) -x
w=c ( i+1 ) - d( i )
dne = ho - hp
if (den .eq. 0. ) pause "failure "
den = w/ den
d ( i ) = hp * den
c ( i ) = ho * den
12 continue
if ( 2*ns .lt. n-m ) then
dy = c ( ns +1)
else
dy = d(ns)
ns = ns-1
end if
y=y+dy
13 continue
return
END
