commons-commits mailing list archives

Site index · List index
Message view « Date » · « Thread »
Top « Date » · « Thread »
From l..@apache.org
Subject svn commit: r1370951 [1/2] - in /commons/proper/math/trunk/src: changes/ main/java/org/apache/commons/math3/analysis/differentiation/ main/java/org/apache/commons/math3/util/ test/java/org/apache/commons/math3/analysis/differentiation/
Date Wed, 08 Aug 2012 20:33:43 GMT
Author: luc
Date: Wed Aug  8 20:33:43 2012
New Revision: 1370951

URL: http://svn.apache.org/viewvc?rev=1370951&view=rev
Log:
Added a new package dealing with differentials.

The package is intended to deals with one or more free parameters and
derivation order 1 or higher.

The core elements are based on Dan Kalman paper "Recursive Multivariate
Automatic Differentiation", Mathematics Magazine, vol. 75, no. 3, June
2002. For efficiency, the recursive structure is compiled as simple
loops once for each pair (number of free parameters, derivation order).

This is work in progress, there are still some features missing even in
the most basic blocks (typically the asin, acos, atan, atant2 and taylor
methods in DSCompiler). There are also still no high level
differentiator implementation.

Added:
    commons/proper/math/trunk/src/main/java/org/apache/commons/math3/analysis/differentiation/
    commons/proper/math/trunk/src/main/java/org/apache/commons/math3/analysis/differentiation/DSCompiler.java   (with props)
    commons/proper/math/trunk/src/main/java/org/apache/commons/math3/analysis/differentiation/DerivativeStructure.java   (with props)
    commons/proper/math/trunk/src/main/java/org/apache/commons/math3/analysis/differentiation/UnivariateDifferential.java   (with props)
    commons/proper/math/trunk/src/main/java/org/apache/commons/math3/analysis/differentiation/UnivariateDifferentiator.java   (with props)
    commons/proper/math/trunk/src/main/java/org/apache/commons/math3/analysis/differentiation/package-info.java   (with props)
    commons/proper/math/trunk/src/test/java/org/apache/commons/math3/analysis/differentiation/
    commons/proper/math/trunk/src/test/java/org/apache/commons/math3/analysis/differentiation/DSCompilerTest.java   (with props)
    commons/proper/math/trunk/src/test/java/org/apache/commons/math3/analysis/differentiation/DerivativeStructureTest.java   (with props)
Modified:
    commons/proper/math/trunk/src/changes/changes.xml
    commons/proper/math/trunk/src/main/java/org/apache/commons/math3/util/FastMath.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=1370951&r1=1370950&r2=1370951&view=diff
==============================================================================
--- commons/proper/math/trunk/src/changes/changes.xml (original)
+++ commons/proper/math/trunk/src/changes/changes.xml Wed Aug  8 20:33:43 2012
@@ -52,6 +52,10 @@ If the output is not quite correct, chec
   <body>
     <release version="3.1" date="TBD" description="
 ">
+      <action dev="luc" type="add" >
+        Added a new package dealing with differentials, for one or more free
+        parameters and derivation order 1 or higher.
+      </action>
       <action dev="tn" type="add" issue="MATH-777" due-to="Reid Hochstedler">
         Added additional crossover policies: "CycleCrossover", "NPointCrossover",
         "OrderedCrossover" and "UniformCrossover".

