View Javadoc
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.math3.analysis.differentiation;
18  
19  import java.io.Serializable;
20  
21  import org.apache.commons.math3.RealFieldElement;
22  import org.apache.commons.math3.Field;
23  import org.apache.commons.math3.FieldElement;
24  import org.apache.commons.math3.exception.DimensionMismatchException;
25  import org.apache.commons.math3.exception.MathArithmeticException;
26  import org.apache.commons.math3.exception.NumberIsTooLargeException;
27  import org.apache.commons.math3.util.FastMath;
28  import org.apache.commons.math3.util.MathArrays;
29  import org.apache.commons.math3.util.MathUtils;
30  
31  /** Class representing both the value and the differentials of a function.
32   * <p>This class is the workhorse of the differentiation package.</p>
33   * <p>This class is an implementation of the extension to Rall's
34   * numbers described in Dan Kalman's paper <a
35   * href="http://www1.american.edu/cas/mathstat/People/kalman/pdffiles/mmgautodiff.pdf">Doubly
36   * Recursive Multivariate Automatic Differentiation</a>, Mathematics Magazine, vol. 75,
37   * no. 3, June 2002.</p>. Rall's numbers are an extension to the real numbers used
38   * throughout mathematical expressions; they hold the derivative together with the
39   * value of a function. Dan Kalman's derivative structures hold all partial derivatives
40   * up to any specified order, with respect to any number of free parameters. Rall's
41   * numbers therefore can be seen as derivative structures for order one derivative and
42   * one free parameter, and real numbers can be seen as derivative structures with zero
43   * order derivative and no free parameters.</p>
44   * <p>{@link DerivativeStructure} instances can be used directly thanks to
45   * the arithmetic operators to the mathematical functions provided as
46   * methods by this class (+, -, *, /, %, sin, cos ...).</p>
47   * <p>Implementing complex expressions by hand using these classes is
48   * a tedious and error-prone task but has the advantage of having no limitation
49   * on the derivation order despite no requiring users to compute the derivatives by
50   * themselves. Implementing complex expression can also be done by developing computation
51   * code using standard primitive double values and to use {@link
52   * UnivariateFunctionDifferentiator differentiators} to create the {@link
53   * DerivativeStructure}-based instances. This method is simpler but may be limited in
54   * the accuracy and derivation orders and may be computationally intensive (this is
55   * typically the case for {@link FiniteDifferencesDifferentiator finite differences
56   * differentiator}.</p>
57   * <p>Instances of this class are guaranteed to be immutable.</p>
58   * @see DSCompiler
59   * @since 3.1
60   */
61  public class DerivativeStructure implements RealFieldElement<DerivativeStructure>, Serializable {
62  
63      /** Serializable UID. */
64      private static final long serialVersionUID = 20120730L;
65  
66      /** Compiler for the current dimensions. */
67      private transient DSCompiler compiler;
68  
69      /** Combined array holding all values. */
70      private final double[] data;
71  
72      /** Build an instance with all values and derivatives set to 0.
73       * @param compiler compiler to use for computation
74       */
75      private DerivativeStructure(final DSCompiler compiler) {
76          this.compiler = compiler;
77          this.data     = new double[compiler.getSize()];
78      }
79  
80      /** Build an instance with all values and derivatives set to 0.
81       * @param parameters number of free parameters
82       * @param order derivation order
83       * @throws NumberIsTooLargeException if order is too large
84       */
85      public DerivativeStructure(final int parameters, final int order)
86          throws NumberIsTooLargeException {
87          this(DSCompiler.getCompiler(parameters, order));
88      }
89  
90      /** Build an instance representing a constant value.
91       * @param parameters number of free parameters
92       * @param order derivation order
93       * @param value value of the constant
94       * @throws NumberIsTooLargeException if order is too large
95       * @see #DerivativeStructure(int, int, int, double)
96       */
97      public DerivativeStructure(final int parameters, final int order, final double value)
98          throws NumberIsTooLargeException {
99          this(parameters, order);
100         this.data[0] = value;
101     }
102 
103     /** Build an instance representing a variable.
104      * <p>Instances built using this constructor are considered
105      * to be the free variables with respect to which differentials
106      * are computed. As such, their differential with respect to
107      * themselves is +1.</p>
108      * @param parameters number of free parameters
109      * @param order derivation order
110      * @param index index of the variable (from 0 to {@code parameters - 1})
111      * @param value value of the variable
112      * @exception NumberIsTooLargeException if {@code index >= parameters}.
113      * @see #DerivativeStructure(int, int, double)
114      */
115     public DerivativeStructure(final int parameters, final int order,
116                                final int index, final double value)
117         throws NumberIsTooLargeException {
118         this(parameters, order, value);
119 
120         if (index >= parameters) {
121             throw new NumberIsTooLargeException(index, parameters, false);
122         }
123 
124         if (order > 0) {
125             // the derivative of the variable with respect to itself is 1.
126             data[DSCompiler.getCompiler(index, order).getSize()] = 1.0;
127         }
128 
129     }
130 
131     /** Linear combination constructor.
132      * The derivative structure built will be a1 * ds1 + a2 * ds2
133      * @param a1 first scale factor
134      * @param ds1 first base (unscaled) derivative structure
135      * @param a2 second scale factor
136      * @param ds2 second base (unscaled) derivative structure
137      * @exception DimensionMismatchException if number of free parameters or orders are inconsistent
138      */
139     public DerivativeStructure(final double a1, final DerivativeStructure ds1,
140                                final double a2, final DerivativeStructure ds2)
141         throws DimensionMismatchException {
142         this(ds1.compiler);
143         compiler.checkCompatibility(ds2.compiler);
144         compiler.linearCombination(a1, ds1.data, 0, a2, ds2.data, 0, data, 0);
145     }
146 
147     /** Linear combination constructor.
148      * The derivative structure built will be a1 * ds1 + a2 * ds2 + a3 * ds3
149      * @param a1 first scale factor
150      * @param ds1 first base (unscaled) derivative structure
151      * @param a2 second scale factor
152      * @param ds2 second base (unscaled) derivative structure
153      * @param a3 third scale factor
154      * @param ds3 third base (unscaled) derivative structure
155      * @exception DimensionMismatchException if number of free parameters or orders are inconsistent
156      */
157     public DerivativeStructure(final double a1, final DerivativeStructure ds1,
158                                final double a2, final DerivativeStructure ds2,
159                                final double a3, final DerivativeStructure ds3)
160         throws DimensionMismatchException {
161         this(ds1.compiler);
162         compiler.checkCompatibility(ds2.compiler);
163         compiler.checkCompatibility(ds3.compiler);
164         compiler.linearCombination(a1, ds1.data, 0, a2, ds2.data, 0, a3, ds3.data, 0, data, 0);
165     }
166 
167     /** Linear combination constructor.
168      * The derivative structure built will be a1 * ds1 + a2 * ds2 + a3 * ds3 + a4 * ds4
169      * @param a1 first scale factor
170      * @param ds1 first base (unscaled) derivative structure
171      * @param a2 second scale factor
172      * @param ds2 second base (unscaled) derivative structure
173      * @param a3 third scale factor
174      * @param ds3 third base (unscaled) derivative structure
175      * @param a4 fourth scale factor
176      * @param ds4 fourth base (unscaled) derivative structure
177      * @exception DimensionMismatchException if number of free parameters or orders are inconsistent
178      */
179     public DerivativeStructure(final double a1, final DerivativeStructure ds1,
180                                final double a2, final DerivativeStructure ds2,
181                                final double a3, final DerivativeStructure ds3,
182                                final double a4, final DerivativeStructure ds4)
183         throws DimensionMismatchException {
184         this(ds1.compiler);
185         compiler.checkCompatibility(ds2.compiler);
186         compiler.checkCompatibility(ds3.compiler);
187         compiler.checkCompatibility(ds4.compiler);
188         compiler.linearCombination(a1, ds1.data, 0, a2, ds2.data, 0,
189                                    a3, ds3.data, 0, a4, ds4.data, 0,
190                                    data, 0);
191     }
192 
193     /** Build an instance from all its derivatives.
194      * @param parameters number of free parameters
195      * @param order derivation order
196      * @param derivatives derivatives sorted according to
197      * {@link DSCompiler#getPartialDerivativeIndex(int...)}
198      * @exception DimensionMismatchException if derivatives array does not match the
199      * {@link DSCompiler#getSize() size} expected by the compiler
200      * @throws NumberIsTooLargeException if order is too large
201      * @see #getAllDerivatives()
202      */
203     public DerivativeStructure(final int parameters, final int order, final double ... derivatives)
204         throws DimensionMismatchException, NumberIsTooLargeException {
205         this(parameters, order);
206         if (derivatives.length != data.length) {
207             throw new DimensionMismatchException(derivatives.length, data.length);
208         }
209         System.arraycopy(derivatives, 0, data, 0, data.length);
210     }
211 
212     /** Copy constructor.
213      * @param ds instance to copy
214      */
215     private DerivativeStructure(final DerivativeStructure ds) {
216         this.compiler = ds.compiler;
217         this.data     = ds.data.clone();
218     }
219 
220     /** Get the number of free parameters.
221      * @return number of free parameters
222      */
223     public int getFreeParameters() {
224         return compiler.getFreeParameters();
225     }
226 
227     /** Get the derivation order.
228      * @return derivation order
229      */
230     public int getOrder() {
231         return compiler.getOrder();
232     }
233 
234     /** Create a constant compatible with instance order and number of parameters.
235      * <p>
236      * This method is a convenience factory method, it simply calls
237      * {@code new DerivativeStructure(getFreeParameters(), getOrder(), c)}
238      * </p>
239      * @param c value of the constant
240      * @return a constant compatible with instance order and number of parameters
241      * @see #DerivativeStructure(int, int, double)
242      * @since 3.3
243      */
244     public DerivativeStructure createConstant(final double c) {
245         return new DerivativeStructure(getFreeParameters(), getOrder(), c);
246     }
247 
248     /** {@inheritDoc}
249      * @since 3.2
250      */
251     public double getReal() {
252         return data[0];
253     }
254 
255     /** Get the value part of the derivative structure.
256      * @return value part of the derivative structure
257      * @see #getPartialDerivative(int...)
258      */
259     public double getValue() {
260         return data[0];
261     }
262 
263     /** Get a partial derivative.
264      * @param orders derivation orders with respect to each variable (if all orders are 0,
265      * the value is returned)
266      * @return partial derivative
267      * @see #getValue()
268      * @exception DimensionMismatchException if the numbers of variables does not
269      * match the instance
270      * @exception NumberIsTooLargeException if sum of derivation orders is larger
271      * than the instance limits
272      */
273     public double getPartialDerivative(final int ... orders)
274         throws DimensionMismatchException, NumberIsTooLargeException {
275         return data[compiler.getPartialDerivativeIndex(orders)];
276     }
277 
278     /** Get all partial derivatives.
279      * @return a fresh copy of partial derivatives, in an array sorted according to
280      * {@link DSCompiler#getPartialDerivativeIndex(int...)}
281      */
282     public double[] getAllDerivatives() {
283         return data.clone();
284     }
285 
286     /** {@inheritDoc}
287      * @since 3.2
288      */
289     public DerivativeStructure add(final double a) {
290         final DerivativeStructure ds = new DerivativeStructure(this);
291         ds.data[0] += a;
292         return ds;
293     }
294 
295     /** {@inheritDoc}
296      * @exception DimensionMismatchException if number of free parameters
297      * or orders do not match
298      */
299     public DerivativeStructure add(final DerivativeStructure a)
300         throws DimensionMismatchException {
301         compiler.checkCompatibility(a.compiler);
302         final DerivativeStructure ds = new DerivativeStructure(this);
303         compiler.add(data, 0, a.data, 0, ds.data, 0);
304         return ds;
305     }
306 
307     /** {@inheritDoc}
308      * @since 3.2
309      */
310     public DerivativeStructure subtract(final double a) {
311         return add(-a);
312     }
313 
314     /** {@inheritDoc}
315      * @exception DimensionMismatchException if number of free parameters
316      * or orders do not match
317      */
318     public DerivativeStructure subtract(final DerivativeStructure a)
319         throws DimensionMismatchException {
320         compiler.checkCompatibility(a.compiler);
321         final DerivativeStructure ds = new DerivativeStructure(this);
322         compiler.subtract(data, 0, a.data, 0, ds.data, 0);
323         return ds;
324     }
325 
326     /** {@inheritDoc} */
327     public DerivativeStructure multiply(final int n) {
328         return multiply((double) n);
329     }
330 
331     /** {@inheritDoc}
332      * @since 3.2
333      */
334     public DerivativeStructure multiply(final double a) {
335         final DerivativeStructure ds = new DerivativeStructure(this);
336         for (int i = 0; i < ds.data.length; ++i) {
337             ds.data[i] *= a;
338         }
339         return ds;
340     }
341 
342     /** {@inheritDoc}
343      * @exception DimensionMismatchException if number of free parameters
344      * or orders do not match
345      */
346     public DerivativeStructure multiply(final DerivativeStructure a)
347         throws DimensionMismatchException {
348         compiler.checkCompatibility(a.compiler);
349         final DerivativeStructure result = new DerivativeStructure(compiler);
350         compiler.multiply(data, 0, a.data, 0, result.data, 0);
351         return result;
352     }
353 
354     /** {@inheritDoc}
355      * @since 3.2
356      */
357     public DerivativeStructure divide(final double a) {
358         final DerivativeStructure ds = new DerivativeStructure(this);
359         for (int i = 0; i < ds.data.length; ++i) {
360             ds.data[i] /= a;
361         }
362         return ds;
363     }
364 
365     /** {@inheritDoc}
366      * @exception DimensionMismatchException if number of free parameters
367      * or orders do not match
368      */
369     public DerivativeStructure divide(final DerivativeStructure a)
370         throws DimensionMismatchException {
371         compiler.checkCompatibility(a.compiler);
372         final DerivativeStructure result = new DerivativeStructure(compiler);
373         compiler.divide(data, 0, a.data, 0, result.data, 0);
374         return result;
375     }
376 
377     /** {@inheritDoc} */
378     public DerivativeStructure remainder(final double a) {
379         final DerivativeStructure ds = new DerivativeStructure(this);
380         ds.data[0] = FastMath.IEEEremainder(ds.data[0], a);
381         return ds;
382     }
383 
384     /** {@inheritDoc}
385      * @exception DimensionMismatchException if number of free parameters
386      * or orders do not match
387      * @since 3.2
388      */
389     public DerivativeStructure remainder(final DerivativeStructure a)
390         throws DimensionMismatchException {
391         compiler.checkCompatibility(a.compiler);
392         final DerivativeStructure result = new DerivativeStructure(compiler);
393         compiler.remainder(data, 0, a.data, 0, result.data, 0);
394         return result;
395     }
396 
397     /** {@inheritDoc} */
398     public DerivativeStructure negate() {
399         final DerivativeStructure ds = new DerivativeStructure(compiler);
400         for (int i = 0; i < ds.data.length; ++i) {
401             ds.data[i] = -data[i];
402         }
403         return ds;
404     }
405 
406     /** {@inheritDoc}
407      * @since 3.2
408      */
409     public DerivativeStructure abs() {
410         if (Double.doubleToLongBits(data[0]) < 0) {
411             // we use the bits representation to also handle -0.0
412             return negate();
413         } else {
414             return this;
415         }
416     }
417 
418     /** {@inheritDoc}
419      * @since 3.2
420      */
421     public DerivativeStructure ceil() {
422         return new DerivativeStructure(compiler.getFreeParameters(),
423                                        compiler.getOrder(),
424                                        FastMath.ceil(data[0]));
425     }
426 
427     /** {@inheritDoc}
428      * @since 3.2
429      */
430     public DerivativeStructure floor() {
431         return new DerivativeStructure(compiler.getFreeParameters(),
432                                        compiler.getOrder(),
433                                        FastMath.floor(data[0]));
434     }
435 
436     /** {@inheritDoc}
437      * @since 3.2
438      */
439     public DerivativeStructure rint() {
440         return new DerivativeStructure(compiler.getFreeParameters(),
441                                        compiler.getOrder(),
442                                        FastMath.rint(data[0]));
443     }
444 
445     /** {@inheritDoc} */
446     public long round() {
447         return FastMath.round(data[0]);
448     }
449 
450     /** {@inheritDoc}
451      * @since 3.2
452      */
453     public DerivativeStructure signum() {
454         return new DerivativeStructure(compiler.getFreeParameters(),
455                                        compiler.getOrder(),
456                                        FastMath.signum(data[0]));
457     }
458 
459     /** {@inheritDoc}
460      * @since 3.2
461      */
462     public DerivativeStructure copySign(final DerivativeStructure sign){
463         long m = Double.doubleToLongBits(data[0]);
464         long s = Double.doubleToLongBits(sign.data[0]);
465         if ((m >= 0 && s >= 0) || (m < 0 && s < 0)) { // Sign is currently OK
466             return this;
467         }
468         return negate(); // flip sign
469     }
470 
471     /** {@inheritDoc}
472      * @since 3.2
473      */
474     public DerivativeStructure copySign(final double sign) {
475         long m = Double.doubleToLongBits(data[0]);
476         long s = Double.doubleToLongBits(sign);
477         if ((m >= 0 && s >= 0) || (m < 0 && s < 0)) { // Sign is currently OK
478             return this;
479         }
480         return negate(); // flip sign
481     }
482 
483     /**
484      * Return the exponent of the instance value, removing the bias.
485      * <p>
486      * For double numbers of the form 2<sup>x</sup>, the unbiased
487      * exponent is exactly x.
488      * </p>
489      * @return exponent for instance in IEEE754 representation, without bias
490      */
491     public int getExponent() {
492         return FastMath.getExponent(data[0]);
493     }
494 
495     /** {@inheritDoc}
496      * @since 3.2
497      */
498     public DerivativeStructure scalb(final int n) {
499         final DerivativeStructure ds = new DerivativeStructure(compiler);
500         for (int i = 0; i < ds.data.length; ++i) {
501             ds.data[i] = FastMath.scalb(data[i], n);
502         }
503         return ds;
504     }
505 
506     /** {@inheritDoc}
507      * @exception DimensionMismatchException if number of free parameters
508      * or orders do not match
509      * @since 3.2
510      */
511     public DerivativeStructure hypot(final DerivativeStructure y)
512         throws DimensionMismatchException {
513 
514         compiler.checkCompatibility(y.compiler);
515 
516         if (Double.isInfinite(data[0]) || Double.isInfinite(y.data[0])) {
517             return new DerivativeStructure(compiler.getFreeParameters(),
518                                            compiler.getFreeParameters(),
519                                            Double.POSITIVE_INFINITY);
520         } else if (Double.isNaN(data[0]) || Double.isNaN(y.data[0])) {
521             return new DerivativeStructure(compiler.getFreeParameters(),
522                                            compiler.getFreeParameters(),
523                                            Double.NaN);
524         } else {
525 
526             final int expX = getExponent();
527             final int expY = y.getExponent();
528             if (expX > expY + 27) {
529                 // y is neglectible with respect to x
530                 return abs();
531             } else if (expY > expX + 27) {
532                 // x is neglectible with respect to y
533                 return y.abs();
534             } else {
535 
536                 // find an intermediate scale to avoid both overflow and underflow
537                 final int middleExp = (expX + expY) / 2;
538 
539                 // scale parameters without losing precision
540                 final DerivativeStructure scaledX = scalb(-middleExp);
541                 final DerivativeStructure scaledY = y.scalb(-middleExp);
542 
543                 // compute scaled hypotenuse
544                 final DerivativeStructure scaledH =
545                         scaledX.multiply(scaledX).add(scaledY.multiply(scaledY)).sqrt();
546 
547                 // remove scaling
548                 return scaledH.scalb(middleExp);
549 
550             }
551 
552         }
553     }
554 
555     /**
556      * Returns the hypotenuse of a triangle with sides {@code x} and {@code y}
557      * - sqrt(<i>x</i><sup>2</sup>&nbsp;+<i>y</i><sup>2</sup>)<br/>
558      * avoiding intermediate overflow or underflow.
559      *
560      * <ul>
561      * <li> If either argument is infinite, then the result is positive infinity.</li>
562      * <li> else, if either argument is NaN then the result is NaN.</li>
563      * </ul>
564      *
565      * @param x a value
566      * @param y a value
567      * @return sqrt(<i>x</i><sup>2</sup>&nbsp;+<i>y</i><sup>2</sup>)
568      * @exception DimensionMismatchException if number of free parameters
569      * or orders do not match
570      * @since 3.2
571      */
572     public static DerivativeStructure hypot(final DerivativeStructure x, final DerivativeStructure y)
573         throws DimensionMismatchException {
574         return x.hypot(y);
575     }
576 
577     /** Compute composition of the instance by a univariate function.
578      * @param f array of value and derivatives of the function at
579      * the current point (i.e. [f({@link #getValue()}),
580      * f'({@link #getValue()}), f''({@link #getValue()})...]).
581      * @return f(this)
582      * @exception DimensionMismatchException if the number of derivatives
583      * in the array is not equal to {@link #getOrder() order} + 1
584      */
585     public DerivativeStructure compose(final double ... f)
586         throws DimensionMismatchException {
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     }
594 
595     /** {@inheritDoc} */
596     public DerivativeStructure reciprocal() {
597         final DerivativeStructure result = new DerivativeStructure(compiler);
598         compiler.pow(data, 0, -1, result.data, 0);
599         return result;
600     }
601 
602     /** {@inheritDoc}
603      * @since 3.2
604      */
605     public DerivativeStructure sqrt() {
606         return rootN(2);
607     }
608 
609     /** {@inheritDoc}
610      * @since 3.2
611      */
612     public DerivativeStructure cbrt() {
613         return rootN(3);
614     }
615 
616     /** {@inheritDoc}
617      * @since 3.2
618      */
619     public DerivativeStructure rootN(final int n) {
620         final DerivativeStructure result = new DerivativeStructure(compiler);
621         compiler.rootN(data, 0, n, result.data, 0);
622         return result;
623     }
624 
625     /** {@inheritDoc} */
626     public Field<DerivativeStructure> getField() {
627         return new Field<DerivativeStructure>() {
628 
629             /** {@inheritDoc} */
630             public DerivativeStructure getZero() {
631                 return new DerivativeStructure(compiler.getFreeParameters(), compiler.getOrder(), 0.0);
632             }
633 
634             /** {@inheritDoc} */
635             public DerivativeStructure getOne() {
636                 return new DerivativeStructure(compiler.getFreeParameters(), compiler.getOrder(), 1.0);
637             }
638 
639             /** {@inheritDoc} */
640             public Class<? extends FieldElement<DerivativeStructure>> getRuntimeClass() {
641                 return DerivativeStructure.class;
642             }
643 
644         };
645     }
646 
647     /** Compute a<sup>x</sup> where a is a double and x a {@link DerivativeStructure}
648      * @param a number to exponentiate
649      * @param x power to apply
650      * @return a<sup>x</sup>
651      * @since 3.3
652      */
653     public static DerivativeStructure pow(final double a, final DerivativeStructure x) {
654         final DerivativeStructure result = new DerivativeStructure(x.compiler);
655         x.compiler.pow(a, x.data, 0, result.data, 0);
656         return result;
657     }
658 
659     /** {@inheritDoc}
660      * @since 3.2
661      */
662     public DerivativeStructure pow(final double p) {
663         final DerivativeStructure result = new DerivativeStructure(compiler);
664         compiler.pow(data, 0, p, result.data, 0);
665         return result;
666     }
667 
668     /** {@inheritDoc}
669      * @since 3.2
670      */
671     public DerivativeStructure pow(final int n) {
672         final DerivativeStructure result = new DerivativeStructure(compiler);
673         compiler.pow(data, 0, n, result.data, 0);
674         return result;
675     }
676 
677     /** {@inheritDoc}
678      * @exception DimensionMismatchException if number of free parameters
679      * or orders do not match
680      * @since 3.2
681      */
682     public DerivativeStructure pow(final DerivativeStructure e)
683         throws DimensionMismatchException {
684         compiler.checkCompatibility(e.compiler);
685         final DerivativeStructure result = new DerivativeStructure(compiler);
686         compiler.pow(data, 0, e.data, 0, result.data, 0);
687         return result;
688     }
689 
690     /** {@inheritDoc}
691      * @since 3.2
692      */
693     public DerivativeStructure exp() {
694         final DerivativeStructure result = new DerivativeStructure(compiler);
695         compiler.exp(data, 0, result.data, 0);
696         return result;
697     }
698 
699     /** {@inheritDoc}
700      * @since 3.2
701      */
702     public DerivativeStructure expm1() {
703         final DerivativeStructure result = new DerivativeStructure(compiler);
704         compiler.expm1(data, 0, result.data, 0);
705         return result;
706     }
707 
708     /** {@inheritDoc}
709      * @since 3.2
710      */
711     public DerivativeStructure log() {
712         final DerivativeStructure result = new DerivativeStructure(compiler);
713         compiler.log(data, 0, result.data, 0);
714         return result;
715     }
716 
717     /** {@inheritDoc}
718      * @since 3.2
719      */
720     public DerivativeStructure log1p() {
721         final DerivativeStructure result = new DerivativeStructure(compiler);
722         compiler.log1p(data, 0, result.data, 0);
723         return result;
724     }
725 
726     /** Base 10 logarithm.
727      * @return base 10 logarithm of the instance
728      */
729     public DerivativeStructure log10() {
730         final DerivativeStructure result = new DerivativeStructure(compiler);
731         compiler.log10(data, 0, result.data, 0);
732         return result;
733     }
734 
735     /** {@inheritDoc}
736      * @since 3.2
737      */
738     public DerivativeStructure cos() {
739         final DerivativeStructure result = new DerivativeStructure(compiler);
740         compiler.cos(data, 0, result.data, 0);
741         return result;
742     }
743 
744     /** {@inheritDoc}
745      * @since 3.2
746      */
747     public DerivativeStructure sin() {
748         final DerivativeStructure result = new DerivativeStructure(compiler);
749         compiler.sin(data, 0, result.data, 0);
750         return result;
751     }
752 
753     /** {@inheritDoc}
754      * @since 3.2
755      */
756     public DerivativeStructure tan() {
757         final DerivativeStructure result = new DerivativeStructure(compiler);
758         compiler.tan(data, 0, result.data, 0);
759         return result;
760     }
761 
762     /** {@inheritDoc}
763      * @since 3.2
764      */
765     public DerivativeStructure acos() {
766         final DerivativeStructure result = new DerivativeStructure(compiler);
767         compiler.acos(data, 0, result.data, 0);
768         return result;
769     }
770 
771     /** {@inheritDoc}
772      * @since 3.2
773      */
774     public DerivativeStructure asin() {
775         final DerivativeStructure result = new DerivativeStructure(compiler);
776         compiler.asin(data, 0, result.data, 0);
777         return result;
778     }
779 
780     /** {@inheritDoc}
781      * @since 3.2
782      */
783     public DerivativeStructure atan() {
784         final DerivativeStructure result = new DerivativeStructure(compiler);
785         compiler.atan(data, 0, result.data, 0);
786         return result;
787     }
788 
789     /** {@inheritDoc}
790      * @since 3.2
791      */
792     public DerivativeStructure atan2(final DerivativeStructure x)
793         throws DimensionMismatchException {
794         compiler.checkCompatibility(x.compiler);
795         final DerivativeStructure result = new DerivativeStructure(compiler);
796         compiler.atan2(data, 0, x.data, 0, result.data, 0);
797         return result;
798     }
799 
800     /** Two arguments arc tangent operation.
801      * @param y first argument of the arc tangent
802      * @param x second argument of the arc tangent
803      * @return atan2(y, x)
804      * @exception DimensionMismatchException if number of free parameters
805      * or orders do not match
806      * @since 3.2
807      */
808     public static DerivativeStructure atan2(final DerivativeStructure y, final DerivativeStructure x)
809         throws DimensionMismatchException {
810         return y.atan2(x);
811     }
812 
813     /** {@inheritDoc}
814      * @since 3.2
815      */
816     public DerivativeStructure cosh() {
817         final DerivativeStructure result = new DerivativeStructure(compiler);
818         compiler.cosh(data, 0, result.data, 0);
819         return result;
820     }
821 
822     /** {@inheritDoc}
823      * @since 3.2
824      */
825     public DerivativeStructure sinh() {
826         final DerivativeStructure result = new DerivativeStructure(compiler);
827         compiler.sinh(data, 0, result.data, 0);
828         return result;
829     }
830 
831     /** {@inheritDoc}
832      * @since 3.2
833      */
834     public DerivativeStructure tanh() {
835         final DerivativeStructure result = new DerivativeStructure(compiler);
836         compiler.tanh(data, 0, result.data, 0);
837         return result;
838     }
839 
840     /** {@inheritDoc}
841      * @since 3.2
842      */
843     public DerivativeStructure acosh() {
844         final DerivativeStructure result = new DerivativeStructure(compiler);
845         compiler.acosh(data, 0, result.data, 0);
846         return result;
847     }
848 
849     /** {@inheritDoc}
850      * @since 3.2
851      */
852     public DerivativeStructure asinh() {
853         final DerivativeStructure result = new DerivativeStructure(compiler);
854         compiler.asinh(data, 0, result.data, 0);
855         return result;
856     }
857 
858     /** {@inheritDoc}
859      * @since 3.2
860      */
861     public DerivativeStructure atanh() {
862         final DerivativeStructure result = new DerivativeStructure(compiler);
863         compiler.atanh(data, 0, result.data, 0);
864         return result;
865     }
866 
867     /** Convert radians to degrees, with error of less than 0.5 ULP
868      *  @return instance converted into degrees
869      */
870     public DerivativeStructure toDegrees() {
871         final DerivativeStructure ds = new DerivativeStructure(compiler);
872         for (int i = 0; i < ds.data.length; ++i) {
873             ds.data[i] = FastMath.toDegrees(data[i]);
874         }
875         return ds;
876     }
877 
878     /** Convert degrees to radians, with error of less than 0.5 ULP
879      *  @return instance converted into radians
880      */
881     public DerivativeStructure toRadians() {
882         final DerivativeStructure ds = new DerivativeStructure(compiler);
883         for (int i = 0; i < ds.data.length; ++i) {
884             ds.data[i] = FastMath.toRadians(data[i]);
885         }
886         return ds;
887     }
888 
889     /** Evaluate Taylor expansion a derivative structure.
890      * @param delta parameters offsets (&Delta;x, &Delta;y, ...)
891      * @return value of the Taylor expansion at x + &Delta;x, y + &Delta;y, ...
892      * @throws MathArithmeticException if factorials becomes too large
893      */
894     public double taylor(final double ... delta) throws MathArithmeticException {
895         return compiler.taylor(data, 0, delta);
896     }
897 
898     /** {@inheritDoc}
899      * @exception DimensionMismatchException if number of free parameters
900      * or orders do not match
901      * @since 3.2
902      */
903     public DerivativeStructure linearCombination(final DerivativeStructure[] a, final DerivativeStructure[] b)
904         throws DimensionMismatchException {
905 
906         // compute an accurate value, taking care of cancellations
907         final double[] aDouble = new double[a.length];
908         for (int i = 0; i < a.length; ++i) {
909             aDouble[i] = a[i].getValue();
910         }
911         final double[] bDouble = new double[b.length];
912         for (int i = 0; i < b.length; ++i) {
913             bDouble[i] = b[i].getValue();
914         }
915         final double accurateValue = MathArrays.linearCombination(aDouble, bDouble);
916 
917         // compute a simple value, with all partial derivatives
918         DerivativeStructure simpleValue = a[0].getField().getZero();
919         for (int i = 0; i < a.length; ++i) {
920             simpleValue = simpleValue.add(a[i].multiply(b[i]));
921         }
922 
923         // create a result with accurate value and all derivatives (not necessarily as accurate as the value)
924         final double[] all = simpleValue.getAllDerivatives();
925         all[0] = accurateValue;
926         return new DerivativeStructure(simpleValue.getFreeParameters(), simpleValue.getOrder(), all);
927 
928     }
929 
930     /** {@inheritDoc}
931      * @exception DimensionMismatchException if number of free parameters
932      * or orders do not match
933      * @since 3.2
934      */
935     public DerivativeStructure linearCombination(final double[] a, final DerivativeStructure[] b)
936         throws DimensionMismatchException {
937 
938         // compute an accurate value, taking care of cancellations
939         final double[] bDouble = new double[b.length];
940         for (int i = 0; i < b.length; ++i) {
941             bDouble[i] = b[i].getValue();
942         }
943         final double accurateValue = MathArrays.linearCombination(a, bDouble);
944 
945         // compute a simple value, with all partial derivatives
946         DerivativeStructure simpleValue = b[0].getField().getZero();
947         for (int i = 0; i < a.length; ++i) {
948             simpleValue = simpleValue.add(b[i].multiply(a[i]));
949         }
950 
951         // create a result with accurate value and all derivatives (not necessarily as accurate as the value)
952         final double[] all = simpleValue.getAllDerivatives();
953         all[0] = accurateValue;
954         return new DerivativeStructure(simpleValue.getFreeParameters(), simpleValue.getOrder(), all);
955 
956     }
957 
958     /** {@inheritDoc}
959      * @exception DimensionMismatchException if number of free parameters
960      * or orders do not match
961      * @since 3.2
962      */
963     public DerivativeStructure linearCombination(final DerivativeStructure a1, final DerivativeStructure b1,
964                                                  final DerivativeStructure a2, final DerivativeStructure b2)
965         throws DimensionMismatchException {
966 
967         // compute an accurate value, taking care of cancellations
968         final double accurateValue = MathArrays.linearCombination(a1.getValue(), b1.getValue(),
969                                                                   a2.getValue(), b2.getValue());
970 
971         // compute a simple value, with all partial derivatives
972         final DerivativeStructure simpleValue = a1.multiply(b1).add(a2.multiply(b2));
973 
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(getFreeParameters(), getOrder(), all);
978 
979     }
980 
981     /** {@inheritDoc}
982      * @exception DimensionMismatchException if number of free parameters
983      * or orders do not match
984      * @since 3.2
985      */
986     public DerivativeStructure linearCombination(final double a1, final DerivativeStructure b1,
987                                                  final double a2, final DerivativeStructure b2)
988         throws DimensionMismatchException {
989 
990         // compute an accurate value, taking care of cancellations
991         final double accurateValue = MathArrays.linearCombination(a1, b1.getValue(),
992                                                                   a2, b2.getValue());
993 
994         // compute a simple value, with all partial derivatives
995         final DerivativeStructure simpleValue = b1.multiply(a1).add(b2.multiply(a2));
996 
997         // create a result with accurate value and all derivatives (not necessarily as accurate as the value)
998         final double[] all = simpleValue.getAllDerivatives();
999         all[0] = accurateValue;
1000         return new DerivativeStructure(getFreeParameters(), getOrder(), all);
1001 
1002     }
1003 
1004     /** {@inheritDoc}
1005      * @exception DimensionMismatchException if number of free parameters
1006      * or orders do not match
1007      * @since 3.2
1008      */
1009     public DerivativeStructure linearCombination(final DerivativeStructure a1, final DerivativeStructure b1,
1010                                                  final DerivativeStructure a2, final DerivativeStructure b2,
1011                                                  final DerivativeStructure a3, final DerivativeStructure b3)
1012         throws DimensionMismatchException {
1013 
1014         // compute an accurate value, taking care of cancellations
1015         final double accurateValue = MathArrays.linearCombination(a1.getValue(), b1.getValue(),
1016                                                                   a2.getValue(), b2.getValue(),
1017                                                                   a3.getValue(), b3.getValue());
1018 
1019         // compute a simple value, with all partial derivatives
1020         final DerivativeStructure simpleValue = a1.multiply(b1).add(a2.multiply(b2)).add(a3.multiply(b3));
1021 
1022         // create a result with accurate value and all derivatives (not necessarily as accurate as the value)
1023         final double[] all = simpleValue.getAllDerivatives();
1024         all[0] = accurateValue;
1025         return new DerivativeStructure(getFreeParameters(), getOrder(), all);
1026 
1027     }
1028 
1029     /** {@inheritDoc}
1030      * @exception DimensionMismatchException if number of free parameters
1031      * or orders do not match
1032      * @since 3.2
1033      */
1034     public DerivativeStructure linearCombination(final double a1, final DerivativeStructure b1,
1035                                                  final double a2, final DerivativeStructure b2,
1036                                                  final double a3, final DerivativeStructure b3)
1037         throws DimensionMismatchException {
1038 
1039         // compute an accurate value, taking care of cancellations
1040         final double accurateValue = MathArrays.linearCombination(a1, b1.getValue(),
1041                                                                   a2, b2.getValue(),
1042                                                                   a3, b3.getValue());
1043 
1044         // compute a simple value, with all partial derivatives
1045         final DerivativeStructure simpleValue = b1.multiply(a1).add(b2.multiply(a2)).add(b3.multiply(a3));
1046 
1047         // create a result with accurate value and all derivatives (not necessarily as accurate as the value)
1048         final double[] all = simpleValue.getAllDerivatives();
1049         all[0] = accurateValue;
1050         return new DerivativeStructure(getFreeParameters(), getOrder(), all);
1051 
1052     }
1053 
1054     /** {@inheritDoc}
1055      * @exception DimensionMismatchException if number of free parameters
1056      * or orders do not match
1057      * @since 3.2
1058      */
1059     public DerivativeStructure linearCombination(final DerivativeStructure a1, final DerivativeStructure b1,
1060                                                  final DerivativeStructure a2, final DerivativeStructure b2,
1061                                                  final DerivativeStructure a3, final DerivativeStructure b3,
1062                                                  final DerivativeStructure a4, final DerivativeStructure b4)
1063         throws DimensionMismatchException {
1064 
1065         // compute an accurate value, taking care of cancellations
1066         final double accurateValue = MathArrays.linearCombination(a1.getValue(), b1.getValue(),
1067                                                                   a2.getValue(), b2.getValue(),
1068                                                                   a3.getValue(), b3.getValue(),
1069                                                                   a4.getValue(), b4.getValue());
1070 
1071         // compute a simple value, with all partial derivatives
1072         final DerivativeStructure simpleValue = a1.multiply(b1).add(a2.multiply(b2)).add(a3.multiply(b3)).add(a4.multiply(b4));
1073 
1074         // create a result with accurate value and all derivatives (not necessarily as accurate as the value)
1075         final double[] all = simpleValue.getAllDerivatives();
1076         all[0] = accurateValue;
1077         return new DerivativeStructure(getFreeParameters(), getOrder(), all);
1078 
1079     }
1080 
1081     /** {@inheritDoc}
1082      * @exception DimensionMismatchException if number of free parameters
1083      * or orders do not match
1084      * @since 3.2
1085      */
1086     public DerivativeStructure linearCombination(final double a1, final DerivativeStructure b1,
1087                                                  final double a2, final DerivativeStructure b2,
1088                                                  final double a3, final DerivativeStructure b3,
1089                                                  final double a4, final DerivativeStructure b4)
1090         throws DimensionMismatchException {
1091 
1092         // compute an accurate value, taking care of cancellations
1093         final double accurateValue = MathArrays.linearCombination(a1, b1.getValue(),
1094                                                                   a2, b2.getValue(),
1095                                                                   a3, b3.getValue(),
1096                                                                   a4, b4.getValue());
1097 
1098         // compute a simple value, with all partial derivatives
1099         final DerivativeStructure simpleValue = b1.multiply(a1).add(b2.multiply(a2)).add(b3.multiply(a3)).add(b4.multiply(a4));
1100 
1101         // create a result with accurate value and all derivatives (not necessarily as accurate as the value)
1102         final double[] all = simpleValue.getAllDerivatives();
1103         all[0] = accurateValue;
1104         return new DerivativeStructure(getFreeParameters(), getOrder(), all);
1105 
1106     }
1107 
1108     /**
1109      * Test for the equality of two derivative structures.
1110      * <p>
1111      * Derivative structures are considered equal if they have the same number
1112      * of free parameters, the same derivation order, and the same derivatives.
1113      * </p>
1114      * @param other Object to test for equality to this
1115      * @return true if two derivative structures are equal
1116      * @since 3.2
1117      */
1118     @Override
1119     public boolean equals(Object other) {
1120 
1121         if (this == other) {
1122             return true;
1123         }
1124 
1125         if (other instanceof DerivativeStructure) {
1126             final DerivativeStructure rhs = (DerivativeStructure)other;
1127             return (getFreeParameters() == rhs.getFreeParameters()) &&
1128                    (getOrder() == rhs.getOrder()) &&
1129                    MathArrays.equals(data, rhs.data);
1130         }
1131 
1132         return false;
1133 
1134     }
1135 
1136     /**
1137      * Get a hashCode for the derivative structure.
1138      * @return a hash code value for this object
1139      * @since 3.2
1140      */
1141     @Override
1142     public int hashCode() {
1143         return 227 + 229 * getFreeParameters() + 233 * getOrder() + 239 * MathUtils.hash(data);
1144     }
1145 
1146     /**
1147      * Replace the instance with a data transfer object for serialization.
1148      * @return data transfer object that will be serialized
1149      */
1150     private Object writeReplace() {
1151         return new DataTransferObject(compiler.getFreeParameters(), compiler.getOrder(), data);
1152     }
1153 
1154     /** Internal class used only for serialization. */
1155     private static class DataTransferObject implements Serializable {
1156 
1157         /** Serializable UID. */
1158         private static final long serialVersionUID = 20120730L;
1159 
1160         /** Number of variables.
1161          * @serial
1162          */
1163         private final int variables;
1164 
1165         /** Derivation order.
1166          * @serial
1167          */
1168         private final int order;
1169 
1170         /** Partial derivatives.
1171          * @serial
1172          */
1173         private final double[] data;
1174 
1175         /** Simple constructor.
1176          * @param variables number of variables
1177          * @param order derivation order
1178          * @param data partial derivatives
1179          */
1180         public DataTransferObject(final int variables, final int order, final double[] data) {
1181             this.variables = variables;
1182             this.order     = order;
1183             this.data      = data;
1184         }
1185 
1186         /** Replace the deserialized data transfer object with a {@link DerivativeStructure}.
1187          * @return replacement {@link DerivativeStructure}
1188          */
1189         private Object readResolve() {
1190             return new DerivativeStructure(variables, order, data);
1191         }
1192 
1193     }
1194 
1195 }