Prev Next

dcsevl

Prototype
doublereal f2cad::dcsevl_(doublereal *x, doublereal *a, integer *n);

Fortran Source
dcsevl.f

Description
The routine dcsevl computes values for a Chebyshev Series; i.e., the return value is  \[
     f(x) = a_0 / 2 + \sum_{i=1}^{n-1} a_i T_i (x)
\]  where  T_i (x) is the i-th order Chebyshev polynomial of the first kind; i.e., it is defined by the recurrence relation  \[
\begin{array}{rcl}
     T_0 (x)     & = & 1 \\
     T_1 (x)     & = & x \\
     T_{i+1} (x) & = & 2 x T_i (x) - T_{i-1} (x)
\end{array}
\]  It follows that the derivatives  T_i^{(1)} (x) satisfy the recurrence relation  \[
\begin{array}{rcl}
     T_0^{(1)} (x)     & = & 0 \\
     T_1^{(1)} (x)     & = & 1 \\
     T_{i+1}^{(1)} (x) & = & 2 T_i (x) + 2 x T_i^{(1)} (x) - T_{i-1}^{(1)} (x)
\end{array}
\]  and the derivative of  f(x) is given by  \[
     f^{(1)} (x) = \sum_{i=1}^{n-1} a_i T_i^{(1)} (x)
\]  
 

# include <f2cad/dcsevl.hpp>

test_result dcsevl(void)
{	bool ok = true;
	// independent variable
	integer  N = 1;
	double   x = .5;
	doublereal X[1];
	X[0]       = x;
	// dependent variable
	integer  M = 1;
	doublereal F[1];
	// polynomial coefficients
	integer n  = 3;
	double     a[3];
	doublereal A[3];
	for(integer i = 0; i < n; i++)
	{	a[i] = 1. / double(i+1);
		A[i] = a[i];
	}

	// declare independent variables
	f2cad::Independent(N, X);

	F[0] = f2cad::dcsevl_(X, A, &n);

	// declare dependent variables
	f2cad::Dependent(M, F);

	// evaluate the Chebyshev polynomials
	double t_0 = 1.;
	double t_1 = x;
	double t_2 = 2. * x * t_1 - t_0;

	// check function value
	double value = f2cad::Value( F[0] );
	double check = a[0] * t_0 / 2. + a[1] * t_1 + a[2] * t_2;
	ok          &= f2cad::near_equal(check, value, 1e-10, 1e-10);

	// evaluate the Derivatives
	double dt_0 = 0.;
	double dt_1 = 1.;
	double dt_2 = 2. * t_1 + 2. * x * dt_1 - dt_0; 
	value       = f2cad::Partial<doublereal>(0, 0);
	check       = a[1] * dt_1 + a[2] * dt_2;
	ok          &= f2cad::near_equal(check, value, 1e-10, 1e-10);

	if( ok )
		return test_pass;
	return test_fail;
}
Input File: example/dcsevl.cpp