Added: 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=1370951&view=auto
==============================================================================
--- commons/proper/math/trunk/src/main/java/org/apache/commons/math3/analysis/differentiation/DSCompiler.java (added)
+++ commons/proper/math/trunk/src/main/java/org/apache/commons/math3/analysis/differentiation/DSCompiler.java Wed Aug  8 20:33:43 2012
@@ -0,0 +1,1246 @@
+/*
+ * Licensed to the Apache Software Foundation (ASF) under one or more
+ * contributor license agreements.  See the NOTICE file distributed with
+ * this work for additional information regarding copyright ownership.
+ * The ASF licenses this file to You under the Apache License, Version 2.0
+ * (the "License"); you may not use this file except in compliance with
+ * the License.  You may obtain a copy of the License at
+ *
+ *      http://www.apache.org/licenses/LICENSE-2.0
+ *
+ * Unless required by applicable law or agreed to in writing, software
+ * distributed under the License is distributed on an "AS IS" BASIS,
+ * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
+ * See the License for the specific language governing permissions and
+ * limitations under the License.
+ */
+package org.apache.commons.math3.analysis.differentiation;
+
+import java.util.ArrayList;
+import java.util.Arrays;
+import java.util.List;
+import java.util.concurrent.atomic.AtomicReference;
+
+import org.apache.commons.math3.exception.DimensionMismatchException;
+import org.apache.commons.math3.exception.NumberIsTooLargeException;
+import org.apache.commons.math3.util.ArithmeticUtils;
+import org.apache.commons.math3.util.FastMath;
+import org.apache.commons.math3.util.MathArrays;
+
+/** Class holding "compiled" computation rules for derivative structures.
+ * <p>This class implements the computation rules described in Dan Kalman's paper <a
+ * href="http://www.math.american.edu/People/kalman/pdffiles/mmgautodiff.pdf">Doubly
+ * Recursive Multivariate Automatic Differentiation</a>, Mathematics Magazine, vol. 75,
+ * no. 3, June 2002. However, in order to avoid performances bottlenecks, the recursive
+ * rules are "compiled" once in an unfold form. This class does this recursion unrolling
+ * and stores the computation rules as simple loops with pre-computed indirection arrays.</p>
+ * <p>
+ * This class maps all derivative computation into single dimension arrays that hold the
+ * value and partial derivatives. The class does not hold these arrays, which remains under
+ * the responsibility of the caller. For each combination of number of free parameters and
+ * derivation order, only one compiler is necessary, and this compiler will be used to
+ * perform computations on all arrays provided to it, which can represent hundreds or
+ * thousands of different parameters kept together with all theur partial derivatives.
+ * </p>
+ * <p>
+ * The arrays on which compilers operate contain only the partial derivatives together
+ * with the 0<sup>th</sup> derivative, i.e. the value. The partial derivatives are stored in
+ * a compiler-specific order, which can be retrieved using methods {@link
+ * #getPartialDerivativeIndex(int...) getPartialDerivativeIndex} and {@link
+ * #getPartialDerivativeOrders(int)}. The value is guaranteed to be stored as the first element
+ * (i.e. the {@link #getPartialDerivativeIndex(int...) getPartialDerivativeIndex} method returns
+ * 0 when called with 0 for all derivation orders and {@link #getPartialDerivativeOrders(int)
+ * getPartialDerivativeOrders} returns an array filled with 0 when called with 0 as the index).
+ * </p>
+ * <p>
+ * Note that the ordering changes with number of parameters and derivation order. For example
+ * given 2 parameters x and y, df/dy is stored at index 2 when derivation order is set to 1 (in
+ * this case the array has three elements: f, df/dx and df/dy). If derivation order is set to
+ * 2, then df/dy will be stored at index 3 (in this case the array has six elements: f, df/dx,
+ * df/dxdx, df/dy, df/dxdy and df/dydy).
+ * </p>
+ * <p>
+ * Given this structure, users can perform some simple operations like adding, subtracting
+ * or multiplying constants and negating the elements by themselves, knowing if they want to
+ * mutate their array or create a new array. These simple operations are not provided by
+ * the compiler. The compiler provides only the more complex operations between several arrays.
+ * </p>
+ * <p>This class is mainly used as the engine for scalar variable {@link DerivativeStructure}.
+ * It can also be used directly to hold several variables in arrays for more complex data
+ * structures. User can for example store a vector of n variables depending on three x, y
+ * and z free parameters in one array as follows:
+ * <pre>
+ *   // parameter 0 is x, parameter 1 is y, parameter 3 is z
+ *   int parameters = 3;
+ *   DSCompiler compiler = DSCompiler.getCompiler(parameters, order);
+ *   int size = compiler.getSize();
+ *
+ *   // pack all elements in a single array
+ *   double[] array = new double[n * size];
+ *   for (int i = 0; i < n; ++i) {
+ *
+ *     // we know value is guaranteed to be the first element
+ *     array[i * size] = v[i];
+ *
+ *     // we don't know where first derivatives are stored, so we ask the compiler
+ *     array[i * size + compiler.getPartialDerivativeIndex(1, 0, 0) = dvOnDx[i][0];
+ *     array[i * size + compiler.getPartialDerivativeIndex(0, 1, 0) = dvOnDy[i][0];
+ *     array[i * size + compiler.getPartialDerivativeIndex(0, 0, 1) = dvOnDz[i][0];
+ *
+ *     // we let all higher order derivatives set to 0
+ *
+ *   }
+ * </pre>
+ * Then in another function, user can perform some operations on all elements stored
+ * in the single array, such as a simple product of all variables:
+ * <pre>
+ *   // compute the product of all elements
+ *   double[] product = new double[size];
+ *   prod[0] = 1.0;
+ *   for (int i = 0; i < n; ++i) {
+ *     double[] tmp = product.clone();
+ *     compiler.multiply(tmp, 0, array, i * size, product, 0);
+ *   }
+ *
+ *   // value
+ *   double p = product[0];
+ *
+ *   // first derivatives
+ *   double dPdX = product[compiler.getPartialDerivativeIndex(1, 0, 0)];
+ *   double dPdY = product[compiler.getPartialDerivativeIndex(0, 1, 0)];
+ *   double dPdZ = product[compiler.getPartialDerivativeIndex(0, 0, 1)];
+ *
+ *   // cross derivatives (assuming order was at least 2)
+ *   double dPdXdX = product[compiler.getPartialDerivativeIndex(2, 0, 0)];
+ *   double dPdXdY = product[compiler.getPartialDerivativeIndex(1, 1, 0)];
+ *   double dPdXdZ = product[compiler.getPartialDerivativeIndex(1, 0, 1)];
+ *   double dPdYdY = product[compiler.getPartialDerivativeIndex(0, 2, 0)];
+ *   double dPdYdZ = product[compiler.getPartialDerivativeIndex(0, 1, 1)];
+ *   double dPdZdZ = product[compiler.getPartialDerivativeIndex(0, 0, 2)];
+ * </p>
+ * @see DerivativeStructure
+ * @version $Id$
+ * @since 3.1
+ */
+public class DSCompiler {
+
+    /** Array of all compilers created so far. */
+    private static AtomicReference<DSCompiler[][]> compilers =
+            new AtomicReference<DSCompiler[][]>(null);
+
+    /** Number of free parameters. */
+    private final int parameters;
+
+    /** Derivation order. */
+    private final int order;
+
+    /** Number of partial derivatives (including the single 0 order derivative element). */
+    private final int[][] sizes;
+
+    /** Indirection array for partial derivatives. */
+    private final int[][] derivativesIndirection;
+
+    /** Indirection array of the lower derivative elements. */
+    private final int[] lowerIndirection;
+
+    /** Indirection arrays for multiplication. */
+    private final int[][][] multIndirection;
+
+    /** Indirection arrays for function composition. */
+    private final int[][][] compIndirection;
+
+    /** Get the compiler for number of free parameters and order.
+     * @param parameters number of free parameters
+     * @param order derivation order
+     * @return cached rules set
+     */
+    public static DSCompiler getCompiler(int parameters, int order) {
+
+        // get the cached compilers
+        final DSCompiler[][] cache = compilers.get();
+        if (cache != null && cache.length > parameters && cache[parameters].length > order) {
+            // the compiler has already been created
+            return cache[parameters][order];
+        }
+
+        // we need to create more compilers
+        final int maxParameters = FastMath.max(parameters, cache == null ? 0 : cache.length);
+        final int maxOrder      = FastMath.max(order,     cache == null ? 0 : cache[0].length);
+        final DSCompiler[][] newCache = new DSCompiler[maxParameters + 1][maxOrder + 1];
+
+        if (cache != null) {
+            // preserve the already created compilers
+            for (int i = 0; i < cache.length; ++i) {
+                System.arraycopy(cache[i], 0, newCache[i], 0, cache[i].length);
+            }
+        }
+
+        // create the array in increasing diagonal order
+        for (int diag = 0; diag <= maxParameters + maxOrder; ++diag) {
+            for (int o = FastMath.max(0, diag - maxParameters); o <= FastMath.min(maxOrder, diag); ++o) {
+                final int p = diag - o;
+                if (newCache[p][o] == null) {
+                    final DSCompiler valueCompiler      = (p == 0) ? null : newCache[p - 1][o];
+                    final DSCompiler derivativeCompiler = (o == 0) ? null : newCache[p][o - 1];
+                    newCache[p][o] = new DSCompiler(p, o, valueCompiler, derivativeCompiler);
+                }
+            }
+        }
+
+        // atomically reset the cached compilers array
+        compilers.compareAndSet(cache, newCache);
+
+        return newCache[parameters][order];
+
+    }
+
+    /** Private constructor, reserved for the factory method {@link #getCompiler(int, int)}.
+     * @param parameters number of free parameters
+     * @param order derivation order
+     * @param valueCompiler compiler for the value part
+     * @param derivativeCompiler compiler for the derivative part
+     */
+    private DSCompiler(final int parameters, final int order,
+                       final DSCompiler valueCompiler, final DSCompiler derivativeCompiler) {
+
+        this.parameters = parameters;
+        this.order      = order;
+        this.sizes      = compileSizes(parameters, order, valueCompiler, derivativeCompiler);
+        this.derivativesIndirection =
+                compileDerivativesIndirection(parameters, order,
+                                              valueCompiler, derivativeCompiler);
+        this.lowerIndirection =
+                compileLowerIndirection(parameters, order,
+                                        valueCompiler, derivativeCompiler);
+        this.multIndirection =
+                compileMultiplicationIndirection(parameters, order,
+                                                 valueCompiler, derivativeCompiler, lowerIndirection);
+        this.compIndirection =
+                compileCompositionIndirection(parameters, order,
+                                              valueCompiler, derivativeCompiler,
+                                              sizes, derivativesIndirection, lowerIndirection);
+
+
+    }
+
+    /** Compile the sizes array.
+     * @param parameters number of free parameters
+     * @param order derivation order
+     * @param valueCompiler compiler for the value part
+     * @param derivativeCompiler compiler for the derivative part
+     * @return sizes array
+     */
+    private static int[][] compileSizes(final int parameters, final int order,
+                                        final DSCompiler valueCompiler,
+                                        final DSCompiler derivativeCompiler) {
+
+        final int[][] sizes = new int[parameters + 1][order + 1];
+        if (parameters == 0) {
+            Arrays.fill(sizes[0], 1);
+        } else {
+            System.arraycopy(valueCompiler.sizes, 0, sizes, 0, parameters);
+            sizes[parameters][0] = 1;
+            for (int i = 0; i < order; ++i) {
+                sizes[parameters][i + 1] = sizes[parameters][i] + sizes[parameters - 1][i + 1];
+            }
+        }
+
+        return sizes;
+
+    }
+
+    /** Compile the derivatives indirection array.
+     * @param parameters number of free parameters
+     * @param order derivation order
+     * @param valueCompiler compiler for the value part
+     * @param derivativeCompiler compiler for the derivative part
+     * @return derivatives indirection array
+     */
+    private static int[][] compileDerivativesIndirection(final int parameters, final int order,
+                                                      final DSCompiler valueCompiler,
+                                                      final DSCompiler derivativeCompiler) {
+
+        if (parameters == 0 || order == 0) {
+            return new int[1][parameters];
+        }
+
+        final int vSize = valueCompiler.derivativesIndirection.length;
+        final int dSize = derivativeCompiler.derivativesIndirection.length;
+        final int[][] derivativesIndirection = new int[vSize + dSize][parameters];
+
+        // set up the indices for the value part
+        for (int i = 0; i < vSize; ++i) {
+            // copy the first indices, the last one remaining set to 0
+            System.arraycopy(valueCompiler.derivativesIndirection[i], 0,
+                             derivativesIndirection[i], 0,
+                             parameters - 1);
+        }
+
+        // set up the indices for the derivative part
+        for (int i = 0; i < dSize; ++i) {
+
+            // copy the indices
+            System.arraycopy(derivativeCompiler.derivativesIndirection[i], 0,
+                             derivativesIndirection[vSize + i], 0,
+                             parameters);
+
+            // increment the derivation order for the last parameter
+            derivativesIndirection[vSize + i][parameters - 1]++;
+
+        }
+
+        return derivativesIndirection;
+
+    }
+
+    /** Compile the lower derivatives indirection array.
+     * <p>
+     * This indirection array contains the indices of all elements
+     * except derivatives for last derivation order.
+     * </p>
+     * @param parameters number of free parameters
+     * @param order derivation order
+     * @param valueCompiler compiler for the value part
+     * @param derivativeCompiler compiler for the derivative part
+     * @return lower derivatives indirection array
+     */
+    private static int[] compileLowerIndirection(final int parameters, final int order,
+                                              final DSCompiler valueCompiler,
+                                              final DSCompiler derivativeCompiler) {
+
+        if (parameters == 0 || order <= 1) {
+            return new int[] { 0 };
+        }
+
+        // this is an implementation of definition 6 in Dan Kalman's paper.
+        final int vSize = valueCompiler.lowerIndirection.length;
+        final int dSize = derivativeCompiler.lowerIndirection.length;
+        final int[] lowerIndirection = new int[vSize + dSize];
+        System.arraycopy(valueCompiler.lowerIndirection, 0, lowerIndirection, 0, vSize);
+        for (int i = 0; i < dSize; ++i) {
+            lowerIndirection[vSize + i] = valueCompiler.getSize() + derivativeCompiler.lowerIndirection[i];
+        }
+
+        return lowerIndirection;
+
+    }
+
+    /** Compile the multiplication indirection array.
+     * <p>
+     * This indirection array contains the indices of all pairs of elements
+     * involved when computing a multiplication. This allows a straightforward
+     * loop-based multiplication (see {@link #multiply(double[], int, double[], int, double[], int)}).
+     * </p>
+     * @param parameters number of free parameters
+     * @param order derivation order
+     * @param valueCompiler compiler for the value part
+     * @param derivativeCompiler compiler for the derivative part
+     * @param lowerIndirection lower derivatives indirection array
+     * @return multiplication indirection array
+     */
+    private static int[][][] compileMultiplicationIndirection(final int parameters, final int order,
+                                                           final DSCompiler valueCompiler,
+                                                           final DSCompiler derivativeCompiler,
+                                                           final int[] lowerIndirection) {
+
+        if ((parameters == 0) || (order == 0)) {
+            return new int[][][] { { { 1, 0, 0 } } };
+        }
+
+        // this is an implementation of definition 3 in Dan Kalman's paper.
+        final int vSize = valueCompiler.multIndirection.length;
+        final int dSize = derivativeCompiler.multIndirection.length;
+        final int[][][] multIndirection = new int[vSize + dSize][][];
+
+        System.arraycopy(valueCompiler.multIndirection, 0, multIndirection, 0, vSize);
+
+        for (int i = 0; i < dSize; ++i) {
+            final int[][] dRow = derivativeCompiler.multIndirection[i];
+            List<int[]> row = new ArrayList<int[]>();
+            for (int j = 0; j < dRow.length; ++j) {
+                row.add(new int[] { dRow[j][0], lowerIndirection[dRow[j][1]], vSize + dRow[j][2] });
+                row.add(new int[] { dRow[j][0], vSize + dRow[j][1], lowerIndirection[dRow[j][2]] });
+            }
+
+            // combine terms with similar derivation orders
+            final List<int[]> combined = new ArrayList<int[]>(row.size());
+            for (int j = 0; j < row.size(); ++j) {
+                final int[] termJ = row.get(j);
+                if (termJ[0] > 0) {
+                    for (int k = j + 1; k < row.size(); ++k) {
+                        final int[] termK = row.get(k);
+                        if (termJ[1] == termK[1] && termJ[2] == termK[2]) {
+                            // combine termJ and termK
+                            termJ[0] += termK[0];
+                            // make sure we will skip termK later on in the outer loop
+                            termK[0] = 0;
+                        }
+                    }
+                    combined.add(termJ);
+                }
+            }
+
+            multIndirection[vSize + i] = combined.toArray(new int[combined.size()][]);
+
+        }
+
+        return multIndirection;
+
+    }
+
+    /** Compile the function composition indirection array.
+     * <p>
+     * This indirection array contains the indices of all sets of elements
+     * involved when computing a composition. This allows a straightforward
+     * loop-based composition (see {@link #compose(double[], int, double[], double[], int)}).
+     * </p>
+     * @param parameters number of free parameters
+     * @param order derivation order
+     * @param valueCompiler compiler for the value part
+     * @param derivativeCompiler compiler for the derivative part
+     * @param sizes sizes array
+     * @param derivativesIndirection derivatives indirection array
+     * @param lowerIndirection lower derivatives indirection array
+     * @return multiplication indirection array
+     */
+    private static int[][][] compileCompositionIndirection(final int parameters, final int order,
+                                                        final DSCompiler valueCompiler,
+                                                        final DSCompiler derivativeCompiler,
+                                                        final int[][] sizes,
+                                                        final int[][] derivativesIndirection,
+                                                        final int[] lowerIndirection) {
+
+        if ((parameters == 0) || (order == 0)) {
+            return new int[][][] { { { 1, 0 } } };
+        }
+
+        final int vSize = valueCompiler.compIndirection.length;
+        final int dSize = derivativeCompiler.compIndirection.length;
+        final int[][][] compIndirection = new int[vSize + dSize][][];
+
+        // the composition rules from the value part can be reused as is
+        System.arraycopy(valueCompiler.compIndirection, 0, compIndirection, 0, vSize);
+
+        // the composition rules for the derivative part are deduced by
+        // differentiation the rules from the underlying compiler once
+        // with respect to the parameter this compiler handles and the
+        // underlying one did not handle
+        for (int i = 0; i < dSize; ++i) {
+            List<int[]> row = new ArrayList<int[]>();
+            for (int[] term : derivativeCompiler.compIndirection[i]) {
+
+                // handle term p * f_k(g(x)) * g_l1(x) * g_l2(x) * ... * g_lp(x)
+
+                // derive the first factor in the term: f_k with respect to new parameter
+                int[] derivedTermF = new int[term.length + 1];
+                derivedTermF[0] = term[0];     // p
+                derivedTermF[1] = term[1] + 1; // f_(k+1)
+                int[] orders = new int[parameters];
+                orders[parameters - 1] = 1;
+                derivedTermF[term.length] = getPartialDerivativeIndex(parameters, order, sizes, orders);  // g_1
+                for (int j = 2; j < term.length; ++j) {
+                    // convert the indices as the mapping for the current order
+                    // is different from the mapping with one less order
+                    derivedTermF[j] = convertIndex(term[j], parameters,
+                                                   derivativeCompiler.derivativesIndirection,
+                                                   parameters, order, sizes);
+                }
+                Arrays.sort(derivedTermF, 2, derivedTermF.length);
+                row.add(derivedTermF);
+
+                // derive the various g_l
+                for (int l = 2; l < term.length; ++l) {
+                    int[] derivedTermG = new int[term.length];
+                    derivedTermG[0] = term[0];
+                    derivedTermG[1] = term[1];
+                    for (int j = 2; j < term.length; ++j) {
+                        // convert the indices as the mapping for the current order
+                        // is different from the mapping with one less order
+                        derivedTermG[j] = convertIndex(term[j], parameters,
+                                                       derivativeCompiler.derivativesIndirection,
+                                                       parameters, order, sizes);
+                        if (j == l) {
+                            // derive this term
+                            System.arraycopy(derivativesIndirection[derivedTermG[j]], 0, orders, 0, parameters);
+                            orders[parameters - 1]++;
+                            derivedTermG[j] = getPartialDerivativeIndex(parameters, order, sizes, orders);
+                        }
+                    }
+                    Arrays.sort(derivedTermG, 2, derivedTermG.length);
+                    row.add(derivedTermG);
+                }
+
+            }
+
+            // combine terms with similar derivation orders
+            final List<int[]> combined = new ArrayList<int[]>(row.size());
+            for (int j = 0; j < row.size(); ++j) {
+                final int[] termJ = row.get(j);
+                if (termJ[0] > 0) {
+                    for (int k = j + 1; k < row.size(); ++k) {
+                        final int[] termK = row.get(k);
+                        boolean equals = termJ.length == termK.length;
+                        for (int l = 1; equals && l < termJ.length; ++l) {
+                            equals &= termJ[l] == termK[l];
+                        }
+                        if (equals) {
+                            // combine termJ and termK
+                            termJ[0] += termK[0];
+                            // make sure we will skip termK later on in the outer loop
+                            termK[0] = 0;
+                        }
+                    }
+                    combined.add(termJ);
+                }
+            }
+
+            compIndirection[vSize + i] = combined.toArray(new int[combined.size()][]);
+
+        }
+
+        return compIndirection;
+
+    }
+
+    /** Get the index of a partial derivative in the array.
+     * <p>
+     * If all orders are set to 0, then the 0<sup>th</sup> order derivative
+     * is returned, which is the value of the function. The index for this
+     * 0<sup>th</sup> order derivative is always 0. the indices of higher
+     * order derivatives is between 1 and {@link #getSize() - 1)}.
+     * </p>
+     * <p>
+     * This method is the inverse of method {@link #getPartialDerivativeOrders(int)}
+     * </p>
+     * @param orders derivation orders with respect to each parameter
+     * @return index of the partial derivative
+     * @exception DimensionMismatchException if the numbers of parameters does not
+     * match the instance
+     * @exception NumberIsTooLargeException if sum of derivation orders is larger
+     * than the instance limits
+     * @see #getPartialDerivativeOrders(int)
+     */
+    public int getPartialDerivativeIndex(final int ... orders)
+            throws DimensionMismatchException, NumberIsTooLargeException {
+
+        // safety check
+        if (orders.length != getFreeParameters()) {
+            throw new DimensionMismatchException(orders.length, getFreeParameters());
+        }
+
+        return getPartialDerivativeIndex(parameters, order, sizes, orders);
+
+    }
+
+    /** Get the index of a partial derivative in an array.
+     * @param parameters number of free parameters
+     * @param order derivation order
+     * @param sizes sizes array
+     * @param orders derivation orders with respect to each parameter
+     * (the lenght of this array must match the number of parameters)
+     * @return index of the partial derivative
+     * @exception NumberIsTooLargeException if sum of derivation orders is larger
+     * than the instance limits
+     */
+    private static int getPartialDerivativeIndex(final int parameters, final int order,
+                                                 final int[][] sizes, final int ... orders)
+        throws NumberIsTooLargeException {
+
+        // the value is obtained by diving into the recursive Dan Kalman's structure
+        // this is theorem 2 of his paper, with recursion replaced by iteration
+        int index     = 0;
+        int m         = order;
+        int ordersSum = 0;
+        for (int i = parameters - 1; i >= 0; --i) {
+
+            // derivative order for current free parameter
+            int derivativeOrder = orders[i];
+
+            // safety check
+            ordersSum += derivativeOrder;
+            if (ordersSum > order) {
+                throw new NumberIsTooLargeException(ordersSum, order, true);
+            }
+
+            while (derivativeOrder-- > 0) {
+                // as long as we differentiate according to current free parameter,
+                // we have to skip the value part and dive into the derivative part
+                // so we add the size of the value part to the base index
+                index += sizes[i][m--];
+            }
+
+        }
+
+        return index;
+
+    }
+
+    /** Convert an index from one (parameters, order) structure to another.
+     * @param index index of a partial derivative in source derivative structure
+     * @param srcP number of free parameters in source derivative structure
+     * @param srcDerivativesIndirection derivatives indirection array for the source
+     * derivative structure
+     * @param destP number of free parameters in destination derivative structure
+     * @param destO derivation order in destination derivative structure
+     * @param destSizes sizes array for the destination derivative structure
+     * @return index of the partial derivative with the <em>same</em> characteristics
+     * in destination derivative structure
+     */
+    private static int convertIndex(final int index,
+                                    final int srcP, final int[][] srcDerivativesIndirection,
+                                    final int destP, final int destO, final int[][] destSizes) {
+        int[] orders = new int[destP];
+        System.arraycopy(srcDerivativesIndirection[index], 0, orders, 0, FastMath.min(srcP, destP));
+        return getPartialDerivativeIndex(destP, destO, destSizes, orders);
+    }
+
+    /** Get the derivation orders for a specific index in the array.
+     * <p>
+     * This method is the inverse of {@link #getPartialDerivativeIndex(int...)}.
+     * </p>
+     * @param index of the partial derivative
+     * @return orders derivation orders with respect to each parameter
+     * @see #getPartialDerivativeIndex(int...)
+     */
+    public int[] getPartialDerivativeOrders(final int index) {
+        return derivativesIndirection[index];
+    }
+
+    /** Get the number of free parameters.
+     * @return number of free parameters
+     */
+    public int getFreeParameters() {
+        return parameters;
+    }
+
+    /** Get the derivation order.
+     * @return derivation order
+     */
+    public int getOrder() {
+        return order;
+    }
+
+    /** Get the array size required for holding partial derivatives data.
+     * <p>
+     * This number includes the single 0 order derivative element, which is
+     * guaranteed to be stored in the first element of the array.
+     * </p>
+     * @return array size required for holding partial derivatives data
+     */
+    public int getSize() {
+        return sizes[parameters][order];
+    }
+
+    /** Compute linear combination.
+     * The derivative structure built will be a1 * ds1 + a2 * ds2
+     * @param a1 first scale factor
+     * @param c1 first base (unscaled) component
+     * @param offset1 offset of first operand in its array
+     * @param a2 second scale factor
+     * @param c2 second base (unscaled) component
+     * @param offset2 offset of second operand in its array
+     * @param result array where result must be stored (it may be
+     * one of the input arrays)
+     * @param resultOffset offset of the result in its array
+     */
+    public void linearCombination(final double a1, final double[] c1, final int offset1,
+                                  final double a2, final double[] c2, final int offset2,
+                                  final double[] result, final int resultOffset) {
+        for (int i = 0; i < getSize(); ++i) {
+            result[resultOffset + i] =
+                    MathArrays.linearCombination(a1, c1[offset1 + i], a2, c2[offset2 + i]);
+        }
+    }
+
+    /** Compute linear combination.
+     * The derivative structure built will be a1 * ds1 + a2 * ds2 + a3 * ds3 + a4 * ds4
+     * @param a1 first scale factor
+     * @param c1 first base (unscaled) component
+     * @param offset1 offset of first operand in its array
+     * @param a2 second scale factor
+     * @param c2 second base (unscaled) component
+     * @param offset2 offset of second operand in its array
+     * @param a3 third scale factor
+     * @param c3 third base (unscaled) component
+     * @param offset3 offset of third operand in its array
+     * @param result array where result must be stored (it may be
+     * one of the input arrays)
+     * @param resultOffset offset of the result in its array
+     */
+    public void linearCombination(final double a1, final double[] c1, final int offset1,
+                                  final double a2, final double[] c2, final int offset2,
+                                  final double a3, final double[] c3, final int offset3,
+                                  final double[] result, final int resultOffset) {
+        for (int i = 0; i < getSize(); ++i) {
+            result[resultOffset + i] =
+                    MathArrays.linearCombination(a1, c1[offset1 + i],
+                                                 a2, c2[offset2 + i],
+                                                 a3, c3[offset3 + i]);
+        }
+    }
+
+    /** Compute linear combination.
+     * The derivative structure built will be a1 * ds1 + a2 * ds2 + a3 * ds3 + a4 * ds4
+     * @param a1 first scale factor
+     * @param c1 first base (unscaled) component
+     * @param offset1 offset of first operand in its array
+     * @param a2 second scale factor
+     * @param c2 second base (unscaled) component
+     * @param offset2 offset of second operand in its array
+     * @param a3 third scale factor
+     * @param c3 third base (unscaled) component
+     * @param offset3 offset of third operand in its array
+     * @param a4 fourth scale factor
+     * @param c4 fourth base (unscaled) component
+     * @param offset4 offset of fourth operand in its array
+     * @param result array where result must be stored (it may be
+     * one of the input arrays)
+     * @param resultOffset offset of the result in its array
+     */
+    public void linearCombination(final double a1, final double[] c1, final int offset1,
+                                  final double a2, final double[] c2, final int offset2,
+                                  final double a3, final double[] c3, final int offset3,
+                                  final double a4, final double[] c4, final int offset4,
+                                  final double[] result, final int resultOffset) {
+        for (int i = 0; i < getSize(); ++i) {
+            result[resultOffset + i] =
+                    MathArrays.linearCombination(a1, c1[offset1 + i],
+                                                 a2, c2[offset2 + i],
+                                                 a3, c3[offset3 + i],
+                                                 a4, c4[offset4 + i]);
+        }
+    }
+
+    /** Perform addition of two derivative structures.
+     * @param lhs array holding left hand side of addition
+     * @param lhsOffset offset of the left hand side in its array
+     * @param rhs array right hand side of addition
+     * @param rhsOffset offset of the right hand side in its array
+     * @param result array where result must be stored (it may be
+     * one of the input arrays)
+     * @param resultOffset offset of the result in its array
+     */
+    public void add(final double[] lhs, final int lhsOffset,
+                    final double[] rhs, final int rhsOffset,
+                    final double[] result, final int resultOffset) {
+        for (int i = 0; i < getSize(); ++i) {
+            result[resultOffset + i] = lhs[lhsOffset + i] + rhs[rhsOffset + i];
+        }
+    }
+    /** Perform subtraction of two derivative structures.
+     * @param lhs array holding left hand side of subtraction
+     * @param lhsOffset offset of the left hand side in its array
+     * @param rhs array right hand side of subtraction
+     * @param rhsOffset offset of the right hand side in its array
+     * @param result array where result must be stored (it may be
+     * one of the input arrays)
+     * @param resultOffset offset of the result in its array
+     */
+    public void subtract(final double[] lhs, final int lhsOffset,
+                         final double[] rhs, final int rhsOffset,
+                         final double[] result, final int resultOffset) {
+        for (int i = 0; i < getSize(); ++i) {
+            result[resultOffset + i] = lhs[lhsOffset + i] - rhs[rhsOffset + i];
+        }
+    }
+
+    /** Perform multiplication of two derivative structures.
+     * @param lhs array holding left hand side of multiplication
+     * @param lhsOffset offset of the left hand side in its array
+     * @param rhs array right hand side of multiplication
+     * @param rhsOffset offset of the right hand side in its array
+     * @param result array where result must be stored (for
+     * multiplication the result array <em>cannot</em> be one of
+     * the input arrays)
+     * @param resultOffset offset of the result in its array
+     */
+    public void multiply(final double[] lhs, final int lhsOffset,
+                         final double[] rhs, final int rhsOffset,
+                         final double[] result, final int resultOffset) {
+        for (int i = 0; i < multIndirection.length; ++i) {
+            final int[][] mappingI = multIndirection[i];
+            double r = 0;
+            for (int j = 0; j < mappingI.length; ++j) {
+                r += mappingI[j][0] *
+                     lhs[lhsOffset + mappingI[j][1]] *
+                     rhs[rhsOffset + mappingI[j][2]];
+            }
+            result[resultOffset + i] = r;
+        }
+    }
+
+    /** Perform division of two derivative structures.
+     * @param lhs array holding left hand side of division
+     * @param lhsOffset offset of the left hand side in its array
+     * @param rhs array right hand side of division
+     * @param rhsOffset offset of the right hand side in its array
+     * @param result array where result must be stored (for
+     * division the result array <em>cannot</em> be one of
+     * the input arrays)
+     * @param resultOffset offset of the result in its array
+     */
+    public void divide(final double[] lhs, final int lhsOffset,
+                       final double[] rhs, final int rhsOffset,
+                       final double[] result, final int resultOffset) {
+        final double[] reciprocal = new double[getSize()];
+        pow(rhs, lhsOffset, -1, reciprocal, 0);
+        multiply(lhs, lhsOffset, reciprocal, 0, result, resultOffset);
+    }
+
+    /** Perform remainder of two derivative structures.
+     * @param lhs array holding left hand side of remainder
+     * @param lhsOffset offset of the left hand side in its array
+     * @param rhs array right hand side of remainder
+     * @param rhsOffset offset of the right hand side in its array
+     * @param result array where result must be stored (it may be
+     * one of the input arrays)
+     * @param resultOffset offset of the result in its array
+     */
+    public void remainder(final double[] lhs, final int lhsOffset,
+                          final double[] rhs, final int rhsOffset,
+                          final double[] result, final int resultOffset) {
+
+        // compute k such that lhs % rhs = lhs - k rhs
+        final double rem = lhs[lhsOffset] % rhs[rhsOffset];
+        final double k   = FastMath.rint((lhs[lhsOffset] - rem) / rhs[rhsOffset]);
+
+        // set up value
+        result[resultOffset] = rem;
+
+        // set up partial derivatives
+        for (int i = 1; i < getSize(); ++i) {
+            result[resultOffset + i] = lhs[lhsOffset + i] - k * rhs[rhsOffset + i];
+        }
+
+    }
+
+    /** Compute power of a derivative structure.
+     * @param operand array holding the operand
+     * @param operandOffset offset of the operand in its array
+     * @param p power to apply
+     * @param result array where result must be stored (for
+     * power the result array <em>cannot</em> be the input
+     * array)
+     * @param resultOffset offset of the result in its array
+     */
+    public void pow(final double[] operand, final int operandOffset, final double p,
+                    final double[] result, final int resultOffset) {
+
+        // create the function value and derivatives
+        // [x^p, px^(p-1), p(p-1)x^(p-2), ... ]
+        double[] function = new double[1 + order];
+        double xk = FastMath.pow(operand[operandOffset], p - order);
+        for (int i = order; i > 0; --i) {
+            function[i] = xk;
+            xk *= operand[operandOffset];
+        }
+        function[0] = xk;
+        double coefficient = p;
+        for (int i = 1; i <= order; ++i) {
+            function[i] *= coefficient;
+            coefficient *= p - i;
+        }
+
+        // apply function composition
+        compose(operand, operandOffset, function, result, resultOffset);
+
+    }
+
+    /** Compute integer power of a derivative structure.
+     * @param operand array holding the operand
+     * @param operandOffset offset of the operand in its array
+     * @param n power to apply
+     * @param result array where result must be stored (for
+     * power the result array <em>cannot</em> be the input
+     * array)
+     * @param resultOffset offset of the result in its array
+     */
+    public void pow(final double[] operand, final int operandOffset, final int n,
+                    final double[] result, final int resultOffset) {
+
+        if (n == 0) {
+            // special case, x^0 = 1 for all x
+            result[resultOffset] = 1.0;
+            Arrays.fill(result, resultOffset + 1, resultOffset + getSize(), 0);
+            return;
+        }
+
+        // create the power function value and derivatives
+        // [x^n, nx^(n-1), n(n-1)x^(n-2), ... ]
+        double[] function = new double[1 + order];
+
+        if (n > 0) {
+            // strictly positive power
+            final int maxOrder = FastMath.min(order, n);
+            double xk = FastMath.pow(operand[operandOffset], n - maxOrder);
+            for (int i = maxOrder; i > 0; --i) {
+                function[i] = xk;
+                xk *= operand[operandOffset];
+            }
+            function[0] = xk;
+        } else {
+            // strictly negative power
+            final double inv = 1.0 / operand[operandOffset];
+            double xk = FastMath.pow(inv, -n);
+            for (int i = 0; i <= order; ++i) {
+                function[i] = xk;
+                xk *= inv;
+            }
+        }
+
+        double coefficient = n;
+        for (int i = 1; i <= order; ++i) {
+            function[i] *= coefficient;
+            coefficient *= n - i;
+        }
+
+        // apply function composition
+        compose(operand, operandOffset, function, result, resultOffset);
+
+    }
+
+    /** Compute n<sup>th</sup> root of a derivative structure.
+     * @param operand array holding the operand
+     * @param operandOffset offset of the operand in its array
+     * @param n order of the root
+     * @param result array where result must be stored (for
+     * n<sup>th</sup> root the result array <em>cannot</em> be the input
+     * array)
+     * @param resultOffset offset of the result in its array
+     */
+    public void rootN(final double[] operand, final int operandOffset, final int n,
+                      final double[] result, final int resultOffset) {
+
+        // create the function value and derivatives
+        // [x^(1/n), (1/n)x^((1/n)-1), (1-n)/n^2x^((1/n)-2), ... ]
+        double[] function = new double[1 + order];
+        double xk;
+        if (n == 2) {
+            xk = FastMath.sqrt(operand[operandOffset]);
+        } else if (n == 3) {
+            xk = FastMath.cbrt(operand[operandOffset]);
+        } else {
+            xk = FastMath.pow(operand[operandOffset], 1.0 / n);
+        }
+        final double nReciprocal = 1.0 / n;
+        final double xReciprocal = 1.0 / operand[operandOffset];
+        for (int i = 0; i <= order; ++i) {
+            function[i] = xk;
+            xk *= xReciprocal * (nReciprocal - i);
+        }
+
+        // apply function composition
+        compose(operand, operandOffset, function, result, resultOffset);
+
+    }
+
+    /** Compute exponential of a derivative structure.
+     * @param operand array holding the operand
+     * @param operandOffset offset of the operand in its array
+     * @param result array where result must be stored (for
+     * exponential the result array <em>cannot</em> be the input
+     * array)
+     * @param resultOffset offset of the result in its array
+     */
+    public void exp(final double[] operand, final int operandOffset,
+                    final double[] result, final int resultOffset) {
+
+        // create the function value and derivatives
+        double[] function = new double[1 + order];
+        Arrays.fill(function, FastMath.exp(operand[operandOffset]));
+
+        // apply function composition
+        compose(operand, operandOffset, function, result, resultOffset);
+
+    }
+
+    /** Compute natural logarithm of a derivative structure.
+     * @param operand array holding the operand
+     * @param operandOffset offset of the operand in its array
+     * @param result array where result must be stored (for
+     * logarithm the result array <em>cannot</em> be the input
+     * array)
+     * @param resultOffset offset of the result in its array
+     */
+    public void log(final double[] operand, final int operandOffset,
+                    final double[] result, final int resultOffset) {
+
+        // create the function value and derivatives
+        double[] function = new double[1 + order];
+        function[0] = FastMath.log(operand[operandOffset]);
+        if (order > 0) {
+            double inv = 1.0 / operand[operandOffset];
+            double xk  = inv;
+            for (int i = 1; i <= order; ++i) {
+                function[i] = xk;
+                xk *= -i * inv;
+            }
+        }
+
+        // apply function composition
+        compose(operand, operandOffset, function, result, resultOffset);
+
+    }
+
+    /** Compute cosine of a derivative structure.
+     * @param operand array holding the operand
+     * @param operandOffset offset of the operand in its array
+     * @param result array where result must be stored (for
+     * cosine the result array <em>cannot</em> be the input
+     * array)
+     * @param resultOffset offset of the result in its array
+     */
+    public void cos(final double[] operand, final int operandOffset,
+                    final double[] result, final int resultOffset) {
+
+        // create the function value and derivatives
+        double[] function = new double[1 + order];
+        function[0] = FastMath.cos(operand[operandOffset]);
+        if (order > 0) {
+            function[1] = -FastMath.sin(operand[operandOffset]);
+            for (int i = 2; i <= order; ++i) {
+                function[i] = -function[i - 2];
+            }
+        }
+
+        // apply function composition
+        compose(operand, operandOffset, function, result, resultOffset);
+
+    }
+
+    /** Compute sine of a derivative structure.
+     * @param operand array holding the operand
+     * @param operandOffset offset of the operand in its array
+     * @param result array where result must be stored (for
+     * sine the result array <em>cannot</em> be the input
+     * array)
+     * @param resultOffset offset of the result in its array
+     */
+    public void sin(final double[] operand, final int operandOffset,
+                    final double[] result, final int resultOffset) {
+
+        // create the function value and derivatives
+        double[] function = new double[1 + order];
+        function[0] = FastMath.sin(operand[operandOffset]);
+        if (order > 0) {
+            function[1] = FastMath.cos(operand[operandOffset]);
+            for (int i = 2; i <= order; ++i) {
+                function[i] = -function[i - 2];
+            }
+        }
+
+        // apply function composition
+        compose(operand, operandOffset, function, result, resultOffset);
+
+    }
+
+    /** Compute tangent of a derivative structure.
+     * @param operand array holding the operand
+     * @param operandOffset offset of the operand in its array
+     * @param result array where result must be stored (for
+     * tangent the result array <em>cannot</em> be the input
+     * array)
+     * @param resultOffset offset of the result in its array
+     */
+    public void tan(final double[] operand, final int operandOffset,
+                    final double[] result, final int resultOffset) {
+
+        // create the function value and derivatives
+        double[] function = new double[1 + order];
+        final double z = operand[operandOffset];
+        final double t = FastMath.tan(z);
+        function[0] = t;
+        if (order > 0) {
+            final double secant2 = 1 + t * t;
+            function[1] = secant2;
+            for (int n = 2; n <= order; ++n) {
+                final int signN4 = ((n & 0x02) == 0) ? 1 : -1;
+                double outerSum = 0;
+                int sign = 1;
+                double secant2Kp2 = secant2;
+                for (int k = 0; k < n; ++k) {
+                    double innerSum = 0;
+                    for (int j = 0; j < k; ++j) {
+                        final double alpha = 2 * (k - j) * z;
+                        final double sc  = ((n & 0x01) == 0) ? FastMath.sin(alpha) : FastMath.cos(alpha);
+                         innerSum += sc * signN4 * ArithmeticUtils.pow(k - j, n - 1) *
+                                    ArithmeticUtils.binomialCoefficient(2 * k, j);
+                    }
+                    double twoNm2K = (n >= 2 * k) ? (1 << (n - 2 * k)) : (1.0 / (1 << (2 * k - n)));
+                    outerSum  += sign * innerSum * ArithmeticUtils.binomialCoefficient(n - 1, k) *
+                                 twoNm2K * secant2Kp2 / (k + 1);
+                    sign       = -sign;
+                    secant2Kp2 *= secant2;
+                }
+                function[n] = n * outerSum;
+            }
+        }
+
+        // apply function composition
+        compose(operand, operandOffset, function, result, resultOffset);
+
+    }
+
+    /** Compute arc cosine of a derivative structure.
+     * @param operand array holding the operand
+     * @param operandOffset offset of the operand in its array
+     * @param result array where result must be stored (for
+     * arc cosine the result array <em>cannot</em> be the input
+     * array)
+     * @param resultOffset offset of the result in its array
+     */
+    public void acos(final double[] operand, final int operandOffset,
+                    final double[] result, final int resultOffset) {
+
+        // create the function value and derivatives
+        double[] function = new double[1 + order];
+        final double x = operand[operandOffset];
+        function[0] = FastMath.acos(x);
+        if (order > 0) {
+            function[1] = -1.0 / FastMath.sqrt(1 - x * x);
+            for (int i = 2; i <= order; ++i) {
+                // TODO compute higher order derivatives
+                function[i] = Double.NaN;
+            }
+        }
+
+        // apply function composition
+        compose(operand, operandOffset, function, result, resultOffset);
+
+    }
+
+    /** Compute arc sine of a derivative structure.
+     * @param operand array holding the operand
+     * @param operandOffset offset of the operand in its array
+     * @param result array where result must be stored (for
+     * arc sine the result array <em>cannot</em> be the input
+     * array)
+     * @param resultOffset offset of the result in its array
+     */
+    public void asin(final double[] operand, final int operandOffset,
+                    final double[] result, final int resultOffset) {
+
+        // create the function value and derivatives
+        double[] function = new double[1 + order];
+        final double x = operand[operandOffset];
+        function[0] = FastMath.asin(x);
+        if (order > 0) {
+            function[1] = 1.0 / FastMath.sqrt(1 - x * x);
+            for (int i = 2; i <= order; ++i) {
+                // TODO compute higher order derivatives
+                function[i] = Double.NaN;
+            }
+        }
+
+        // apply function composition
+        compose(operand, operandOffset, function, result, resultOffset);
+
+    }
+
+    /** Compute arc tangent of a derivative structure.
+     * @param operand array holding the operand
+     * @param operandOffset offset of the operand in its array
+     * @param result array where result must be stored (for
+     * arc tangent the result array <em>cannot</em> be the input
+     * array)
+     * @param resultOffset offset of the result in its array
+     */
+    public void atan(final double[] operand, final int operandOffset,
+                     final double[] result, final int resultOffset) {
+
+        // create the function value and derivatives
+        double[] function = new double[1 + order];
+        final double x = operand[operandOffset];
+        function[0] = FastMath.atan(x);
+        if (order > 0) {
+            function[1] = 1.0 / (1 + x * x);
+            for (int i = 2; i <= order; ++i) {
+                // TODO compute higher order derivatives
+                function[i] = Double.NaN;
+            }
+        }
+
+        // apply function composition
+        compose(operand, operandOffset, function, result, resultOffset);
+
+    }
+
+    /** Compute two arguments arc tangent of a derivative structure.
+     * @param y array holding the first operand
+     * @param yOffset offset of the first operand in its array
+     * @param x array holding the second operand
+     * @param xOffset offset of the second operand in its array
+     * @param result array where result must be stored (for
+     * two arguments arc tangent the result array <em>cannot</em>
+     * be the input array)
+     * @param resultOffset offset of the result in its array
+     */
+    public void atan2(final double[] y, final int yOffset,
+                      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 composition of a derivative structure by a function.
+     * @param operand array holding the operand
+     * @param operandOffset offset of the operand in its array
+     * @param f array of value and derivatives of the function at
+     * the current point (i.e. at {@code operand[operandOffset]}).
+     * @param result array where result must be stored (for
+     * composition the result array <em>cannot</em> be the input
+     * array)
+     * @param resultOffset offset of the result in its array
+     */
+    public void compose(final double[] operand, final int operandOffset, final double[] f,
+                        final double[] result, final int resultOffset) {
+        for (int i = 0; i < compIndirection.length; ++i) {
+            final int[][] mappingI = compIndirection[i];
+            double r = 0;
+            for (int j = 0; j < mappingI.length; ++j) {
+                final int[] mappingIJ = mappingI[j];
+                double product = mappingIJ[0] * f[mappingIJ[1]];
+                for (int k = 2; k < mappingIJ.length; ++k) {
+                    product *= operand[operandOffset + mappingIJ[k]];
+                }
+                r += product;
+            }
+            result[resultOffset + i] = r;
+        }
+    }
+
+    /** Evaluate Taylor expansion a derivative structure.
+     * @param ds array holding the derivative structure 
+     * @param dsOffset offset of the derivative structure in its array
+     * @param offsets parameters offsets (dx, dy, ...)
+     * @return value of the Taylor expansion at x+dx, y.dy, ...
+     */
+    public double taylor(final double[] ds, final int dsOffset, final double ... offsets) {
+        // TODO
+        return Double.NaN;
+    }
+
+    /** Check rules set compatibility.
+     * @param compiler other compiler to check against instance
+     * @exception DimensionMismatchException if number of free parameters or orders are inconsistent
+     */
+    public void checkCompatibility(final DSCompiler compiler)
+            throws DimensionMismatchException {
+        if (parameters != compiler.parameters) {
+            throw new DimensionMismatchException(parameters, compiler.parameters);
+        }
+        if (order != compiler.order) {
+            throw new DimensionMismatchException(order, compiler.order);
+        }
+    }
+
+}

Propchange: commons/proper/math/trunk/src/main/java/org/apache/commons/math3/analysis/differentiation/DSCompiler.java
------------------------------------------------------------------------------
    svn:eol-style = native

Propchange: commons/proper/math/trunk/src/main/java/org/apache/commons/math3/analysis/differentiation/DSCompiler.java
------------------------------------------------------------------------------
    svn:keywords = "Author Date Id Revision"

Added: commons/proper/math/trunk/src/main/java/org/apache/commons/math3/analysis/differentiation/DerivativeStructure.java
URL: http://svn.apache.org/viewvc/commons/proper/math/trunk/src/main/java/org/apache/commons/math3/analysis/differentiation/DerivativeStructure.java?rev=1370951&view=auto
==============================================================================
--- commons/proper/math/trunk/src/main/java/org/apache/commons/math3/analysis/differentiation/DerivativeStructure.java (added)
+++ commons/proper/math/trunk/src/main/java/org/apache/commons/math3/analysis/differentiation/DerivativeStructure.java Wed Aug  8 20:33:43 2012
@@ -0,0 +1,576 @@
+/*
+ * Licensed to the Apache Software Foundation (ASF) under one or more
+ * contributor license agreements.  See the NOTICE file distributed with
+ * this work for additional information regarding copyright ownership.
+ * The ASF licenses this file to You under the Apache License, Version 2.0
+ * (the "License"); you may not use this file except in compliance with
+ * the License.  You may obtain a copy of the License at
+ *
+ *      http://www.apache.org/licenses/LICENSE-2.0
+ *
+ * Unless required by applicable law or agreed to in writing, software
+ * distributed under the License is distributed on an "AS IS" BASIS,
+ * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
+ * See the License for the specific language governing permissions and
+ * limitations under the License.
+ */
+package org.apache.commons.math3.analysis.differentiation;
+
+import java.io.Serializable;
+
+import org.apache.commons.math3.Field;
+import org.apache.commons.math3.FieldElement;
+import org.apache.commons.math3.exception.DimensionMismatchException;
+import org.apache.commons.math3.exception.NumberIsTooLargeException;
+
+/** Class representing both the value and the differentials of a function.
+ * <p>This class is the workhorse of the differentiation package.</p>
+ * <p>This class is an implementation of the extension to Rall's
+ * numbers described in Dan Kalman's paper <a
+ * href="http://www.math.american.edu/People/kalman/pdffiles/mmgautodiff.pdf">Doubly
+ * Recursive Multivariate Automatic Differentiation</a>, Mathematics Magazine, vol. 75,
+ * no. 3, June 2002.</p>. Rall's numbers are an extension to the real numbers used
+ * throughout mathematical expressions; they hold the derivative together with the
+ * value of a function. Dan Kalman's derivative structures holds all partial derivatives
+ * up to any specified order, with respect to any number of free variables. Rall's
+ * number therefore can be seen as derivative structures for order one derivative and
+ * one free variable, and real numbers can be seen as derivative structures with zero
+ * order derivative and no free variables.</p>
+ * <p>{@link DerivativeStructure} instances can be used directly thanks to
+ * the arithmetic operators to the mathematical functions provided as static
+ * methods by this class (+, -, *, /, %, sin, cos ...).</p>
+ * <p>Implementing complex expressions by hand using these classes is
+ * however a complex and error-prone task, so the classical use is
+ * simply to develop computation code using standard primitive double
+ * values and to use {@link UnivariateDifferentiator differentiators} to create
+ * the {@link DerivativeStructure}-based instances.</p>
+ * <p>Instances of this class are guaranteed to be immutable.</p>
+ * @see DSCompiler
+ * @version $Id$
+ * @since 3.1
+ */
+public class DerivativeStructure implements FieldElement<DerivativeStructure>, Serializable {
+
+    /** Serializable UID. */
+    private static final long serialVersionUID = 20120730L;
+
+    /** Compiler for the current dimensions. */
+    private transient DSCompiler compiler;
+
+    /** Combined array holding all values. */
+    private final double[] data;
+
+    /** Build an instance with all values and derivatives set to 0.
+     * @param compiler compiler to use for computation
+     */
+    private DerivativeStructure(final DSCompiler compiler) {
+        this.compiler = compiler;
+        this.data     = new double[compiler.getSize()];
+    }
+
+    /** Build an instance with all values and derivatives set to 0.
+     * @param variables number of variables
+     * @param order derivation order
+     */
+    public DerivativeStructure(final int variables, final int order) {
+        this(DSCompiler.getCompiler(variables, order));
+    }
+
+    /** Build an instance representing a constant value.
+     * @param variables number of variables
+     * @param order derivation order
+     * @param value value of the constant
+     * @see #DerivativeStructure(int, int, int, double)
+     */
+    public DerivativeStructure(final int variables, final int order, final double value) {
+        this(variables, order);
+        this.data[0] = value;
+    }
+
+    /** Build an instance representing a variable.
+     * <p>Instances built using this constructor are considered
+     * to be the free variables with respect to which differentials
+     * are computed. As such, their differential with respect to
+     * themselves is +1.</p>
+     * @param variables number of variables
+     * @param order derivation order
+     * @param index index of the variable (from 0 to {@code variables - 1})
+     * @param value value of the variable
+     * @exception NumberIsTooLargeException if index is equal to variables or larger
+     * @see #DerivativeStructure(int, int, double)
+     */
+    public DerivativeStructure(final int variables, final int order,
+                               final int index, final double value)
+        throws NumberIsTooLargeException {
+        this(variables, order, value);
+
+        if (index >= variables) {
+            throw new NumberIsTooLargeException(index, variables, false);
+        }
+
+        if (order > 0) {
+            // the derivative of the variable with respect to itself is 1.0
+            data[DSCompiler.getCompiler(index, order).getSize()] = 1.0;
+        }
+
+    }
+
+    /** Linear combination constructor.
+     * The derivative structure built will be a1 * ds1 + a2 * ds2
+     * @param a1 first scale factor
+     * @param ds1 first base (unscaled) derivative structure
+     * @param a2 second scale factor
+     * @param ds2 second base (unscaled) derivative structure
+     * @exception DimensionMismatchException if number of free parameters or orders are inconsistent
+     */
+    public DerivativeStructure(final double a1, final DerivativeStructure ds1,
+                               final double a2, final DerivativeStructure ds2)
+        throws DimensionMismatchException {
+        this(ds1.compiler);
+        compiler.checkCompatibility(ds2.compiler);
+        compiler.linearCombination(a1, ds1.data, 0, a2, ds2.data, 0, data, 0);
+    }
+
+    /** Linear combination constructor.
+     * The derivative structure built will be a1 * ds1 + a2 * ds2 + a3 * ds3
+     * @param a1 first scale factor
+     * @param ds1 first base (unscaled) derivative structure
+     * @param a2 second scale factor
+     * @param ds2 second base (unscaled) derivative structure
+     * @param a3 third scale factor
+     * @param ds3 third base (unscaled) derivative structure
+     * @exception DimensionMismatchException if number of free parameters or orders are inconsistent
+     */
+    public DerivativeStructure(final double a1, final DerivativeStructure ds1,
+                               final double a2, final DerivativeStructure ds2,
+                               final double a3, final DerivativeStructure ds3)
+        throws DimensionMismatchException {
+        this(ds1.compiler);
+        compiler.checkCompatibility(ds2.compiler);
+        compiler.checkCompatibility(ds3.compiler);
+        compiler.linearCombination(a1, ds1.data, 0, a2, ds2.data, 0, a3, ds3.data, 0, data, 0);
+    }
+
+    /** Linear combination constructor.
+     * The derivative structure built will be a1 * ds1 + a2 * ds2 + a3 * ds3 + a4 * ds4
+     * @param a1 first scale factor
+     * @param ds1 first base (unscaled) derivative structure
+     * @param a2 second scale factor
+     * @param ds2 second base (unscaled) derivative structure
+     * @param a3 third scale factor
+     * @param ds3 third base (unscaled) derivative structure
+     * @param a4 fourth scale factor
+     * @param ds4 fourth base (unscaled) derivative structure
+     * @exception DimensionMismatchException if number of free parameters or orders are inconsistent
+     */
+    public DerivativeStructure(final double a1, final DerivativeStructure ds1,
+                               final double a2, final DerivativeStructure ds2,
+                               final double a3, final DerivativeStructure ds3,
+                               final double a4, final DerivativeStructure ds4)
+        throws DimensionMismatchException {
+        this(ds1.compiler);
+        compiler.checkCompatibility(ds2.compiler);
+        compiler.checkCompatibility(ds3.compiler);
+        compiler.checkCompatibility(ds4.compiler);
+        compiler.linearCombination(a1, ds1.data, 0, a2, ds2.data, 0,
+                                   a3, ds3.data, 0, a4, ds4.data, 0,
+                                   data, 0);
+    }
+
+    /** Copy constructor.
+     * @param instance to copy
+     */
+    private DerivativeStructure(final DerivativeStructure ds) {
+        this.compiler = ds.compiler;
+        this.data     = ds.data.clone();
+    }
+
+    /** Get the number of free parameters.
+     * @return number of free parameters
+     */
+    public int getFreeParameters() {
+        return compiler.getFreeParameters();
+    }
+
+    /** Get the derivation order.
+     * @return derivation order
+     */
+    public int getOrder() {
+        return compiler.getOrder();
+    }
+
+    /** Get the value part of the derivative structure.
+     * @return value part of the derivative structure
+     * @see #getPartialDerivative(int...)
+     */
+    public double getValue() {
+        return data[0];
+    }
+
+    /** Get a partial derivative.
+     * @param orders derivation orders with respect to each variable (if all orders are 0,
+     * the value is returned)
+     * @return partial derivative
+     * @see #getValue()
+     * @exception DimensionMismatchException if the numbers of variables does not
+     * match the instance
+     * @exception NumberIsTooLargeException if sum of derivation orders is larger
+     * than the instance limits
+     */
+    public double getPartialDerivative(final int ... orders)
+        throws DimensionMismatchException, NumberIsTooLargeException {
+        return data[compiler.getPartialDerivativeIndex(orders)];
+    }
+
+    /** '+' operator.
+     * @param a right hand side parameter of the operator
+     * @return this+a
+     */
+    public DerivativeStructure add(final double a) {
+        final DerivativeStructure ds = new DerivativeStructure(this);
+        ds.data[0] += a;
+        return ds;
+    }
+
+    /** '+' operator.
+     * @param a right hand side parameter of the operator
+     * @return this+a
+     * @exception DimensionMismatchException if number of free parameters or orders are inconsistent
+     */
+    public DerivativeStructure add(final DerivativeStructure a)
+        throws DimensionMismatchException {
+        compiler.checkCompatibility(a.compiler);
+        final DerivativeStructure ds = new DerivativeStructure(this);
+        compiler.add(data, 0, a.data, 0, ds.data, 0);
+        return ds;
+    }
+
+    /** '-' operator.
+     * @param a right hand side parameter of the operator
+     * @return this-a
+     */
+    public DerivativeStructure subtract(final double a) {
+        return add(-a);
+    }
+
+    /** '-' operator.
+     * @param a right hand side parameter of the operator
+     * @return this-a
+     * @exception DimensionMismatchException if number of free parameters or orders are inconsistent
+     */
+    public DerivativeStructure subtract(final DerivativeStructure a)
+        throws DimensionMismatchException {
+        compiler.checkCompatibility(a.compiler);
+        final DerivativeStructure ds = new DerivativeStructure(this);
+        compiler.subtract(data, 0, a.data, 0, ds.data, 0);
+        return ds;
+    }
+
+    /** {@inheritDoc} */
+    public DerivativeStructure multiply(final int n) {
+        return multiply((double) n);
+    }
+
+    /** '&times;' operator.
+     * @param a right hand side parameter of the operator
+     * @return this&times;a
+     */
+    public DerivativeStructure multiply(final double a) {
+        final DerivativeStructure ds = new DerivativeStructure(this);
+        for (int i = 0; i < ds.data.length; ++i) {
+            ds.data[i] *= a;
+        }
+        return ds;
+    }
+
+    /** '&times;' operator.
+     * @param a right hand side parameter of the operator
+     * @return this&times;a
+     * @exception DimensionMismatchException if number of free parameters or orders are inconsistent
+     */
+    public DerivativeStructure multiply(final DerivativeStructure a)
+        throws DimensionMismatchException {
+        compiler.checkCompatibility(a.compiler);
+        final DerivativeStructure result = new DerivativeStructure(compiler);
+        compiler.multiply(data, 0, a.data, 0, result.data, 0);
+        return result;
+    }
+
+    /** '&divides;' operator.
+     * @param a right hand side parameter of the operator
+     * @return this&divides;a
+     */
+    public DerivativeStructure divide(final double a) {
+        final DerivativeStructure ds = new DerivativeStructure(this);
+        for (int i = 0; i < ds.data.length; ++i) {
+            ds.data[i] /= a;
+        }
+        return ds;
+    }
+
+    /** '&divides;' operator.
+     * @param a right hand side parameter of the operator
+     * @return this&divides;a
+     * @exception DimensionMismatchException if number of free parameters or orders are inconsistent
+     */
+    public DerivativeStructure divide(final DerivativeStructure a)
+        throws DimensionMismatchException {
+        compiler.checkCompatibility(a.compiler);
+        final DerivativeStructure result = new DerivativeStructure(compiler);
+        compiler.divide(data, 0, a.data, 0, result.data, 0);
+        return result;
+    }
+
+    /** '%' operator.
+     * @param a right hand side parameter of the operator
+     * @return this%a
+     */
+    public DerivativeStructure remainder(final double a) {
+        final DerivativeStructure ds = new DerivativeStructure(this);
+        ds.data[0] = ds.data[0] % a;
+        return ds;
+    }
+
+    /** '%' operator.
+     * @param a right hand side parameter of the operator
+     * @return this%a
+     * @exception DimensionMismatchException if number of free parameters or orders are inconsistent
+     */
+    public DerivativeStructure remainder(final DerivativeStructure a)
+        throws DimensionMismatchException {
+        compiler.checkCompatibility(a.compiler);
+        final DerivativeStructure result = new DerivativeStructure(compiler);
+        compiler.remainder(data, 0, a.data, 0, result.data, 0);
+        return result;
+    }
+
+    /** unary '-' operator.
+     * @return -this
+     */
+    public DerivativeStructure negate() {
+        final DerivativeStructure ds = new DerivativeStructure(compiler);
+        for (int i = 0; i < ds.data.length; ++i) {
+            ds.data[i] = -data[i];
+        }
+        return ds;
+    }
+
+    /** {@inheritDoc} */
+    public DerivativeStructure reciprocal() {
+        final DerivativeStructure result = new DerivativeStructure(compiler);
+        compiler.pow(data, 0, -1, result.data, 0);
+        return result;
+    }
+
+    /** Square root.
+     * @return square root of the instance
+     */
+    public DerivativeStructure sqrt() {
+        return rootN(2);
+    }
+
+    /** Cubic root.
+     * @return cubic root of the instance
+     */
+    public DerivativeStructure cbrt() {
+        return rootN(3);
+    }
+
+    /** N<sup>th</sup> root.
+     * @param n order of the root
+     * @return n<sup>th</sup> root of the instance
+     */
+    public DerivativeStructure rootN(final int n) {
+        final DerivativeStructure result = new DerivativeStructure(compiler);
+        compiler.rootN(data, 0, n, result.data, 0);
+        return result;
+    }
+
+    /** {@inheritDoc} */
+    public Field<DerivativeStructure> getField() {
+        return new Field<DerivativeStructure>() {
+
+            /** {@inheritDoc} */
+            public DerivativeStructure getZero() {
+                return new DerivativeStructure(compiler.getFreeParameters(), compiler.getOrder(), 0.0);
+            }
+
+            /** {@inheritDoc} */
+            public DerivativeStructure getOne() {
+                return new DerivativeStructure(compiler.getFreeParameters(), compiler.getOrder(), 1.0);
+            }
+
+            /** {@inheritDoc} */
+            public Class<? extends FieldElement<DerivativeStructure>> getRuntimeClass() {
+                return DerivativeStructure.class;
+            }
+
+        };
+    }
+
+    /** Power operation.
+     * @param p power to apply
+     * @return this<sup>p</sup>
+     */
+    public DerivativeStructure pow(final double p) {
+        final DerivativeStructure result = new DerivativeStructure(compiler);
+        compiler.pow(data, 0, p, result.data, 0);
+        return result;
+    }
+
+    /** Integer power operation.
+     * @param n power to apply
+     * @return this<sup>n</sup>
+     */
+    public DerivativeStructure pow(final int n) {
+        final DerivativeStructure result = new DerivativeStructure(compiler);
+        compiler.pow(data, 0, n, result.data, 0);
+        return result;
+    }
+
+    /** Exponential.
+     * @return exponential of the instance
+     */
+    public DerivativeStructure exp() {
+        final DerivativeStructure result = new DerivativeStructure(compiler);
+        compiler.exp(data, 0, result.data, 0);
+        return result;
+    }
+
+    /** Natural logarithm.
+     * @return logarithm of the instance
+     */
+    public DerivativeStructure log() {
+        final DerivativeStructure result = new DerivativeStructure(compiler);
+        compiler.log(data, 0, result.data, 0);
+        return result;
+    }
+
+    /** Cosine operation.
+     * @return cos(this)
+     */
+    public DerivativeStructure cos() {
+        final DerivativeStructure result = new DerivativeStructure(compiler);
+        compiler.cos(data, 0, result.data, 0);
+        return result;
+    }
+
+    /** Sine operation.
+     * @return sin(this)
+     */
+    public DerivativeStructure sin() {
+        final DerivativeStructure result = new DerivativeStructure(compiler);
+        compiler.sin(data, 0, result.data, 0);
+        return result;
+    }
+
+    /** Tangent operation.
+     * @return tan(this)
+     */
+    public DerivativeStructure tan() {
+        final DerivativeStructure result = new DerivativeStructure(compiler);
+        compiler.tan(data, 0, result.data, 0);
+        return result;
+    }
+
+    /** Arc cosine operation.
+     * @return acos(this)
+     */
+    public DerivativeStructure acos() {
+        final DerivativeStructure result = new DerivativeStructure(compiler);
+        compiler.acos(data, 0, result.data, 0);
+        return result;
+    }
+
+    /** Arc sine operation.
+     * @return asin(this)
+     */
+    public DerivativeStructure asin() {
+        final DerivativeStructure result = new DerivativeStructure(compiler);
+        compiler.asin(data, 0, result.data, 0);
+        return result;
+    }
+
+    /** Arc tangent operation.
+     * @return tan(this)
+     */
+    public DerivativeStructure atan() {
+        final DerivativeStructure result = new DerivativeStructure(compiler);
+        compiler.atan(data, 0, result.data, 0);
+        return result;
+    }
+
+    /** Two arguments arc tangent operation.
+     * @param y first argument of the arc tangent
+     * @param x second argument of the arc tangent
+     * @return atan2(y, x)
+     * @exception DimensionMismatchException if number of free parameters or orders are inconsistent
+     */
+    public static DerivativeStructure atan2(final DerivativeStructure y, final DerivativeStructure x)
+        throws DimensionMismatchException {
+        y.compiler.checkCompatibility(x.compiler);
+        final DerivativeStructure result = new DerivativeStructure(y.compiler);
+        y.compiler.atan2(y.data, 0, x.data, 0, result.data, 0);
+        return result;
+    }
+
+    /** Evaluate Taylor expansion a derivative structure.
+     * @param offsets parameters offsets (dx, dy, ...)
+     * @return value of the Taylor expansion at x+dx, y.dy, ...
+     */
+    public double taylor(final double ... offsets) {
+        return compiler.taylor(data, 0, offsets);
+    }
+
+    /**
+     * Replace the instance with a data transfer object for serialization.
+     * @return data transfer object that will be serialized
+     */
+    private Object writeReplace() {
+        return new DataTransferObject(compiler.getFreeParameters(), compiler.getOrder(), data);
+    }
+
+    /** Internal class used only for serialization. */
+    private static class DataTransferObject implements Serializable {
+
+        /** Serializable UID. */
+        private static final long serialVersionUID = 20120730L;
+
+        /** Number of variables.
+         * @Serial
+         */
+        private final int variables;
+
+        /** Derivation order.
+         * @Serial
+         */
+        private final int order;
+
+        /** Partial derivatives.
+         * @Serial
+         */
+        private final double[] data;
+
+        /** Simple constructor.
+         * @param variables number of variables
+         * @param order derivation order
+         * @param data partial derivatives
+         */
+        public DataTransferObject(final int variables, final int order, final double[] data) {
+            this.variables = variables;
+            this.order     = order;
+            this.data      = data;
+        }
+
+        /** Replace the deserialized data transfer object with a {@link DerivativeStructure}.
+         * @return replacement {@link DerivativeStructure}
+         */
+        private Object readResolve() {
+            final DerivativeStructure ds = new DerivativeStructure(variables, order);
+            System.arraycopy(data, 0, ds.data, 0, ds.data.length);
+            return ds;
+        }
+
+    }
+
+}

Propchange: commons/proper/math/trunk/src/main/java/org/apache/commons/math3/analysis/differentiation/DerivativeStructure.java
------------------------------------------------------------------------------
    svn:eol-style = native

Propchange: commons/proper/math/trunk/src/main/java/org/apache/commons/math3/analysis/differentiation/DerivativeStructure.java
------------------------------------------------------------------------------
    svn:keywords = "Author Date Id Revision"

Added: commons/proper/math/trunk/src/main/java/org/apache/commons/math3/analysis/differentiation/UnivariateDifferential.java
URL: http://svn.apache.org/viewvc/commons/proper/math/trunk/src/main/java/org/apache/commons/math3/analysis/differentiation/UnivariateDifferential.java?rev=1370951&view=auto
==============================================================================
--- commons/proper/math/trunk/src/main/java/org/apache/commons/math3/analysis/differentiation/UnivariateDifferential.java (added)
+++ commons/proper/math/trunk/src/main/java/org/apache/commons/math3/analysis/differentiation/UnivariateDifferential.java Wed Aug  8 20:33:43 2012
@@ -0,0 +1,61 @@
+/*
+ * Licensed to the Apache Software Foundation (ASF) under one or more
+ * contributor license agreements.  See the NOTICE file distributed with
+ * this work for additional information regarding copyright ownership.
+ * The ASF licenses this file to You under the Apache License, Version 2.0
+ * (the "License"); you may not use this file except in compliance with
+ * the License.  You may obtain a copy of the License at
+ *
+ *      http://www.apache.org/licenses/LICENSE-2.0
+ *
+ * Unless required by applicable law or agreed to in writing, software
+ * distributed under the License is distributed on an "AS IS" BASIS,
+ * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
+ * See the License for the specific language governing permissions and
+ * limitations under the License.
+ */
+package org.apache.commons.math3.analysis.differentiation;
+
+import org.apache.commons.math3.analysis.UnivariateFunction;
+
+/** Interface for univariate functions derivatives.
+ * <p>This interface represents a simple function which computes
+ * both the value and the first derivative of a mathematical function.
+ * The derivative is computed with respect to the input variable.</p>
+ * @see UnivariateDifferentiable
+ * @see UnivariateDifferentiator
+ * @since 3.1
+ * @version $Id$
+ */
+public interface UnivariateDifferential {
+
+    /** Get the primitive function associated with this differential.
+     * <p>Each {@link UnivariateDifferential} instance is tightly bound
+     * to an {@link UnivariateDifferentiable} instance. If the state of
+     * the primitive instance changes in any way that affects the
+     * differential computation, this binding allows this change to
+     * be immediately seen by the derivative instance, there is no need
+     * to differentiate the primitive again. The existing instance is aware
+     * of the primitive changes.</p>
+     * <p>In other words in the following code snippet, the three values
+     * f1, f2 and f3 should be equal (at least at machine tolerance level)</p>
+     * <pre>
+     *    UnivariateDifferential derivative = differentiator.differentiate(derivable);
+     *    derivable.someFunctionThatMutatesHeavilyTheInstance();
+     *    double f1 = derivable.f(t);
+     *    double f2 = derivative.getPrimitive().f(t);
+     *    double f3 = derivative.f(new DerivativeStructure(variables, order, index, t)).getValue();
+     * </pre>
+     * @return primitive function bound to this derivative
+     */
+    UnivariateFunction getPrimitive();
+
+    /** Simple mathematical function.
+     * <p>{@link UnivariateDifferential} classes compute both the
+     * value and the first derivative of the function.</p>
+     * @param t function input value
+     * @return function result
+     */
+    DerivativeStructure f(DerivativeStructure t);
+
+}

Propchange: commons/proper/math/trunk/src/main/java/org/apache/commons/math3/analysis/differentiation/UnivariateDifferential.java
------------------------------------------------------------------------------
    svn:eol-style = native

Propchange: commons/proper/math/trunk/src/main/java/org/apache/commons/math3/analysis/differentiation/UnivariateDifferential.java
------------------------------------------------------------------------------
    svn:keywords = "Author Date Id Revision"

Added: commons/proper/math/trunk/src/main/java/org/apache/commons/math3/analysis/differentiation/UnivariateDifferentiator.java
URL: http://svn.apache.org/viewvc/commons/proper/math/trunk/src/main/java/org/apache/commons/math3/analysis/differentiation/UnivariateDifferentiator.java?rev=1370951&view=auto
==============================================================================
--- commons/proper/math/trunk/src/main/java/org/apache/commons/math3/analysis/differentiation/UnivariateDifferentiator.java (added)
+++ commons/proper/math/trunk/src/main/java/org/apache/commons/math3/analysis/differentiation/UnivariateDifferentiator.java Wed Aug  8 20:33:43 2012
@@ -0,0 +1,34 @@
+/*
+ * Licensed to the Apache Software Foundation (ASF) under one or more
+ * contributor license agreements.  See the NOTICE file distributed with
+ * this work for additional information regarding copyright ownership.
+ * The ASF licenses this file to You under the Apache License, Version 2.0
+ * (the "License"); you may not use this file except in compliance with
+ * the License.  You may obtain a copy of the License at
+ *
+ *      http://www.apache.org/licenses/LICENSE-2.0
+ *
+ * Unless required by applicable law or agreed to in writing, software
+ * distributed under the License is distributed on an "AS IS" BASIS,
+ * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
+ * See the License for the specific language governing permissions and
+ * limitations under the License.
+ */
+package org.apache.commons.math3.analysis.differentiation;
+
+import org.apache.commons.math3.analysis.UnivariateFunction;
+
+/** Interface defining the function differentiation operation.
+ * @version $Id$
+ * @since 3.1
+ */
+public interface UnivariateDifferentiator {
+
+    /** Create an implementation of a differential for a
+     * {@link UnivariateDifferentiable differentiable function}.
+     * @param function function to differentiate
+     * @return differential function
+     */
+    UnivariateDifferential differentiate(UnivariateFunction function);
+
+}

Propchange: commons/proper/math/trunk/src/main/java/org/apache/commons/math3/analysis/differentiation/UnivariateDifferentiator.java
------------------------------------------------------------------------------
    svn:eol-style = native

Propchange: commons/proper/math/trunk/src/main/java/org/apache/commons/math3/analysis/differentiation/UnivariateDifferentiator.java
------------------------------------------------------------------------------
    svn:keywords = "Author Date Id Revision"



Mime
View raw message