Re: FORTRAN TYPE construct?



I am in fact using the F77 SPLINE and SPLINT. My code only populates
the splinegroups and then calls SPLINE, verifies that there are the
correct values in the arrays, then calls SPLINT.

I added the same code to output the array contents into the SPLINT
subroutine and found that upon calling SPLINT, all the arrays were
empty (though their contents printed fine just before calling SPLINT),
and after exiting the SPLINT subroutine, the arrays were empty.


glen herrmannsfeldt wrote:
michael@xxxxxxxxxxxxxxxx wrote:
Andrew:
Cubic Spline interpolation routines (it is better to use Shape
preserving Cubic Splines ), in general require two arrays a(n) and
b(4,n) that hold the calculated coefficients of the splines. N is the
number of points you wish to interpoalte.
So make sure that you allocate enough space to hold these coefficients

That was my thought, too, so I looked in my Numerical Recipes.

It seems that NR has routines named SPLINE and SPLINT, though I
only have the original (Fortran 77) version of the book.
Instead of returning the four spline coefficients, they return the
second derivative of the spline function at each point.
It seems to be easy to interpolate the cubic spline value from
the Y and Y'' values at the two nearest points.

The Fortran 77 version has a call like:

CALL SPLINE(X,Y,N,YP1,YPN,Y2)

with X, Y, and Y2 dimensioned N. YP1 and YPN are the input
first derivative at the end points, or 1e30 to set the second
derivative to zero at the end for a natural spline.

A newer version presumably doesn't need N, but should still need
YP1 and YPN.

-- glen

.



Relevant Pages


Loading