DerivativeStructure.java

  1. /*
  2.  * Licensed to the Apache Software Foundation (ASF) under one or more
  3.  * contributor license agreements.  See the NOTICE file distributed with
  4.  * this work for additional information regarding copyright ownership.
  5.  * The ASF licenses this file to You under the Apache License, Version 2.0
  6.  * (the "License"); you may not use this file except in compliance with
  7.  * the License.  You may obtain a copy of the License at
  8.  *
  9.  *      http://www.apache.org/licenses/LICENSE-2.0
  10.  *
  11.  * Unless required by applicable law or agreed to in writing, software
  12.  * distributed under the License is distributed on an "AS IS" BASIS,
  13.  * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
  14.  * See the License for the specific language governing permissions and
  15.  * limitations under the License.
  16.  */
  17. package org.apache.commons.math4.legacy.analysis.differentiation;

  19. import java.util.Arrays;

  21. import org.apache.commons.numbers.core.Sum;
  22. import org.apache.commons.math4.legacy.core.Field;
  23. import org.apache.commons.math4.legacy.core.FieldElement;
  24. import org.apache.commons.math4.legacy.core.RealFieldElement;
  25. import org.apache.commons.math4.legacy.exception.DimensionMismatchException;
  26. import org.apache.commons.math4.legacy.exception.NumberIsTooLargeException;
  27. import org.apache.commons.math4.core.jdkmath.JdkMath;
  28. import org.apache.commons.math4.legacy.core.MathArrays;

  30. /** Class representing both the value and the differentials of a function.
  31.  * <p>This class is the workhorse of the differentiation package.</p>
  32.  * <p>This class is an implementation of the extension to Rall's
  33.  * numbers described in Dan Kalman's paper <a
  34.  * href="http://www1.american.edu/cas/mathstat/People/kalman/pdffiles/mmgautodiff.pdf">Doubly
  35.  * Recursive Multivariate Automatic Differentiation</a>, Mathematics Magazine, vol. 75,
  36.  * no. 3, June 2002. Rall's numbers are an extension to the real numbers used
  37.  * throughout mathematical expressions; they hold the derivative together with the
  38.  * value of a function. Dan Kalman's derivative structures hold all partial derivatives
  39.  * up to any specified order, with respect to any number of free parameters. Rall's
  40.  * numbers therefore can be seen as derivative structures for order one derivative and
  41.  * one free parameter, and real numbers can be seen as derivative structures with zero
  42.  * order derivative and no free parameters.</p>
  43.  * <p>{@link DerivativeStructure} instances can be used directly thanks to
  44.  * the arithmetic operators to the mathematical functions provided as
  45.  * methods by this class (+, -, *, /, %, sin, cos ...).</p>
  46.  * <p>Implementing complex expressions by hand using these classes is
  47.  * a tedious and error-prone task but has the advantage of having no limitation
  48.  * on the derivation order despite not requiring users to compute the derivatives by
  49.  * themselves. Implementing complex expression can also be done by developing computation
  50.  * code using standard primitive double values and to use {@link
  51.  * UnivariateFunctionDifferentiator differentiators} to create the {@link
  52.  * DerivativeStructure}-based instances. This method is simpler but may be limited in
  53.  * the accuracy and derivation orders and may be computationally intensive (this is
  54.  * typically the case for {@link FiniteDifferencesDifferentiator finite differences
  55.  * differentiator}.</p>
  56.  * <p>Instances of this class are guaranteed to be immutable.</p>
  57.  * @see DSCompiler
  58.  * @since 3.1
  59.  */
  60. public class DerivativeStructure implements RealFieldElement<DerivativeStructure> {
  61.     /** Compiler for the current dimensions. */
  62.     private DSCompiler compiler;

  64.     /** Combined array holding all values. */
  65.     private final double[] data;

  67.     /** Build an instance with all values and derivatives set to 0.
  68.      * @param compiler compiler to use for computation
  69.      */
  70.     private DerivativeStructure(final DSCompiler compiler) {
  71.         this.compiler = compiler;
  72.         this.data     = new double[compiler.getSize()];
  73.     }

  75.     /** Build an instance with all values and derivatives set to 0.
  76.      * @param parameters number of free parameters
  77.      * @param order derivation order
  78.      * @throws NumberIsTooLargeException if order is too large.
  79.      */
  80.     public DerivativeStructure(final int parameters, final int order) {
  81.         this(DSCompiler.getCompiler(parameters, order));
  82.     }

  84.     /** Build an instance representing a constant value.
  85.      * @param parameters number of free parameters
  86.      * @param order derivation order
  87.      * @param value value of the constant
  88.      * @throws NumberIsTooLargeException if order is too large.
  89.      * @see #DerivativeStructure(int, int, int, double)
  90.      */
  91.     public DerivativeStructure(final int parameters, final int order, final double value) {
  92.         this(parameters, order);
  93.         this.data[0] = value;
  94.     }

  96.     /** Build an instance representing a variable.
  97.      * <p>Instances built using this constructor are considered
  98.      * to be the free variables with respect to which differentials
  99.      * are computed. As such, their differential with respect to
  100.      * themselves is +1.</p>
  101.      * @param parameters number of free parameters
  102.      * @param order derivation order
  103.      * @param index index of the variable (from 0 to {@code parameters - 1})
  104.      * @param value value of the variable
  105.      * @throws NumberIsTooLargeException if {@code index &ge; parameters}.
  106.      * @see #DerivativeStructure(int, int, double)
  107.      */
  108.     public DerivativeStructure(final int parameters, final int order,
  109.                                final int index, final double value) {
  110.         this(parameters, order, value);

  112.         if (index >= parameters) {
  113.             throw new NumberIsTooLargeException(index, parameters, false);
  114.         }

  116.         if (order > 0) {
  117.             // the derivative of the variable with respect to itself is 1.
  118.             data[DSCompiler.getCompiler(index, order).getSize()] = 1.0;
  119.         }
  120.     }

  122.     /** Linear combination constructor.
  123.      * The derivative structure built will be a1 * ds1 + a2 * ds2
  124.      * @param a1 first scale factor
  125.      * @param ds1 first base (unscaled) derivative structure
  126.      * @param a2 second scale factor
  127.      * @param ds2 second base (unscaled) derivative structure
  128.      * @throws org.apache.commons.math4.legacy.exception.DimensionMismatchException
  129.      * if number of free parameters or orders are inconsistent
  130.      */
  131.     public DerivativeStructure(final double a1, final DerivativeStructure ds1,
  132.                                final double a2, final DerivativeStructure ds2) {
  133.         this(ds1.compiler);
  134.         compiler.checkCompatibility(ds2.compiler);
  135.         compiler.linearCombination(a1, ds1.data, 0, a2, ds2.data, 0, data, 0);
  136.     }

  138.     /** Linear combination constructor.
  139.      * The derivative structure built will be a1 * ds1 + a2 * ds2 + a3 * ds3
  140.      * @param a1 first scale factor
  141.      * @param ds1 first base (unscaled) derivative structure
  142.      * @param a2 second scale factor
  143.      * @param ds2 second base (unscaled) derivative structure
  144.      * @param a3 third scale factor
  145.      * @param ds3 third base (unscaled) derivative structure
  146.      * @throws org.apache.commons.math4.legacy.exception.DimensionMismatchException
  147.      * if number of free parameters or orders are inconsistent.
  148.      */
  149.     public DerivativeStructure(final double a1, final DerivativeStructure ds1,
  150.                                final double a2, final DerivativeStructure ds2,
  151.                                final double a3, final DerivativeStructure ds3) {
  152.         this(ds1.compiler);
  153.         compiler.checkCompatibility(ds2.compiler);
  154.         compiler.checkCompatibility(ds3.compiler);
  155.         compiler.linearCombination(a1, ds1.data, 0, a2, ds2.data, 0, a3, ds3.data, 0, data, 0);
  156.     }

  158.     /** Linear combination constructor.
  159.      * The derivative structure built will be a1 * ds1 + a2 * ds2 + a3 * ds3 + a4 * ds4
  160.      * @param a1 first scale factor
  161.      * @param ds1 first base (unscaled) derivative structure
  162.      * @param a2 second scale factor
  163.      * @param ds2 second base (unscaled) derivative structure
  164.      * @param a3 third scale factor
  165.      * @param ds3 third base (unscaled) derivative structure
  166.      * @param a4 fourth scale factor
  167.      * @param ds4 fourth base (unscaled) derivative structure
  168.      * @throws DimensionMismatchException if number of free parameters or orders are inconsistent.
  169.      */
  170.     public DerivativeStructure(final double a1, final DerivativeStructure ds1,
  171.                                final double a2, final DerivativeStructure ds2,
  172.                                final double a3, final DerivativeStructure ds3,
  173.                                final double a4, final DerivativeStructure ds4) {
  174.         this(ds1.compiler);
  175.         compiler.checkCompatibility(ds2.compiler);
  176.         compiler.checkCompatibility(ds3.compiler);
  177.         compiler.checkCompatibility(ds4.compiler);
  178.         compiler.linearCombination(a1, ds1.data, 0, a2, ds2.data, 0,
  179.                                    a3, ds3.data, 0, a4, ds4.data, 0,
  180.                                    data, 0);
  181.     }

  183.     /** Build an instance from all its derivatives.
  184.      * @param parameters number of free parameters
  185.      * @param order derivation order
  186.      * @param derivatives derivatives sorted according to
  187.      * {@link DSCompiler#getPartialDerivativeIndex(int...)}
  188.      * @throws DimensionMismatchException if derivatives array does not match the
  189.      * {@link DSCompiler#getSize() size} expected by the compiler.
  190.      * @throws NumberIsTooLargeException if order is too large.
  191.      * @see #getAllDerivatives()
  192.      */
  193.     public DerivativeStructure(final int parameters, final int order, final double ... derivatives) {
  194.         this(parameters, order);
  195.         if (derivatives.length != data.length) {
  196.             throw new DimensionMismatchException(derivatives.length, data.length);
  197.         }
  198.         System.arraycopy(derivatives, 0, data, 0, data.length);
  199.     }

  201.     /** Copy constructor.
  202.      * @param ds instance to copy
  203.      */
  204.     private DerivativeStructure(final DerivativeStructure ds) {
  205.         this.compiler = ds.compiler;
  206.         this.data     = ds.data.clone();
  207.     }

  209.     /** Get the number of free parameters.
  210.      * @return number of free parameters
  211.      */
  212.     public int getFreeParameters() {
  213.         return compiler.getFreeParameters();
  214.     }

  216.     /** Get the derivation order.
  217.      * @return derivation order
  218.      */
  219.     public int getOrder() {
  220.         return compiler.getOrder();
  221.     }

  223.     /** Create a constant compatible with instance order and number of parameters.
  224.      * <p>
  225.      * This method is a convenience factory method, it simply calls
  226.      * {@code new DerivativeStructure(getFreeParameters(), getOrder(), c)}
  227.      * </p>
  228.      * @param c value of the constant
  229.      * @return a constant compatible with instance order and number of parameters
  230.      * @see #DerivativeStructure(int, int, double)
  231.      * @since 3.3
  232.      */
  233.     public DerivativeStructure createConstant(final double c) {
  234.         return new DerivativeStructure(getFreeParameters(), getOrder(), c);
  235.     }

  237.     /** {@inheritDoc}
  238.      * @since 3.2
  239.      */
  240.     @Override
  241.     public double getReal() {
  242.         return data[0];
  243.     }

  245.     /** Get the value part of the derivative structure.
  246.      * @return value part of the derivative structure
  247.      * @see #getPartialDerivative(int...)
  248.      */
  249.     public double getValue() {
  250.         return data[0];
  251.     }

  253.     /** Get a partial derivative.
  254.      * @param orders derivation orders with respect to each variable (if all orders are 0,
  255.      * the value is returned)
  256.      * @return partial derivative
  257.      * @see #getValue()
  258.      * @throws DimensionMismatchException if the numbers of variables does not
  259.      * match the instance.
  260.      * @throws NumberIsTooLargeException if the sum of derivation orders is larger than
  261.      * the instance limits.
  262.      */
  263.     public double getPartialDerivative(final int ... orders) {
  264.         return data[compiler.getPartialDerivativeIndex(orders)];
  265.     }

  267.     /** Get all partial derivatives.
  268.      * @return a fresh copy of partial derivatives, in an array sorted according to
  269.      * {@link DSCompiler#getPartialDerivativeIndex(int...)}
  270.      */
  271.     public double[] getAllDerivatives() {
  272.         return data.clone();
  273.     }

  275.     /** {@inheritDoc}
  276.      * @since 3.2
  277.      */
  278.     @Override
  279.     public DerivativeStructure add(final double a) {
  280.         final DerivativeStructure ds = new DerivativeStructure(this);
  281.         ds.data[0] += a;
  282.         return ds;
  283.     }

  285.     /** {@inheritDoc}
  286.      * @throws DimensionMismatchException if number of free parameters
  287.      * or orders do not match.
  288.      */
  289.     @Override
  290.     public DerivativeStructure add(final DerivativeStructure a) {
  291.         compiler.checkCompatibility(a.compiler);
  292.         final DerivativeStructure ds = new DerivativeStructure(this);
  293.         compiler.add(data, 0, a.data, 0, ds.data, 0);
  294.         return ds;
  295.     }

  297.     /** {@inheritDoc}
  298.      * @since 3.2
  299.      */
  300.     @Override
  301.     public DerivativeStructure subtract(final double a) {
  302.         return add(-a);
  303.     }

  305.     /** {@inheritDoc}
  306.      * @throws DimensionMismatchException if number of free parameters
  307.      * or orders do not match
  308.      */
  309.     @Override
  310.     public DerivativeStructure subtract(final DerivativeStructure a) {
  311.         compiler.checkCompatibility(a.compiler);
  312.         final DerivativeStructure ds = new DerivativeStructure(this);
  313.         compiler.subtract(data, 0, a.data, 0, ds.data, 0);
  314.         return ds;
  315.     }

  317.     /** {@inheritDoc} */
  318.     @Override
  319.     public DerivativeStructure multiply(final int n) {
  320.         return multiply((double) n);
  321.     }

  323.     /** {@inheritDoc}
  324.      * @since 3.2
  325.      */
  326.     @Override
  327.     public DerivativeStructure multiply(final double a) {
  328.         final DerivativeStructure ds = new DerivativeStructure(this);
  329.         for (int i = 0; i < ds.data.length; ++i) {
  330.             ds.data[i] *= a;
  331.         }
  332.         return ds;
  333.     }

  335.     /** {@inheritDoc}
  336.      * @throws DimensionMismatchException if number of free parameters
  337.      * or orders do not match
  338.      */
  339.     @Override
  340.     public DerivativeStructure multiply(final DerivativeStructure a) {
  341.         compiler.checkCompatibility(a.compiler);
  342.         final DerivativeStructure result = new DerivativeStructure(compiler);
  343.         compiler.multiply(data, 0, a.data, 0, result.data, 0);
  344.         return result;
  345.     }

  347.     /** {@inheritDoc}
  348.      * @since 3.2
  349.      */
  350.     @Override
  351.     public DerivativeStructure divide(final double a) {
  352.         final DerivativeStructure ds = new DerivativeStructure(this);
  353.         for (int i = 0; i < ds.data.length; ++i) {
  354.             ds.data[i] /= a;
  355.         }
  356.         return ds;
  357.     }

  359.     /** {@inheritDoc}
  360.      * @throws DimensionMismatchException if number of free parameters
  361.      * or orders do not match
  362.      */
  363.     @Override
  364.     public DerivativeStructure divide(final DerivativeStructure a) {
  365.         compiler.checkCompatibility(a.compiler);
  366.         final DerivativeStructure result = new DerivativeStructure(compiler);
  367.         compiler.divide(data, 0, a.data, 0, result.data, 0);
  368.         return result;
  369.     }

  371.     /** {@inheritDoc} */
  372.     @Override
  373.     public DerivativeStructure remainder(final double a) {
  374.         final DerivativeStructure ds = new DerivativeStructure(this);
  375.         ds.data[0] = JdkMath.IEEEremainder(ds.data[0], a);
  376.         return ds;
  377.     }

  379.     /** {@inheritDoc}
  380.      * @throws DimensionMismatchException if number of free parameters
  381.      * or orders do not match
  382.      * @since 3.2
  383.      */
  384.     @Override
  385.     public DerivativeStructure remainder(final DerivativeStructure a) {
  386.         compiler.checkCompatibility(a.compiler);
  387.         final DerivativeStructure result = new DerivativeStructure(compiler);
  388.         compiler.remainder(data, 0, a.data, 0, result.data, 0);
  389.         return result;
  390.     }

  392.     /** {@inheritDoc} */
  393.     @Override
  394.     public DerivativeStructure negate() {
  395.         final DerivativeStructure ds = new DerivativeStructure(compiler);
  396.         for (int i = 0; i < ds.data.length; ++i) {
  397.             ds.data[i] = -data[i];
  398.         }
  399.         return ds;
  400.     }

  402.     /** {@inheritDoc}
  403.      * @since 3.2
  404.      */
  405.     @Override
  406.     public DerivativeStructure abs() {
  407.         if (Double.doubleToLongBits(data[0]) < 0) {
  408.             // we use the bits representation to also handle -0.0
  409.             return negate();
  410.         } else {
  411.             return this;
  412.         }
  413.     }

  415.     /** {@inheritDoc}
  416.      * @since 3.2
  417.      */
  418.     @Override
  419.     public DerivativeStructure ceil() {
  420.         return new DerivativeStructure(compiler.getFreeParameters(),
  421.                                        compiler.getOrder(),
  422.                                        JdkMath.ceil(data[0]));
  423.     }

  425.     /** {@inheritDoc}
  426.      * @since 3.2
  427.      */
  428.     @Override
  429.     public DerivativeStructure floor() {
  430.         return new DerivativeStructure(compiler.getFreeParameters(),
  431.                                        compiler.getOrder(),
  432.                                        JdkMath.floor(data[0]));
  433.     }

  435.     /** {@inheritDoc}
  436.      * @since 3.2
  437.      */
  438.     @Override
  439.     public DerivativeStructure rint() {
  440.         return new DerivativeStructure(compiler.getFreeParameters(),
  441.                                        compiler.getOrder(),
  442.                                        JdkMath.rint(data[0]));
  443.     }

  445.     /** {@inheritDoc} */
  446.     @Override
  447.     public long round() {
  448.         return JdkMath.round(data[0]);
  449.     }

  451.     /** {@inheritDoc}
  452.      * @since 3.2
  453.      */
  454.     @Override
  455.     public DerivativeStructure signum() {
  456.         return new DerivativeStructure(compiler.getFreeParameters(),
  457.                                        compiler.getOrder(),
  458.                                        JdkMath.signum(data[0]));
  459.     }

  461.     /** {@inheritDoc}
  462.      * @since 3.2
  463.      */
  464.     @Override
  465.     public DerivativeStructure copySign(final DerivativeStructure sign){
  466.         long m = Double.doubleToLongBits(data[0]);
  467.         long s = Double.doubleToLongBits(sign.data[0]);
  468.         if ((m >= 0 && s >= 0) || (m < 0 && s < 0)) { // Sign is currently OK
  469.             return this;
  470.         }
  471.         return negate(); // flip sign
  472.     }

  474.     /** {@inheritDoc}
  475.      * @since 3.2
  476.      */
  477.     @Override
  478.     public DerivativeStructure copySign(final double sign) {
  479.         long m = Double.doubleToLongBits(data[0]);
  480.         long s = Double.doubleToLongBits(sign);
  481.         if ((m >= 0 && s >= 0) || (m < 0 && s < 0)) { // Sign is currently OK
  482.             return this;
  483.         }
  484.         return negate(); // flip sign
  485.     }

  487.     /**
  488.      * Return the exponent of the instance value, removing the bias.
  489.      * <p>
  490.      * For double numbers of the form 2<sup>x</sup>, the unbiased
  491.      * exponent is exactly x.
  492.      * </p>
  493.      * @return exponent for instance in IEEE754 representation, without bias
  494.      */
  495.     public int getExponent() {
  496.         return JdkMath.getExponent(data[0]);
  497.     }

  499.     /** {@inheritDoc}
  500.      * @since 3.2
  501.      */
  502.     @Override
  503.     public DerivativeStructure scalb(final int n) {
  504.         final DerivativeStructure ds = new DerivativeStructure(compiler);
  505.         for (int i = 0; i < ds.data.length; ++i) {
  506.             ds.data[i] = JdkMath.scalb(data[i], n);
  507.         }
  508.         return ds;
  509.     }

  511.     /** {@inheritDoc}
  512.      * @throws DimensionMismatchException if number of free parameters
  513.      * or orders do not match
  514.      * @since 3.2
  515.      */
  516.     @Override
  517.     public DerivativeStructure hypot(final DerivativeStructure y) {
  518.         compiler.checkCompatibility(y.compiler);

  520.         if (Double.isInfinite(data[0]) || Double.isInfinite(y.data[0])) {
  521.             return new DerivativeStructure(compiler.getFreeParameters(),
  522.                                            compiler.getFreeParameters(),
  523.                                            Double.POSITIVE_INFINITY);
  524.         } else if (Double.isNaN(data[0]) || Double.isNaN(y.data[0])) {
  525.             return new DerivativeStructure(compiler.getFreeParameters(),
  526.                                            compiler.getFreeParameters(),
  527.                                            Double.NaN);
  528.         } else {

  530.             final int expX = getExponent();
  531.             final int expY = y.getExponent();
  532.             if (expX > expY + 27) {
  533.                 // y is neglectible with respect to x
  534.                 return abs();
  535.             } else if (expY > expX + 27) {
  536.                 // x is neglectible with respect to y
  537.                 return y.abs();
  538.             } else {

  540.                 // find an intermediate scale to avoid both overflow and underflow
  541.                 final int middleExp = (expX + expY) / 2;

  543.                 // scale parameters without losing precision
  544.                 final DerivativeStructure scaledX = scalb(-middleExp);
  545.                 final DerivativeStructure scaledY = y.scalb(-middleExp);

  547.                 // compute scaled hypotenuse
  548.                 final DerivativeStructure scaledH =
  549.                         scaledX.multiply(scaledX).add(scaledY.multiply(scaledY)).sqrt();

  551.                 // remove scaling
  552.                 return scaledH.scalb(middleExp);
  553.             }
  554.         }
  555.     }

  557.     /**
  558.      * Returns the hypotenuse of a triangle with sides {@code x} and {@code y}
  559.      * - sqrt(<i>x</i><sup>2</sup>&nbsp;+<i>y</i><sup>2</sup>)
  560.      * avoiding intermediate overflow or underflow.
  561.      *
  562.      * <ul>
  563.      * <li> If either argument is infinite, then the result is positive infinity.</li>
  564.      * <li> else, if either argument is NaN then the result is NaN.</li>
  565.      * </ul>
  566.      *
  567.      * @param x a value
  568.      * @param y a value
  569.      * @return sqrt(<i>x</i><sup>2</sup>&nbsp;+<i>y</i><sup>2</sup>)
  570.      * @throws DimensionMismatchException if number of free parameters
  571.      * or orders do not match
  572.      * @since 3.2
  573.      */
  574.     public static DerivativeStructure hypot(final DerivativeStructure x, final DerivativeStructure y) {
  575.         return x.hypot(y);
  576.     }

  578.     /** Compute composition of the instance by a univariate function.
  579.      * @param f array of value and derivatives of the function at
  580.      * the current point (i.e. [f({@link #getValue()}),
  581.      * f'({@link #getValue()}), f''({@link #getValue()})...]).
  582.      * @return f(this)
  583.      * @throws DimensionMismatchException if the number of derivatives
  584.      * in the array is not equal to {@link #getOrder() order} + 1
  585.      */
  586.     public DerivativeStructure compose(final double ... f) {
  587.         if (f.length != getOrder() + 1) {
  588.             throw new DimensionMismatchException(f.length, getOrder() + 1);
  589.         }
  590.         final DerivativeStructure result = new DerivativeStructure(compiler);
  591.         compiler.compose(data, 0, f, result.data, 0);
  592.         return result;
  593.     }

  595.     /** {@inheritDoc} */
  596.     @Override
  597.     public DerivativeStructure reciprocal() {
  598.         final DerivativeStructure result = new DerivativeStructure(compiler);
  599.         compiler.pow(data, 0, -1, result.data, 0);
  600.         return result;
  601.     }

  603.     /** {@inheritDoc}
  604.      * @since 3.2
  605.      */
  606.     @Override
  607.     public DerivativeStructure sqrt() {
  608.         return rootN(2);
  609.     }

  611.     /** {@inheritDoc}
  612.      * @since 3.2
  613.      */
  614.     @Override
  615.     public DerivativeStructure cbrt() {
  616.         return rootN(3);
  617.     }

  619.     /** {@inheritDoc}
  620.      * @since 3.2
  621.      */
  622.     @Override
  623.     public DerivativeStructure rootN(final int n) {
  624.         final DerivativeStructure result = new DerivativeStructure(compiler);
  625.         compiler.rootN(data, 0, n, result.data, 0);
  626.         return result;
  627.     }

  629.     /** {@inheritDoc} */
  630.     @Override
  631.     public Field<DerivativeStructure> getField() {
  632.         return new Field<DerivativeStructure>() {

  634.             /** {@inheritDoc} */
  635.             @Override
  636.             public DerivativeStructure getZero() {
  637.                 return new DerivativeStructure(compiler.getFreeParameters(), compiler.getOrder(), 0.0);
  638.             }

  640.             /** {@inheritDoc} */
  641.             @Override
  642.             public DerivativeStructure getOne() {
  643.                 return new DerivativeStructure(compiler.getFreeParameters(), compiler.getOrder(), 1.0);
  644.             }

  646.             /** {@inheritDoc} */
  647.             @Override
  648.             public Class<? extends FieldElement<DerivativeStructure>> getRuntimeClass() {
  649.                 return DerivativeStructure.class;
  650.             }
  651.         };
  652.     }

  654.     /** Compute a<sup>x</sup> where a is a double and x a {@link DerivativeStructure}.
  655.      * @param a number to exponentiate
  656.      * @param x power to apply
  657.      * @return a<sup>x</sup>
  658.      * @since 3.3
  659.      */
  660.     public static DerivativeStructure pow(final double a, final DerivativeStructure x) {
  661.         final DerivativeStructure result = new DerivativeStructure(x.compiler);
  662.         x.compiler.pow(a, x.data, 0, result.data, 0);
  663.         return result;
  664.     }

  666.     /** {@inheritDoc}
  667.      * @since 3.2
  668.      */
  669.     @Override
  670.     public DerivativeStructure pow(final double p) {
  671.         final DerivativeStructure result = new DerivativeStructure(compiler);
  672.         compiler.pow(data, 0, p, result.data, 0);
  673.         return result;
  674.     }

  676.     /** {@inheritDoc}
  677.      * @since 3.2
  678.      */
  679.     @Override
  680.     public DerivativeStructure pow(final int n) {
  681.         final DerivativeStructure result = new DerivativeStructure(compiler);
  682.         compiler.pow(data, 0, n, result.data, 0);
  683.         return result;
  684.     }

  686.     /** {@inheritDoc}
  687.      * @throws DimensionMismatchException if number of free parameters
  688.      * or orders do not match
  689.      * @since 3.2
  690.      */
  691.     @Override
  692.     public DerivativeStructure pow(final DerivativeStructure e) {
  693.         compiler.checkCompatibility(e.compiler);
  694.         final DerivativeStructure result = new DerivativeStructure(compiler);
  695.         compiler.pow(data, 0, e.data, 0, result.data, 0);
  696.         return result;
  697.     }

  699.     /** {@inheritDoc}
  700.      * @since 3.2
  701.      */
  702.     @Override
  703.     public DerivativeStructure exp() {
  704.         final DerivativeStructure result = new DerivativeStructure(compiler);
  705.         compiler.exp(data, 0, result.data, 0);
  706.         return result;
  707.     }

  709.     /** {@inheritDoc}
  710.      * @since 3.2
  711.      */
  712.     @Override
  713.     public DerivativeStructure expm1() {
  714.         final DerivativeStructure result = new DerivativeStructure(compiler);
  715.         compiler.expm1(data, 0, result.data, 0);
  716.         return result;
  717.     }

  719.     /** {@inheritDoc}
  720.      * @since 3.2
  721.      */
  722.     @Override
  723.     public DerivativeStructure log() {
  724.         final DerivativeStructure result = new DerivativeStructure(compiler);
  725.         compiler.log(data, 0, result.data, 0);
  726.         return result;
  727.     }

  729.     /** {@inheritDoc}
  730.      * @since 3.2
  731.      */
  732.     @Override
  733.     public DerivativeStructure log1p() {
  734.         final DerivativeStructure result = new DerivativeStructure(compiler);
  735.         compiler.log1p(data, 0, result.data, 0);
  736.         return result;
  737.     }

  739.     /** Base 10 logarithm.
  740.      * @return base 10 logarithm of the instance
  741.      */
  742.     @Override
  743.     public DerivativeStructure log10() {
  744.         final DerivativeStructure result = new DerivativeStructure(compiler);
  745.         compiler.log10(data, 0, result.data, 0);
  746.         return result;
  747.     }

  749.     /** {@inheritDoc}
  750.      * @since 3.2
  751.      */
  752.     @Override
  753.     public DerivativeStructure cos() {
  754.         final DerivativeStructure result = new DerivativeStructure(compiler);
  755.         compiler.cos(data, 0, result.data, 0);
  756.         return result;
  757.     }

  759.     /** {@inheritDoc}
  760.      * @since 3.2
  761.      */
  762.     @Override
  763.     public DerivativeStructure sin() {
  764.         final DerivativeStructure result = new DerivativeStructure(compiler);
  765.         compiler.sin(data, 0, result.data, 0);
  766.         return result;
  767.     }

  769.     /** {@inheritDoc}
  770.      * @since 3.2
  771.      */
  772.     @Override
  773.     public DerivativeStructure tan() {
  774.         final DerivativeStructure result = new DerivativeStructure(compiler);
  775.         compiler.tan(data, 0, result.data, 0);
  776.         return result;
  777.     }

  779.     /** {@inheritDoc}
  780.      * @since 3.2
  781.      */
  782.     @Override
  783.     public DerivativeStructure acos() {
  784.         final DerivativeStructure result = new DerivativeStructure(compiler);
  785.         compiler.acos(data, 0, result.data, 0);
  786.         return result;
  787.     }

  789.     /** {@inheritDoc}
  790.      * @since 3.2
  791.      */
  792.     @Override
  793.     public DerivativeStructure asin() {
  794.         final DerivativeStructure result = new DerivativeStructure(compiler);
  795.         compiler.asin(data, 0, result.data, 0);
  796.         return result;
  797.     }

  799.     /** {@inheritDoc}
  800.      * @since 3.2
  801.      */
  802.     @Override
  803.     public DerivativeStructure atan() {
  804.         final DerivativeStructure result = new DerivativeStructure(compiler);
  805.         compiler.atan(data, 0, result.data, 0);
  806.         return result;
  807.     }

  809.     /** {@inheritDoc}
  810.      * @since 3.2
  811.      */
  812.     @Override
  813.     public DerivativeStructure atan2(final DerivativeStructure x) {
  814.         compiler.checkCompatibility(x.compiler);
  815.         final DerivativeStructure result = new DerivativeStructure(compiler);
  816.         compiler.atan2(data, 0, x.data, 0, result.data, 0);
  817.         return result;
  818.     }

  820.     /** Two arguments arc tangent operation.
  821.      * @param y first argument of the arc tangent
  822.      * @param x second argument of the arc tangent
  823.      * @return atan2(y, x)
  824.      * @throws DimensionMismatchException if number of free parameters
  825.      * or orders do not match
  826.      * @since 3.2
  827.      */
  828.     public static DerivativeStructure atan2(final DerivativeStructure y, final DerivativeStructure x) {
  829.         return y.atan2(x);
  830.     }

  832.     /** {@inheritDoc}
  833.      * @since 3.2
  834.      */
  835.     @Override
  836.     public DerivativeStructure cosh() {
  837.         final DerivativeStructure result = new DerivativeStructure(compiler);
  838.         compiler.cosh(data, 0, result.data, 0);
  839.         return result;
  840.     }

  842.     /** {@inheritDoc}
  843.      * @since 3.2
  844.      */
  845.     @Override
  846.     public DerivativeStructure sinh() {
  847.         final DerivativeStructure result = new DerivativeStructure(compiler);
  848.         compiler.sinh(data, 0, result.data, 0);
  849.         return result;
  850.     }

  852.     /** {@inheritDoc}
  853.      * @since 3.2
  854.      */
  855.     @Override
  856.     public DerivativeStructure tanh() {
  857.         final DerivativeStructure result = new DerivativeStructure(compiler);
  858.         compiler.tanh(data, 0, result.data, 0);
  859.         return result;
  860.     }

  862.     /** {@inheritDoc}
  863.      * @since 3.2
  864.      */
  865.     @Override
  866.     public DerivativeStructure acosh() {
  867.         final DerivativeStructure result = new DerivativeStructure(compiler);
  868.         compiler.acosh(data, 0, result.data, 0);
  869.         return result;
  870.     }

  872.     /** {@inheritDoc}
  873.      * @since 3.2
  874.      */
  875.     @Override
  876.     public DerivativeStructure asinh() {
  877.         final DerivativeStructure result = new DerivativeStructure(compiler);
  878.         compiler.asinh(data, 0, result.data, 0);
  879.         return result;
  880.     }

  882.     /** {@inheritDoc}
  883.      * @since 3.2
  884.      */
  885.     @Override
  886.     public DerivativeStructure atanh() {
  887.         final DerivativeStructure result = new DerivativeStructure(compiler);
  888.         compiler.atanh(data, 0, result.data, 0);
  889.         return result;
  890.     }

  892.     /** Convert radians to degrees, with error of less than 0.5 ULP.
  893.      *  @return instance converted into degrees
  894.      */
  895.     public DerivativeStructure toDegrees() {
  896.         final DerivativeStructure ds = new DerivativeStructure(compiler);
  897.         for (int i = 0; i < ds.data.length; ++i) {
  898.             ds.data[i] = JdkMath.toDegrees(data[i]);
  899.         }
  900.         return ds;
  901.     }

  903.     /** Convert degrees to radians, with error of less than 0.5 ULP.
  904.      *  @return instance converted into radians
  905.      */
  906.     public DerivativeStructure toRadians() {
  907.         final DerivativeStructure ds = new DerivativeStructure(compiler);
  908.         for (int i = 0; i < ds.data.length; ++i) {
  909.             ds.data[i] = JdkMath.toRadians(data[i]);
  910.         }
  911.         return ds;
  912.     }

  914.     /** Evaluate Taylor expansion a derivative structure.
  915.      * @param delta parameters offsets (&Delta;x, &Delta;y, ...)
  916.      * @return value of the Taylor expansion at x + &Delta;x, y + &Delta;y, ...
  917.      * @throws org.apache.commons.math4.legacy.exception.MathArithmeticException
  918.      * if factorials becomes too large
  919.      */
  920.     public double taylor(final double ... delta) {
  921.         return compiler.taylor(data, 0, delta);
  922.     }

  924.     /** {@inheritDoc}
  925.      * @throws DimensionMismatchException if number of free parameters
  926.      * or orders do not match
  927.      * @since 3.2
  928.      */
  929.     @Override
  930.     public DerivativeStructure linearCombination(final DerivativeStructure[] a, final DerivativeStructure[] b) {
  931.         // compute an accurate value, taking care of cancellations
  932.         final double[] aDouble = new double[a.length];
  933.         for (int i = 0; i < a.length; ++i) {
  934.             aDouble[i] = a[i].getValue();
  935.         }
  936.         final double[] bDouble = new double[b.length];
  937.         for (int i = 0; i < b.length; ++i) {
  938.             bDouble[i] = b[i].getValue();
  939.         }
  940.         final double accurateValue = Sum.ofProducts(aDouble, bDouble).getAsDouble();

  942.         // compute a simple value, with all partial derivatives
  943.         DerivativeStructure simpleValue = a[0].getField().getZero();
  944.         for (int i = 0; i < a.length; ++i) {
  945.             simpleValue = simpleValue.add(a[i].multiply(b[i]));
  946.         }

  948.         // create a result with accurate value and all derivatives (not necessarily as accurate as the value)
  949.         final double[] all = simpleValue.getAllDerivatives();
  950.         all[0] = accurateValue;
  951.         return new DerivativeStructure(simpleValue.getFreeParameters(), simpleValue.getOrder(), all);
  952.     }

  954.     /** {@inheritDoc}
  955.      * @throws DimensionMismatchException if number of free parameters
  956.      * or orders do not match
  957.      * @since 3.2
  958.      */
  959.     @Override
  960.     public DerivativeStructure linearCombination(final double[] a, final DerivativeStructure[] b) {
  961.         // compute an accurate value, taking care of cancellations
  962.         final double[] bDouble = new double[b.length];
  963.         for (int i = 0; i < b.length; ++i) {
  964.             bDouble[i] = b[i].getValue();
  965.         }
  966.         final double accurateValue = Sum.ofProducts(a, bDouble).getAsDouble();

  968.         // compute a simple value, with all partial derivatives
  969.         DerivativeStructure simpleValue = b[0].getField().getZero();
  970.         for (int i = 0; i < a.length; ++i) {
  971.             simpleValue = simpleValue.add(b[i].multiply(a[i]));
  972.         }

  974.         // create a result with accurate value and all derivatives (not necessarily as accurate as the value)
  975.         final double[] all = simpleValue.getAllDerivatives();
  976.         all[0] = accurateValue;
  977.         return new DerivativeStructure(simpleValue.getFreeParameters(), simpleValue.getOrder(), all);
  978.     }

  980.     /** {@inheritDoc}
  981.      * @throws DimensionMismatchException if number of free parameters
  982.      * or orders do not match
  983.      * @since 3.2
  984.      */
  985.     @Override
  986.     public DerivativeStructure linearCombination(final DerivativeStructure a1, final DerivativeStructure b1,
  987.                                                  final DerivativeStructure a2, final DerivativeStructure b2) {
  988.         // compute an accurate value, taking care of cancellations
  989.         final double accurateValue = Sum.create()
  990.             .addProduct(a1.getValue(), b1.getValue())
  991.             .addProduct(a2.getValue(), b2.getValue()).getAsDouble();

  993.         // compute a simple value, with all partial derivatives
  994.         final DerivativeStructure simpleValue = a1.multiply(b1).add(a2.multiply(b2));

  996.         // create a result with accurate value and all derivatives (not necessarily as accurate as the value)
  997.         final double[] all = simpleValue.getAllDerivatives();
  998.         all[0] = accurateValue;
  999.         return new DerivativeStructure(getFreeParameters(), getOrder(), all);
  1000.     }

  1002.     /** {@inheritDoc}
  1003.      * @throws DimensionMismatchException if number of free parameters
  1004.      * or orders do not match
  1005.      * @since 3.2
  1006.      */
  1007.     @Override
  1008.     public DerivativeStructure linearCombination(final double a1, final DerivativeStructure b1,
  1009.                                                  final double a2, final DerivativeStructure b2) {
  1010.         // compute an accurate value, taking care of cancellations
  1011.         final double accurateValue = Sum.create()
  1012.             .addProduct(a1, b1.getValue())
  1013.             .addProduct(a2, b2.getValue()).getAsDouble();

  1015.         // compute a simple value, with all partial derivatives
  1016.         final DerivativeStructure simpleValue = b1.multiply(a1).add(b2.multiply(a2));

  1018.         // create a result with accurate value and all derivatives (not necessarily as accurate as the value)
  1019.         final double[] all = simpleValue.getAllDerivatives();
  1020.         all[0] = accurateValue;
  1021.         return new DerivativeStructure(getFreeParameters(), getOrder(), all);
  1022.     }

  1024.     /** {@inheritDoc}
  1025.      * @throws DimensionMismatchException if number of free parameters
  1026.      * or orders do not match
  1027.      * @since 3.2
  1028.      */
  1029.     @Override
  1030.     public DerivativeStructure linearCombination(final DerivativeStructure a1, final DerivativeStructure b1,
  1031.                                                  final DerivativeStructure a2, final DerivativeStructure b2,
  1032.                                                  final DerivativeStructure a3, final DerivativeStructure b3) {
  1033.         // compute an accurate value, taking care of cancellations
  1034.         final double accurateValue = Sum.create()
  1035.             .addProduct(a1.getValue(), b1.getValue())
  1036.             .addProduct(a2.getValue(), b2.getValue())
  1037.             .addProduct(a3.getValue(), b3.getValue()).getAsDouble();

  1039.         // compute a simple value, with all partial derivatives
  1040.         final DerivativeStructure simpleValue = a1.multiply(b1).add(a2.multiply(b2)).add(a3.multiply(b3));

  1042.         // create a result with accurate value and all derivatives (not necessarily as accurate as the value)
  1043.         final double[] all = simpleValue.getAllDerivatives();
  1044.         all[0] = accurateValue;
  1045.         return new DerivativeStructure(getFreeParameters(), getOrder(), all);
  1046.     }

  1048.     /** {@inheritDoc}
  1049.      * @throws DimensionMismatchException if number of free parameters
  1050.      * or orders do not match
  1051.      * @since 3.2
  1052.      */
  1053.     @Override
  1054.     public DerivativeStructure linearCombination(final double a1, final DerivativeStructure b1,
  1055.                                                  final double a2, final DerivativeStructure b2,
  1056.                                                  final double a3, final DerivativeStructure b3) {
  1057.         // compute an accurate value, taking care of cancellations
  1058.         final double accurateValue = Sum.create()
  1059.             .addProduct(a1, b1.getValue())
  1060.             .addProduct(a2, b2.getValue())
  1061.             .addProduct(a3, b3.getValue()).getAsDouble();

  1063.         // compute a simple value, with all partial derivatives
  1064.         final DerivativeStructure simpleValue = b1.multiply(a1).add(b2.multiply(a2)).add(b3.multiply(a3));

  1066.         // create a result with accurate value and all derivatives (not necessarily as accurate as the value)
  1067.         final double[] all = simpleValue.getAllDerivatives();
  1068.         all[0] = accurateValue;
  1069.         return new DerivativeStructure(getFreeParameters(), getOrder(), all);
  1070.     }

  1072.     /** {@inheritDoc}
  1073.      * @throws DimensionMismatchException if number of free parameters
  1074.      * or orders do not match
  1075.      * @since 3.2
  1076.      */
  1077.     @Override
  1078.     public DerivativeStructure linearCombination(final DerivativeStructure a1, final DerivativeStructure b1,
  1079.                                                  final DerivativeStructure a2, final DerivativeStructure b2,
  1080.                                                  final DerivativeStructure a3, final DerivativeStructure b3,
  1081.                                                  final DerivativeStructure a4, final DerivativeStructure b4) {
  1082.         // compute an accurate value, taking care of cancellations
  1083.         final double accurateValue = Sum.create()
  1084.             .addProduct(a1.getValue(), b1.getValue())
  1085.             .addProduct(a2.getValue(), b2.getValue())
  1086.             .addProduct(a3.getValue(), b3.getValue())
  1087.             .addProduct(a4.getValue(), b4.getValue()).getAsDouble();

  1089.         // compute a simple value, with all partial derivatives
  1090.         final DerivativeStructure simpleValue = a1.multiply(b1).add(a2.multiply(b2)).add(a3.multiply(b3)).add(a4.multiply(b4));

  1092.         // create a result with accurate value and all derivatives (not necessarily as accurate as the value)
  1093.         final double[] all = simpleValue.getAllDerivatives();
  1094.         all[0] = accurateValue;
  1095.         return new DerivativeStructure(getFreeParameters(), getOrder(), all);
  1096.     }

  1098.     /** {@inheritDoc}
  1099.      * @throws DimensionMismatchException if number of free parameters
  1100.      * or orders do not match
  1101.      * @since 3.2
  1102.      */
  1103.     @Override
  1104.     public DerivativeStructure linearCombination(final double a1, final DerivativeStructure b1,
  1105.                                                  final double a2, final DerivativeStructure b2,
  1106.                                                  final double a3, final DerivativeStructure b3,
  1107.                                                  final double a4, final DerivativeStructure b4) {
  1108.         // compute an accurate value, taking care of cancellations
  1109.         final double accurateValue = Sum.create()
  1110.             .addProduct(a1, b1.getValue())
  1111.             .addProduct(a2, b2.getValue())
  1112.             .addProduct(a3, b3.getValue())
  1113.             .addProduct(a4, b4.getValue()).getAsDouble();

  1115.         // compute a simple value, with all partial derivatives
  1116.         final DerivativeStructure simpleValue = b1.multiply(a1).add(b2.multiply(a2)).add(b3.multiply(a3)).add(b4.multiply(a4));

  1118.         // create a result with accurate value and all derivatives (not necessarily as accurate as the value)
  1119.         final double[] all = simpleValue.getAllDerivatives();
  1120.         all[0] = accurateValue;
  1121.         return new DerivativeStructure(getFreeParameters(), getOrder(), all);
  1122.     }

  1124.     /**
  1125.      * Test for the equality of two derivative structures.
  1126.      * <p>
  1127.      * Derivative structures are considered equal if they have the same number
  1128.      * of free parameters, the same derivation order, and the same derivatives.
  1129.      * </p>
  1130.      * @param other Object to test for equality to this
  1131.      * @return true if two derivative structures are equal
  1132.      * @since 3.2
  1133.      */
  1134.     @Override
  1135.     public boolean equals(Object other) {

  1137.         if (this == other) {
  1138.             return true;
  1139.         }

  1141.         if (other instanceof DerivativeStructure) {
  1142.             final DerivativeStructure rhs = (DerivativeStructure)other;
  1143.             return getFreeParameters() == rhs.getFreeParameters() &&
  1144.                    getOrder() == rhs.getOrder() &&
  1145.                    MathArrays.equals(data, rhs.data);
  1146.         }

  1148.         return false;
  1149.     }

  1151.     /**
  1152.      * Get a hashCode for the derivative structure.
  1153.      * @return a hash code value for this object
  1154.      * @since 3.2
  1155.      */
  1156.     @Override
  1157.     public int hashCode() {
  1158.         return 227 + 229 * getFreeParameters() + 233 * getOrder() + 239 * Arrays.hashCode(data);
  1159.     }
  1160. }
Created with JaCoCo 0.8.10.202304240956Code Coverage Report for Apache Commons Math 4.0-SNAPSHOT