CTaylor – a C++ template library for automatic differentiation

This is a library to calculate with truncated taylor series (instead of double values). It was inspired by Ulf Ekstrӧm’s libtaylor but has been extended to automatically grow the order of involved derivatives. One simply creates an independent variable by passing an enumeration and the maximum order of derivative one wants to calculate. Such an independent variable only carries the value and the first derivative value (1). This library is using meta programming (using the boost::mpl library) to create different types matching the different number of independent variables involved and different order. Certain nonlinear operations increase the order of the involved derivatives. E.g. multiplying two different independent variables (e.g. x and y) with each other will result in a value containing derivatives vs. d/dx and d/dy and d^2/dx/dy. The library does not store derivatives but polynomial coefficients, which differ from derivatives by faculty values.

Here is an example piece of source code and here is the resulting executable (for windows 64). The executable needs some command line parameters (this is not checked for and you will get a crash if they are not provided):

.\cder 3 1e-6 4 1e-6 5 1e-6 16384

This will give the following output:

1:2:1: = 227058, 1:2: = -205341, 2:2: = 16983, 1:1:2: = -43480.7, 1:1:1: = 68352.4, 2:1:1: = -5652.75, 3:1: = 2902.46, 1:1: = 17422.9, 2:1: = -10662.1, 1:0:2: = -5446.79, 2:0:2: = 322.691, 3:0:1: = -483.084, 1:0:1: = -2900.09, 2:0:1: = 1774.62, 1:0:3: = 3016.72, 1:3: = -454722, 4: = 29.1334, 2: = -531.299, 1: = 5099.3,  = 141734, 3: = -11.2339, 0:3: = -1.70336e+006, 0:3:1: = 2.9838e+006, 0:0:3: = 6206.25, 0:0:1: = -23393.4, 0:0:2: = 10542.2, 0:0:4: = -5154.73, 0:1:3: = 105264, 0:1: = 140554, 0:1:1: = -103281, 0:1:2: = -132955, 0:2: = 310248, 0:2:1: = 850492, 0:2:2: = -815717, 0:4: = -4.48177e+006,

The performance of this library compared to libtaylor is considerable higher, as a lot of zero derivatives are not being calculated and stored.

The digits before the equal sign are order of derivatives for the first, second or third variable. To find the value look for an empty field (“,  = 141734,”). Thus “1:2:1: = 227058” is represents d^4/dx/dy^2/dz.

Consider that the output does not represent the plain derivatives for some point, but the sum of such function and their derivatives:

Some hints for source code to be mailed in:

·         Use the c++0x auto and decltype keywords.

·         Use templates for every function parameter, which may be of type CTaylor.

·         Wrong:
static double myFunction(const double _dX, const double _dY)
{              return std::sin(_dX*_dX + _dY*_dY);
}

·         Correct:
template<typename TX, typename TY>
inline auto myFunction(const TX &_dX, const TY &_dY) -> decltype(std::sin(_dX*_dX + _dY*_dY))
{              return std::sin(_dX*_dX + _dY*_dY);
}

·         Wrong:
template<typename TX, typename TY>
inline auto myFunction(const TX &_dX, const TY &_dY) -> decltype(std::sin(_dX*_dX*_dX + _dY*_dY*_dY))
{              TX X3 = _dX*_dX*_dX;
TY Y3 = _dY*_dY*_dY;
return std::sin(X3 + Y3);
}

·         Correct:
template<typename TX, typename TY>
inline auto myFunction(const TX &_dX, const TY &_dY) -> decltype(std::sin(_dX*_dX*_dX + _dY*_dY*_dY))
{              auto X3 = _dX*_dX*_dX;
auto Y3 = _dY*_dY*_dY;
return std::sin(X3 + Y3);
}

·         Reusing variables: Be careful when reusing variables in different assignments, as the type may have to change. In the code below the assignment b=c*c may not compile anymore.
decltype(a*a) b;
if (a < 0)
b = 0;
else
b = a*a;
// use b:
auto c = b*b - d;
if (c < 0)
b = 0;
else
b = c*c;

·         Use the boost::remove_const class to remove const specifiers from types:
const auto b = a*a;
boost::remove_const<decltype(b*b)>::type c = b*b;

·         In case of your code contains if-statements, which select different functions depending on the value of some input variable (or some value influenced by some independent variable), then your code is likely not smooth in all derivatives. It only makes sense to calculate derivatives of an order, which are still smooth. There are some methods to make step-wise equations smooth for all derivatives.

mailto:ctaylor@foelsche.com