You are viewing a plain text version of this content. The canonical link for it is here.
Posted to commits@commons.apache.org by lu...@apache.org on 2012/09/28 15:28:43 UTC

svn commit: r1391451 - in /commons/proper/math/trunk/src: changes/ main/java/org/apache/commons/math3/analysis/differentiation/ test/java/org/apache/commons/math3/analysis/differentiation/ test/java/org/apache/commons/math3/analysis/function/

Author: luc
Date: Fri Sep 28 13:28:42 2012
New Revision: 1391451

URL: http://svn.apache.org/viewvc?rev=1391451&view=rev
Log:
Fixed some issues in nth root derivatives at 0.

The current behavior is correct with respect to infinities and NaN being
appropriately generated, but in some cases counter-intuitive.

Modified:
    commons/proper/math/trunk/src/changes/changes.xml
    commons/proper/math/trunk/src/main/java/org/apache/commons/math3/analysis/differentiation/DSCompiler.java
    commons/proper/math/trunk/src/test/java/org/apache/commons/math3/analysis/differentiation/DerivativeStructureTest.java
    commons/proper/math/trunk/src/test/java/org/apache/commons/math3/analysis/function/SqrtTest.java

Modified: commons/proper/math/trunk/src/changes/changes.xml
URL: http://svn.apache.org/viewvc/commons/proper/math/trunk/src/changes/changes.xml?rev=1391451&r1=1391450&r2=1391451&view=diff
==============================================================================
--- commons/proper/math/trunk/src/changes/changes.xml (original)
+++ commons/proper/math/trunk/src/changes/changes.xml Fri Sep 28 13:28:42 2012
@@ -52,6 +52,9 @@ If the output is not quite correct, chec
   <body>
     <release version="3.1" date="TBD" description="
 ">
+      <action dev="luc" type="fix" >
+        Fixed some issues in nth root derivatives at 0.
+      </action>
       <action dev="tn" type="fix" issue="MATH-848">
         Fixed transformation to a Schur matrix for certain input matrices.
       </action>

