Program 2 : int1d1, d1int1, spd, fd1, fdm

      program prog2 
      implicit none
      
      real cubic,rms
      external cubic,rms  
      real f(10),x(10),h(10),xi(15)
      real a(10),b(10),c(10)
      real fi(15),s(15),sslp(15)
      real sslp1,sslpn
      integer i

# initialize f 
#   define a vector of arbitrary spaced coordinates, x 
#   compute f at these coordinates
#   compute the spacing of coordinates
      do i=1,10
         x(i)=i**2
         f(i)=cubic(x(i))
         print *,'f(',x(i),' )=',f(i)
         if (i.ne.1) h(i-1)=x(i)-x(i-1) 
      enddo
      print *,''

# initialize coordinates to be interpolated, xi,
# and true value of f at these coordinates
      do i=1,15
         xi(i)=15+0.3*i**2
         fi(i)=cubic(xi(i))
      enddo

# compute spline slope at each end
      call fd1(sslp1,f,h)   
      call fdm(sslpn,f,h,10)

      print *,'On an interpolation of 15 values,'
# two ways for interpolating f at coordinates xi:
# first way: with spd and d1int1
#   compute end slopes of the spline estimating f
#   compute slopes at all points
      call spd(sslp,f,10,h,0,sslp1,0,sslpn,a,b,c)      
      call d1int1(s,f,xi,x,sslp,h,10,15)
      print *,'rms error of first way (spd and d1int1) is:',rms(s,fi,15)

# second way: with int1d1
      call int1d1(s,f,xi,x,sslp,h,10,15,0,sslp1,0,sslpn,a,b,c)
      print *,'rms error of the second way (int1d1) is:',rms(s,fi,15)

      stop
      end


      real function cubic(x)
      real x
      cubic=x**3+15*x**2+60*x+6
      return
      end
         
      real function rms(a,b,n)
      real a(n),b(n),sum
      integer i,n
      sum=0
      do i=1,n
         sum=sum+(a(i)-b(i))**2
      enddo
      rms=sqrt(sum/n)
      return
      end