Java Code Examples for org.apache.commons.math3.analysis.differentiation.DerivativeStructure#getValue()
The following examples show how to use
org.apache.commons.math3.analysis.differentiation.DerivativeStructure#getValue() .
You can vote up the ones you like or vote down the ones you don't like,
and go to the original project or source file by following the links above each example. You may check out the related API usage on the sidebar.
Example 1
| Source File: PolynomialSplineFunction.java From astor with GNU General Public License v2.0 | 6 votes |
|
/** {@inheritDoc}
* @since 3.1
*/
public DerivativeStructure value(final DerivativeStructure t) {
final double t0 = t.getValue();
if (t0 < knots[0] || t0 > knots[n]) {
throw new OutOfRangeException(t0, knots[0], knots[n]);
}
int i = Arrays.binarySearch(knots, t0);
if (i < 0) {
i = -i - 2;
}
// This will handle the case where t is the last knot value
// There are only n-1 polynomials, so if t is the last knot
// then we will use the last polynomial to calculate the value.
if ( i >= polynomials.length ) {
i--;
}
return polynomials[i].value(t.subtract(knots[i]));
}
Example 2
| Source File: NewtonRaphsonSolver.java From astor with GNU General Public License v2.0 | 6 votes |
|
/**
* {@inheritDoc}
*/
@Override
protected double doSolve()
throws TooManyEvaluationsException {
final double startValue = getStartValue();
final double absoluteAccuracy = getAbsoluteAccuracy();
double x0 = startValue;
double x1;
while (true) {
final DerivativeStructure y0 = computeObjectiveValueAndDerivative(x0);
x1 = x0 - (y0.getValue() / y0.getPartialDerivative(1));
if (FastMath.abs(x1 - x0) <= absoluteAccuracy) {
return x1;
}
x0 = x1;
}
}
Example 3
| Source File: MinpackTest.java From astor with GNU General Public License v2.0 | 6 votes |
|
@Override
public DerivativeStructure[] value(DerivativeStructure[] variables) {
DerivativeStructure x1 = variables[0];
DerivativeStructure x2 = variables[1];
DerivativeStructure x3 = variables[2];
DerivativeStructure tmp1 = variables[0].getField().getZero();
if (x1.getValue() == 0) {
tmp1 = tmp1.add((x2.getValue() >= 0) ? 0.25 : -0.25);
} else {
tmp1 = x2.divide(x1).atan().divide(twoPi);
if (x1.getValue() < 0) {
tmp1 = tmp1.add(0.5);
}
}
DerivativeStructure tmp2 = x1.multiply(x1).add(x2.multiply(x2)).sqrt();
return new DerivativeStructure[] {
x3.subtract(tmp1.multiply(10)).multiply(10),
tmp2.subtract(1).multiply(10),
x3
};
}
Example 4
| Source File: PolynomialSplineFunction.java From astor with GNU General Public License v2.0 | 6 votes |
|
/** {@inheritDoc}
* @since 3.1
*/
public DerivativeStructure value(final DerivativeStructure t) {
final double t0 = t.getValue();
if (t0 < knots[0] || t0 > knots[n]) {
throw new OutOfRangeException(t0, knots[0], knots[n]);
}
int i = Arrays.binarySearch(knots, t0);
if (i < 0) {
i = -i - 2;
}
// This will handle the case where t is the last knot value
// There are only n-1 polynomials, so if t is the last knot
// then we will use the last polynomial to calculate the value.
if ( i >= polynomials.length ) {
i--;
}
return polynomials[i].value(t.subtract(knots[i]));
}
Example 5
| Source File: MinpackTest.java From astor with GNU General Public License v2.0 | 6 votes |
|
@Override
public DerivativeStructure[] value(DerivativeStructure[] variables) {
DerivativeStructure x1 = variables[0];
DerivativeStructure x2 = variables[1];
DerivativeStructure x3 = variables[2];
DerivativeStructure tmp1 = variables[0].getField().getZero();
if (x1.getValue() == 0) {
tmp1 = tmp1.add((x2.getValue() >= 0) ? 0.25 : -0.25);
} else {
tmp1 = x2.divide(x1).atan().divide(twoPi);
if (x1.getValue() < 0) {
tmp1 = tmp1.add(0.5);
}
}
DerivativeStructure tmp2 = x1.multiply(x1).add(x2.multiply(x2)).sqrt();
return new DerivativeStructure[] {
x3.subtract(tmp1.multiply(10)).multiply(10),
tmp2.subtract(1).multiply(10),
x3
};
}
Example 6
| Source File: NewtonRaphsonSolver.java From astor with GNU General Public License v2.0 | 6 votes |
|
/**
* {@inheritDoc}
*/
@Override
protected double doSolve()
throws TooManyEvaluationsException {
final double startValue = getStartValue();
final double absoluteAccuracy = getAbsoluteAccuracy();
double x0 = startValue;
double x1;
while (true) {
final DerivativeStructure y0 = computeObjectiveValueAndDerivative(x0);
x1 = x0 - (y0.getValue() / y0.getPartialDerivative(1));
if (FastMath.abs(x1 - x0) <= absoluteAccuracy) {
return x1;
}
x0 = x1;
}
}
Example 7
| Source File: HarmonicOscillator.java From astor with GNU General Public License v2.0 | 6 votes |
|
/** {@inheritDoc}
* @since 3.1
*/
public DerivativeStructure value(final DerivativeStructure t) {
final double x = t.getValue();
double[] f = new double[t.getOrder() + 1];
final double alpha = omega * x + phase;
f[0] = amplitude * FastMath.cos(alpha);
if (f.length > 1) {
f[1] = -amplitude * omega * FastMath.sin(alpha);
final double mo2 = - omega * omega;
for (int i = 2; i < f.length; ++i) {
f[i] = mo2 * f[i - 2];
}
}
return t.compose(f);
}
Example 8
| Source File: NewtonRaphsonSolver.java From astor with GNU General Public License v2.0 | 6 votes |
|
/**
* {@inheritDoc}
*/
@Override
protected double doSolve()
throws TooManyEvaluationsException {
final double startValue = getStartValue();
final double absoluteAccuracy = getAbsoluteAccuracy();
double x0 = startValue;
double x1;
while (true) {
final DerivativeStructure y0 = computeObjectiveValueAndDerivative(x0);
x1 = x0 - (y0.getValue() / y0.getPartialDerivative(1));
if (FastMath.abs(x1 - x0) <= absoluteAccuracy) {
return x1;
}
x0 = x1;
}
}
Example 9
| Source File: MinpackTest.java From astor with GNU General Public License v2.0 | 6 votes |
|
@Override
public DerivativeStructure[] value(DerivativeStructure[] variables) {
DerivativeStructure x1 = variables[0];
DerivativeStructure x2 = variables[1];
DerivativeStructure x3 = variables[2];
DerivativeStructure tmp1 = variables[0].getField().getZero();
if (x1.getValue() == 0) {
tmp1 = tmp1.add((x2.getValue() >= 0) ? 0.25 : -0.25);
} else {
tmp1 = x2.divide(x1).atan().divide(twoPi);
if (x1.getValue() < 0) {
tmp1 = tmp1.add(0.5);
}
}
DerivativeStructure tmp2 = x1.multiply(x1).add(x2.multiply(x2)).sqrt();
return new DerivativeStructure[] {
x3.subtract(tmp1.multiply(10)).multiply(10),
tmp2.subtract(1).multiply(10),
x3
};
}
Example 10
| Source File: HarmonicOscillator.java From astor with GNU General Public License v2.0 | 6 votes |
|
/** {@inheritDoc}
* @since 3.1
*/
public DerivativeStructure value(final DerivativeStructure t)
throws DimensionMismatchException {
final double x = t.getValue();
double[] f = new double[t.getOrder() + 1];
final double alpha = omega * x + phase;
f[0] = amplitude * FastMath.cos(alpha);
if (f.length > 1) {
f[1] = -amplitude * omega * FastMath.sin(alpha);
final double mo2 = - omega * omega;
for (int i = 2; i < f.length; ++i) {
f[i] = mo2 * f[i - 2];
}
}
return t.compose(f);
}
Example 11
| Source File: PolynomialSplineFunction.java From astor with GNU General Public License v2.0 | 6 votes |
|
/** {@inheritDoc}
* @since 3.1
*/
public DerivativeStructure value(final DerivativeStructure t) {
final double t0 = t.getValue();
if (t0 < knots[0] || t0 > knots[n]) {
throw new OutOfRangeException(t0, knots[0], knots[n]);
}
int i = Arrays.binarySearch(knots, t0);
if (i < 0) {
i = -i - 2;
}
// This will handle the case where t is the last knot value
// There are only n-1 polynomials, so if t is the last knot
// then we will use the last polynomial to calculate the value.
if ( i >= polynomials.length ) {
i--;
}
return polynomials[i].value(t.subtract(knots[i]));
}
Example 12
| Source File: MinpackTest.java From astor with GNU General Public License v2.0 | 6 votes |
|
@Override
public DerivativeStructure[] value(DerivativeStructure[] variables) {
DerivativeStructure x1 = variables[0];
DerivativeStructure x2 = variables[1];
DerivativeStructure x3 = variables[2];
DerivativeStructure tmp1 = variables[0].getField().getZero();
if (x1.getValue() == 0) {
tmp1 = tmp1.add((x2.getValue() >= 0) ? 0.25 : -0.25);
} else {
tmp1 = x2.divide(x1).atan().divide(twoPi);
if (x1.getValue() < 0) {
tmp1 = tmp1.add(0.5);
}
}
DerivativeStructure tmp2 = x1.multiply(x1).add(x2.multiply(x2)).sqrt();
return new DerivativeStructure[] {
x3.subtract(tmp1.multiply(10)).multiply(10),
tmp2.subtract(1).multiply(10),
x3
};
}
Example 13
| Source File: NewtonRaphsonSolver.java From astor with GNU General Public License v2.0 | 6 votes |
|
/**
* {@inheritDoc}
*/
@Override
protected double doSolve()
throws TooManyEvaluationsException {
final double startValue = getStartValue();
final double absoluteAccuracy = getAbsoluteAccuracy();
double x0 = startValue;
double x1;
while (true) {
final DerivativeStructure y0 = computeObjectiveValueAndDerivative(x0);
x1 = x0 - (y0.getValue() / y0.getPartialDerivative(1));
if (FastMath.abs(x1 - x0) <= absoluteAccuracy) {
return x1;
}
x0 = x1;
}
}
Example 14
| Source File: Sinc.java From astor with GNU General Public License v2.0 | 4 votes |
|
/** {@inheritDoc}
* @since 3.1
*/
public DerivativeStructure value(final DerivativeStructure t)
throws DimensionMismatchException {
final double scaledX = (normalized ? FastMath.PI : 1) * t.getValue();
final double scaledX2 = scaledX * scaledX;
double[] f = new double[t.getOrder() + 1];
if (FastMath.abs(scaledX) <= SHORTCUT) {
for (int i = 0; i < f.length; ++i) {
final int k = i / 2;
if ((i & 0x1) == 0) {
// even derivation order
f[i] = (((k & 0x1) == 0) ? 1 : -1) *
(1.0 / (i + 1) - scaledX2 * (1.0 / (2 * i + 6) - scaledX2 / (24 * i + 120)));
} else {
// odd derivation order
f[i] = (((k & 0x1) == 0) ? -scaledX : scaledX) *
(1.0 / (i + 2) - scaledX2 * (1.0 / (6 * i + 24) - scaledX2 / (120 * i + 720)));
}
}
} else {
final double inv = 1 / scaledX;
final double cos = FastMath.cos(scaledX);
final double sin = FastMath.sin(scaledX);
f[0] = inv * sin;
// the nth order derivative of sinc has the form:
// dn(sinc(x)/dxn = [S_n(x) sin(x) + C_n(x) cos(x)] / x^(n+1)
// where S_n(x) is an even polynomial with degree n-1 or n (depending on parity)
// and C_n(x) is an odd polynomial with degree n-1 or n (depending on parity)
// S_0(x) = 1, S_1(x) = -1, S_2(x) = -x^2 + 2, S_3(x) = 3x^2 - 6...
// C_0(x) = 0, C_1(x) = x, C_2(x) = -2x, C_3(x) = -x^3 + 6x...
// the general recurrence relations for S_n and C_n are:
// S_n(x) = x S_(n-1)'(x) - n S_(n-1)(x) - x C_(n-1)(x)
// C_n(x) = x C_(n-1)'(x) - n C_(n-1)(x) + x S_(n-1)(x)
// as per polynomials parity, we can store both S_n and C_n in the same array
final double[] sc = new double[f.length];
sc[0] = 1;
double coeff = inv;
for (int n = 1; n < f.length; ++n) {
double s = 0;
double c = 0;
// update and evaluate polynomials S_n(x) and C_n(x)
final int kStart;
if ((n & 0x1) == 0) {
// even derivation order, S_n is degree n and C_n is degree n-1
sc[n] = 0;
kStart = n;
} else {
// odd derivation order, S_n is degree n-1 and C_n is degree n
sc[n] = sc[n - 1];
c = sc[n];
kStart = n - 1;
}
// in this loop, k is always even
for (int k = kStart; k > 1; k -= 2) {
// sine part
sc[k] = (k - n) * sc[k] - sc[k - 1];
s = s * scaledX2 + sc[k];
// cosine part
sc[k - 1] = (k - 1 - n) * sc[k - 1] + sc[k -2];
c = c * scaledX2 + sc[k - 1];
}
sc[0] *= -n;
s = s * scaledX2 + sc[0];
coeff *= inv;
f[n] = coeff * (s * sin + c * scaledX * cos);
}
}
if (normalized) {
double scale = FastMath.PI;
for (int i = 1; i < f.length; ++i) {
f[i] *= scale;
scale *= FastMath.PI;
}
}
return t.compose(f);
}
Example 15
| Source File: Gaussian.java From astor with GNU General Public License v2.0 | 4 votes |
|
/** {@inheritDoc}
* @since 3.1
*/
public DerivativeStructure value(final DerivativeStructure t)
throws DimensionMismatchException {
final double u = is * (t.getValue() - mean);
double[] f = new double[t.getOrder() + 1];
// the nth order derivative of the Gaussian has the form:
// dn(g(x)/dxn = (norm / s^n) P_n(u) exp(-u^2/2) with u=(x-m)/s
// where P_n(u) is a degree n polynomial with same parity as n
// P_0(u) = 1, P_1(u) = -u, P_2(u) = u^2 - 1, P_3(u) = -u^3 + 3 u...
// the general recurrence relation for P_n is:
// P_n(u) = P_(n-1)'(u) - u P_(n-1)(u)
// as per polynomial parity, we can store coefficients of both P_(n-1) and P_n in the same array
final double[] p = new double[f.length];
p[0] = 1;
final double u2 = u * u;
double coeff = norm * FastMath.exp(-0.5 * u2);
if (coeff <= Precision.SAFE_MIN) {
Arrays.fill(f, 0.0);
} else {
f[0] = coeff;
for (int n = 1; n < f.length; ++n) {
// update and evaluate polynomial P_n(x)
double v = 0;
p[n] = -p[n - 1];
for (int k = n; k >= 0; k -= 2) {
v = v * u2 + p[k];
if (k > 2) {
p[k - 2] = (k - 1) * p[k - 1] - p[k - 3];
} else if (k == 2) {
p[0] = p[1];
}
}
if ((n & 0x1) == 1) {
v *= u;
}
coeff *= is;
f[n] = coeff * v;
}
}
return t.compose(f);
}
Example 16
| Source File: Logit.java From astor with GNU General Public License v2.0 | 4 votes |
|
/** {@inheritDoc}
* @since 3.1
* @exception OutOfRangeException if parameter is outside of function domain
*/
public DerivativeStructure value(final DerivativeStructure t)
throws OutOfRangeException {
final double x = t.getValue();
if (x < lo || x > hi) {
throw new OutOfRangeException(x, lo, hi);
}
double[] f = new double[t.getOrder() + 1];
// function value
f[0] = FastMath.log((x - lo) / (hi - x));
if (Double.isInfinite(f[0])) {
if (f.length > 1) {
f[1] = Double.POSITIVE_INFINITY;
}
// fill the array with infinities
// (for x close to lo the signs will flip between -inf and +inf,
// for x close to hi the signs will always be +inf)
// this is probably overkill, since the call to compose at the end
// of the method will transform most infinities into NaN ...
for (int i = 2; i < f.length; ++i) {
f[i] = f[i - 2];
}
} else {
// function derivatives
final double invL = 1.0 / (x - lo);
double xL = invL;
final double invH = 1.0 / (hi - x);
double xH = invH;
for (int i = 1; i < f.length; ++i) {
f[i] = xL + xH;
xL *= -i * invL;
xH *= i * invH;
}
}
return t.compose(f);
}
Example 17
| Source File: Logit.java From astor with GNU General Public License v2.0 | 4 votes |
|
/** {@inheritDoc}
* @since 3.1
* @exception OutOfRangeException if parameter is outside of function domain
*/
public DerivativeStructure value(final DerivativeStructure t)
throws OutOfRangeException {
final double x = t.getValue();
if (x < lo || x > hi) {
throw new OutOfRangeException(x, lo, hi);
}
double[] f = new double[t.getOrder() + 1];
// function value
f[0] = FastMath.log((x - lo) / (hi - x));
if (Double.isInfinite(f[0])) {
if (f.length > 1) {
f[1] = Double.POSITIVE_INFINITY;
}
// fill the array with infinities
// (for x close to lo the signs will flip between -inf and +inf,
// for x close to hi the signs will always be +inf)
// this is probably overkill, since the call to compose at the end
// of the method will transform most infinities into NaN ...
for (int i = 2; i < f.length; ++i) {
f[i] = f[i - 2];
}
} else {
// function derivatives
final double invL = 1.0 / (x - lo);
double xL = invL;
final double invH = 1.0 / (hi - x);
double xH = invH;
for (int i = 1; i < f.length; ++i) {
f[i] = xL + xH;
xL *= -i * invL;
xH *= i * invH;
}
}
return t.compose(f);
}
Example 18
| Source File: Gaussian.java From astor with GNU General Public License v2.0 | 4 votes |
|
/** {@inheritDoc}
* @since 3.1
*/
public DerivativeStructure value(final DerivativeStructure t)
throws DimensionMismatchException {
final double u = is * (t.getValue() - mean);
double[] f = new double[t.getOrder() + 1];
// the nth order derivative of the Gaussian has the form:
// dn(g(x)/dxn = (norm / s^n) P_n(u) exp(-u^2/2) with u=(x-m)/s
// where P_n(u) is a degree n polynomial with same parity as n
// P_0(u) = 1, P_1(u) = -u, P_2(u) = u^2 - 1, P_3(u) = -u^3 + 3 u...
// the general recurrence relation for P_n is:
// P_n(u) = P_(n-1)'(u) - u P_(n-1)(u)
// as per polynomial parity, we can store coefficients of both P_(n-1) and P_n in the same array
final double[] p = new double[f.length];
p[0] = 1;
final double u2 = u * u;
double coeff = norm * FastMath.exp(-0.5 * u2);
if (coeff <= Precision.SAFE_MIN) {
Arrays.fill(f, 0.0);
} else {
f[0] = coeff;
for (int n = 1; n < f.length; ++n) {
// update and evaluate polynomial P_n(x)
double v = 0;
p[n] = -p[n - 1];
for (int k = n; k >= 0; k -= 2) {
v = v * u2 + p[k];
if (k > 2) {
p[k - 2] = (k - 1) * p[k - 1] - p[k - 3];
} else if (k == 2) {
p[0] = p[1];
}
}
if ((n & 0x1) == 1) {
v *= u;
}
coeff *= is;
f[n] = coeff * v;
}
}
return t.compose(f);
}
Example 19
| Source File: Sinc.java From astor with GNU General Public License v2.0 | 4 votes |
|
/** {@inheritDoc}
* @since 3.1
*/
public DerivativeStructure value(final DerivativeStructure t)
throws DimensionMismatchException {
final double scaledX = (normalized ? FastMath.PI : 1) * t.getValue();
final double scaledX2 = scaledX * scaledX;
double[] f = new double[t.getOrder() + 1];
if (FastMath.abs(scaledX) <= SHORTCUT) {
for (int i = 0; i < f.length; ++i) {
final int k = i / 2;
if ((i & 0x1) == 0) {
// even derivation order
f[i] = (((k & 0x1) == 0) ? 1 : -1) *
(1.0 / (i + 1) - scaledX2 * (1.0 / (2 * i + 6) - scaledX2 / (24 * i + 120)));
} else {
// odd derivation order
f[i] = (((k & 0x1) == 0) ? -scaledX : scaledX) *
(1.0 / (i + 2) - scaledX2 * (1.0 / (6 * i + 24) - scaledX2 / (120 * i + 720)));
}
}
} else {
final double inv = 1 / scaledX;
final double cos = FastMath.cos(scaledX);
final double sin = FastMath.sin(scaledX);
f[0] = inv * sin;
// the nth order derivative of sinc has the form:
// dn(sinc(x)/dxn = [S_n(x) sin(x) + C_n(x) cos(x)] / x^(n+1)
// where S_n(x) is an even polynomial with degree n-1 or n (depending on parity)
// and C_n(x) is an odd polynomial with degree n-1 or n (depending on parity)
// S_0(x) = 1, S_1(x) = -1, S_2(x) = -x^2 + 2, S_3(x) = 3x^2 - 6...
// C_0(x) = 0, C_1(x) = x, C_2(x) = -2x, C_3(x) = -x^3 + 6x...
// the general recurrence relations for S_n and C_n are:
// S_n(x) = x S_(n-1)'(x) - n S_(n-1)(x) - x C_(n-1)(x)
// C_n(x) = x C_(n-1)'(x) - n C_(n-1)(x) + x S_(n-1)(x)
// as per polynomials parity, we can store both S_n and C_n in the same array
final double[] sc = new double[f.length];
sc[0] = 1;
double coeff = inv;
for (int n = 1; n < f.length; ++n) {
double s = 0;
double c = 0;
// update and evaluate polynomials S_n(x) and C_n(x)
final int kStart;
if ((n & 0x1) == 0) {
// even derivation order, S_n is degree n and C_n is degree n-1
sc[n] = 0;
kStart = n;
} else {
// odd derivation order, S_n is degree n-1 and C_n is degree n
sc[n] = sc[n - 1];
c = sc[n];
kStart = n - 1;
}
// in this loop, k is always even
for (int k = kStart; k > 1; k -= 2) {
// sine part
sc[k] = (k - n) * sc[k] - sc[k - 1];
s = s * scaledX2 + sc[k];
// cosine part
sc[k - 1] = (k - 1 - n) * sc[k - 1] + sc[k -2];
c = c * scaledX2 + sc[k - 1];
}
sc[0] *= -n;
s = s * scaledX2 + sc[0];
coeff *= inv;
f[n] = coeff * (s * sin + c * scaledX * cos);
}
}
if (normalized) {
double scale = FastMath.PI;
for (int i = 1; i < f.length; ++i) {
f[i] *= scale;
scale *= FastMath.PI;
}
}
return t.compose(f);
}
Example 20
| Source File: Sinc.java From astor with GNU General Public License v2.0 | 4 votes |
|
/** {@inheritDoc}
* @since 3.1
*/
public DerivativeStructure value(final DerivativeStructure t)
throws DimensionMismatchException {
final double scaledX = (normalized ? FastMath.PI : 1) * t.getValue();
final double scaledX2 = scaledX * scaledX;
double[] f = new double[t.getOrder() + 1];
if (FastMath.abs(scaledX) <= SHORTCUT) {
for (int i = 0; i < f.length; ++i) {
final int k = i / 2;
if ((i & 0x1) == 0) {
// even derivation order
f[i] = (((k & 0x1) == 0) ? 1 : -1) *
(1.0 / (i + 1) - scaledX2 * (1.0 / (2 * i + 6) - scaledX2 / (24 * i + 120)));
} else {
// odd derivation order
f[i] = (((k & 0x1) == 0) ? -scaledX : scaledX) *
(1.0 / (i + 2) - scaledX2 * (1.0 / (6 * i + 24) - scaledX2 / (120 * i + 720)));
}
}
} else {
final double inv = 1 / scaledX;
final double cos = FastMath.cos(scaledX);
final double sin = FastMath.sin(scaledX);
f[0] = inv * sin;
// the nth order derivative of sinc has the form:
// dn(sinc(x)/dxn = [S_n(x) sin(x) + C_n(x) cos(x)] / x^(n+1)
// where S_n(x) is an even polynomial with degree n-1 or n (depending on parity)
// and C_n(x) is an odd polynomial with degree n-1 or n (depending on parity)
// S_0(x) = 1, S_1(x) = -1, S_2(x) = -x^2 + 2, S_3(x) = 3x^2 - 6...
// C_0(x) = 0, C_1(x) = x, C_2(x) = -2x, C_3(x) = -x^3 + 6x...
// the general recurrence relations for S_n and C_n are:
// S_n(x) = x S_(n-1)'(x) - n S_(n-1)(x) - x C_(n-1)(x)
// C_n(x) = x C_(n-1)'(x) - n C_(n-1)(x) + x S_(n-1)(x)
// as per polynomials parity, we can store both S_n and C_n in the same array
final double[] sc = new double[f.length];
sc[0] = 1;
double coeff = inv;
for (int n = 1; n < f.length; ++n) {
double s = 0;
double c = 0;
// update and evaluate polynomials S_n(x) and C_n(x)
final int kStart;
if ((n & 0x1) == 0) {
// even derivation order, S_n is degree n and C_n is degree n-1
sc[n] = 0;
kStart = n;
} else {
// odd derivation order, S_n is degree n-1 and C_n is degree n
sc[n] = sc[n - 1];
c = sc[n];
kStart = n - 1;
}
// in this loop, k is always even
for (int k = kStart; k > 1; k -= 2) {
// sine part
sc[k] = (k - n) * sc[k] - sc[k - 1];
s = s * scaledX2 + sc[k];
// cosine part
sc[k - 1] = (k - 1 - n) * sc[k - 1] + sc[k -2];
c = c * scaledX2 + sc[k - 1];
}
sc[0] *= -n;
s = s * scaledX2 + sc[0];
coeff *= inv;
f[n] = coeff * (s * sin + c * scaledX * cos);
}
}
if (normalized) {
double scale = FastMath.PI;
for (int i = 1; i < f.length; ++i) {
f[i] *= scale;
scale *= FastMath.PI;
}
}
return t.compose(f);
}