Modified: commons/proper/math/trunk/src/main/java/org/apache/commons/math3/analysis/differentiation/DSCompiler.java
URL: http://svn.apache.org/viewvc/commons/proper/math/trunk/src/main/java/org/apache/commons/math3/analysis/differentiation/DSCompiler.java?rev=1391451&r1=1391450&r2=1391451&view=diff
==============================================================================
--- commons/proper/math/trunk/src/main/java/org/apache/commons/math3/analysis/differentiation/DSCompiler.java (original)
+++ commons/proper/math/trunk/src/main/java/org/apache/commons/math3/analysis/differentiation/DSCompiler.java Fri Sep 28 13:28:42 2012
@@ -952,15 +952,18 @@ public class DSCompiler {
         double[] function = new double[1 + order];
         double xk;
         if (n == 2) {
-            xk = FastMath.sqrt(operand[operandOffset]);
+            function[0] = FastMath.sqrt(operand[operandOffset]);
+            xk          = 0.5 / function[0];
         } else if (n == 3) {
-            xk = FastMath.cbrt(operand[operandOffset]);
+            function[0] = FastMath.cbrt(operand[operandOffset]);
+            xk          = 1.0 / (3.0 * function[0] * function[0]);
         } else {
-            xk = FastMath.pow(operand[operandOffset], 1.0 / n);
+            function[0] = FastMath.pow(operand[operandOffset], 1.0 / n);
+            xk          = 1.0 / (n * FastMath.pow(function[0], n - 1));
         }
         final double nReciprocal = 1.0 / n;
         final double xReciprocal = 1.0 / operand[operandOffset];
-        for (int i = 0; i <= order; ++i) {
+        for (int i = 1; i <= order; ++i) {
             function[i] = xk;
             xk *= xReciprocal * (nReciprocal - i);
         }

Modified: commons/proper/math/trunk/src/test/java/org/apache/commons/math3/analysis/differentiation/DerivativeStructureTest.java
URL: http://svn.apache.org/viewvc/commons/proper/math/trunk/src/test/java/org/apache/commons/math3/analysis/differentiation/DerivativeStructureTest.java?rev=1391451&r1=1391450&r2=1391451&view=diff
==============================================================================
--- commons/proper/math/trunk/src/test/java/org/apache/commons/math3/analysis/differentiation/DerivativeStructureTest.java (original)
+++ commons/proper/math/trunk/src/test/java/org/apache/commons/math3/analysis/differentiation/DerivativeStructureTest.java Fri Sep 28 13:28:42 2012
@@ -338,6 +338,69 @@ public class DerivativeStructureTest {
     }
 
     @Test
+    public void testRootNSingularity() {
+        for (int n = 2; n < 10; ++n) {
+            for (int maxOrder = 0; maxOrder < 12; ++maxOrder) {
+                DerivativeStructure dsZero = new DerivativeStructure(1, maxOrder, 0, 0.0);
+                DerivativeStructure rootN  = dsZero.rootN(n);
+                Assert.assertEquals(0.0, rootN.getValue(), 1.0e-20);
+                if (maxOrder > 0) {
+                    Assert.assertTrue(Double.isInfinite(rootN.getPartialDerivative(1)));
+                    Assert.assertTrue(rootN.getPartialDerivative(1) > 0);
+                    for (int order = 2; order <= maxOrder; ++order) {
+                        // the following checks shows a LIMITATION of the current implementation
+                        // we have no way to tell dsZero is a pure linear variable x = 0
+                        // we only say: "dsZero is a structure with value = 0.0,
+                        // first derivative = 1.0, second and higher derivatives = 0.0".
+                        // Function composition rule for second derivatives is:
+                        // d2[f(g(x))]/dx2 = f''(g(x)) * [g'(x)]^2 + f'(g(x)) * g''(x)
+                        // when function f is the nth root and x = 0 we have:
+                        // f(0) = 0, f'(0) = +infinity, f''(0) = -infinity (and higher
+                        // derivatives keep switching between +infinity and -infinity)
+                        // so given that in our case dsZero represents g, we have g(x) = 0,
+                        // g'(x) = 1 and g''(x) = 0
+                        // applying the composition rules gives:
+                        // d2[f(g(x))]/dx2 = f''(g(x)) * [g'(x)]^2 + f'(g(x)) * g''(x)
+                        //                 = -infinity * 1^2 + +infinity * 0
+                        //                 = -infinity + NaN
+                        //                 = NaN
+                        // if we knew dsZero is really the x variable and not the identity
+                        // function applied to x, we would not have computed f'(g(x)) * g''(x)
+                        // and we would have found that the result was -infinity and not NaN
+                        Assert.assertTrue(Double.isNaN(rootN.getPartialDerivative(order)));
+                    }
+                }
+
+                // the following shows that the limitation explained above is NOT a bug...
+                // if we set up the higher order derivatives for g appropriately, we do
+                // compute the higher order derivatives of the composition correctly
+                double[] gDerivatives = new double[ 1 + maxOrder];
+                gDerivatives[0] = 0.0;
+                for (int k = 1; k <= maxOrder; ++k) {
+                    gDerivatives[k] = FastMath.pow(-1.0, k + 1);
+                }
+                DerivativeStructure correctRoot = new DerivativeStructure(1, maxOrder, gDerivatives).rootN(n);
+                Assert.assertEquals(0.0, correctRoot.getValue(), 1.0e-20);
+                if (maxOrder > 0) {
+                    Assert.assertTrue(Double.isInfinite(correctRoot.getPartialDerivative(1)));
+                    Assert.assertTrue(correctRoot.getPartialDerivative(1) > 0);
+                    for (int order = 2; order <= maxOrder; ++order) {
+                        Assert.assertTrue(Double.isInfinite(correctRoot.getPartialDerivative(order)));
+                        if ((order % 2) == 0) {
+                            Assert.assertTrue(correctRoot.getPartialDerivative(order) < 0);
+                        } else {
+                            Assert.assertTrue(correctRoot.getPartialDerivative(order) > 0);
+                        }
+                    }
+                }
+
+            }
+
+        }
+
+    }
+
+    @Test
     public void testSqrtPow2() {
         double[] epsilon = new double[] { 1.0e-16, 3.0e-16, 2.0e-15, 6.0e-14, 6.0e-12 };
         for (int maxOrder = 0; maxOrder < 5; ++maxOrder) {
@@ -662,7 +725,7 @@ 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 };
+        double[] epsilon = new double[] { 5.0e-16, 3.0e-15, 2.2e-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);

Modified: commons/proper/math/trunk/src/test/java/org/apache/commons/math3/analysis/function/SqrtTest.java
URL: http://svn.apache.org/viewvc/commons/proper/math/trunk/src/test/java/org/apache/commons/math3/analysis/function/SqrtTest.java?rev=1391451&r1=1391450&r2=1391451&view=diff
==============================================================================
--- commons/proper/math/trunk/src/test/java/org/apache/commons/math3/analysis/function/SqrtTest.java (original)
+++ commons/proper/math/trunk/src/test/java/org/apache/commons/math3/analysis/function/SqrtTest.java Fri Sep 28 13:28:42 2012
@@ -65,7 +65,7 @@ public class SqrtTest {
        Assert.assertEquals(-0.1901814435781826783,  s.getPartialDerivative(2), 1.0e-16);
        Assert.assertEquals(0.23772680447272834785,  s.getPartialDerivative(3), 1.0e-16);
        Assert.assertEquals(-0.49526417598485072465,   s.getPartialDerivative(4), 1.0e-16);
-       Assert.assertEquals(1.4445205132891479465,  s.getPartialDerivative(5), 3.0e-16);
+       Assert.assertEquals(1.4445205132891479465,  s.getPartialDerivative(5), 5.0e-16);
    }
 
 }