Added support for atan2 in DSCompiler.

git-svn-id: https://svn.apache.org/repos/asf/commons/proper/math/trunk@1372414 13f79535-47bb-0310-9956-ffa450edef68
This commit is contained in:
Luc Maisonobe 2012-08-13 13:06:26 +00:00
parent 1ec0b8d527
commit 07aa93aae7
2 changed files with 55 additions and 7 deletions

View File

@ -1286,14 +1286,36 @@ public class DSCompiler {
final double[] x, final int xOffset,
final double[] result, final int resultOffset) {
final double y0 = y[yOffset];
final double x0 = x[xOffset];
result[resultOffset] = FastMath.atan2(y0, x0);
if (order > 0) {
for (int i = 1; i <= order; ++i) {
// TODO compute higher order derivatives
result[resultOffset + i] = Double.NaN;
// compute r = sqrt(x^2+y^2)
double[] tmp1 = new double[getSize()];
multiply(x, xOffset, x, xOffset, tmp1, 0); // x^2
double[] tmp2 = new double[getSize()];
multiply(y, yOffset, y, yOffset, tmp2, 0); // y^2
add(tmp1, 0, tmp2, 0, tmp2, 0); // x^2 + y^2
rootN(tmp2, 0, 2, tmp1, 0); // r = sqrt(x^2 + y^2)
if (x[xOffset] >= 0) {
// compute atan2(y, x) = 2 atan(y / (r + x))
add(tmp1, 0, x, xOffset, tmp2, 0); // r + x
divide(y, yOffset, tmp2, 0, tmp1, 0); // y /(r + x)
atan(tmp1, 0, tmp2, 0); // atan(y / (r + x))
for (int i = 0; i < tmp2.length; ++i) {
result[resultOffset + i] = 2 * tmp2[i]; // 2 * atan(y / (r + x))
}
} else {
// compute atan2(y, x) = +/- pi - 2 atan(y / (r - x))
subtract(tmp1, 0, x, xOffset, tmp2, 0); // r - x
divide(y, yOffset, tmp2, 0, tmp1, 0); // y /(r - x)
atan(tmp1, 0, tmp2, 0); // atan(y / (r - x))
result[resultOffset] =
((tmp2[0] <= 0) ? -FastMath.PI : FastMath.PI) - 2 * tmp2[0]; // +/-pi - 2 * atan(y / (r - x))
for (int i = 1; i < tmp2.length; ++i) {
result[resultOffset + i] = -2 * tmp2[i]; // +/-pi - 2 * atan(y / (r - x))
}
}
}

View File

@ -542,6 +542,32 @@ public class DerivativeStructureTest {
}
}
@Test
public void testAtan2() {
double[] epsilon = new double[] { 5.0e-16, 3.0e-15, 2.0e-14, 1.0e-12, 8.0e-11 };
for (int maxOrder = 0; maxOrder < 5; ++maxOrder) {
for (double x = -1.7; x < 2; x += 0.2) {
DerivativeStructure dsX = new DerivativeStructure(2, maxOrder, 0, x);
for (double y = -1.7; y < 2; y += 0.2) {
DerivativeStructure dsY = new DerivativeStructure(2, maxOrder, 1, y);
DerivativeStructure atan2 = DerivativeStructure.atan2(dsY, dsX);
DerivativeStructure ref = dsY.divide(dsX).atan();
if (x < 0) {
ref = (y < 0) ? ref.subtract(FastMath.PI) : ref.add(FastMath.PI);
}
DerivativeStructure zero = atan2.subtract(ref);
for (int n = 0; n <= maxOrder; ++n) {
for (int m = 0; m <= maxOrder; ++m) {
if (n + m <= maxOrder) {
Assert.assertEquals(0, zero.getPartialDerivative(n, m), epsilon[n + m]);
}
}
}
}
}
}
}
@Test
public void testCompositionOneVariableY() {
double epsilon = 1.0e-13;