private void sety2 ( double [ ] x , double [ ] y ) { double p , qn , sig , un ; double [ ] u = new double [ n - 1 ] ; y2 [ 0 ] = u [ 0 ] = 0.0 ; for ( int i = 1 ; i < n - 1 ; i ++ ) { sig = ( x [ i ] - x [ i - 1 ] ) / ( x [ i + 1 ] - x [ i - 1 ] ) ; p = sig * y2 [ i - 1 ] + 2.0 ; y2 [ i ] = ( sig - 1.0 ) / p ; u [ i ] = ( y [ i + 1 ] - y [ i ] ) / ( x [ i + 1 ] - x [ i ] ) - ( y [ i ] - y [ i - 1 ] ) / ( x [ i ] - x [ i - 1 ] ) ; u [ i ] = ( 6.0 * u [ i ] / ( x [ i + 1 ] - x [ i - 1 ] ) - sig * u [ i - 1 ] ) / p ; } qn = un = 0.0 ; y2 [ n - 1 ] = ( un - qn * u [ n - 2 ] ) / ( qn * y2 [ n - 2 ] + 1.0 ) ; for ( int k = n - 2 ; k >= 0 ; k -- ) { y2 [ k ] = y2 [ k ] * y2 [ k + 1 ] + u [ k ] ; } }
Calculate the second derivatives of the interpolating function at the tabulated points . At the endpoints we use a natural spline with zero second derivative on that boundary .