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