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.math4.legacy.analysis.differentiation;
18  
19  import java.util.Collections;
20  import java.util.HashMap;
21  import java.util.Map;
22  
23  import org.apache.commons.numbers.core.Sum;
24  import org.apache.commons.numbers.core.Precision;
25  import org.apache.commons.math4.legacy.core.Field;
26  import org.apache.commons.math4.legacy.core.FieldElement;
27  import org.apache.commons.math4.legacy.core.RealFieldElement;
28  import org.apache.commons.math4.core.jdkmath.JdkMath;
29  
30  /**
31   * First derivative computation with large number of variables.
32   * <p>
33   * This class plays a similar role to {@link DerivativeStructure}, with
34   * a focus on efficiency when dealing with large number of independent variables
35   * and most computation depend only on a few of them, and when only first derivative
36   * is desired. When these conditions are met, this class should be much faster than
37   * {@link DerivativeStructure} and use less memory.
38   * </p>
39   *
40   * @since 3.3
41   */
42  public final class SparseGradient implements RealFieldElement<SparseGradient> {
43      /** Value of the calculation. */
44      private double value;
45  
46      /** Stored derivative, each key representing a different independent variable. */
47      private final Map<Integer, Double> derivatives;
48  
49      /** Internal constructor.
50       * @param value value of the function
51       * @param derivatives derivatives map, a deep copy will be performed,
52       * so the map given here will remain safe from changes in the new instance,
53       * may be null to create an empty derivatives map, i.e. a constant value
54       */
55      private SparseGradient(final double value, final Map<Integer, Double> derivatives) {
56          this.value = value;
57          this.derivatives = new HashMap<>();
58          if (derivatives != null) {
59              this.derivatives.putAll(derivatives);
60          }
61      }
62  
63      /** Internal constructor.
64       * @param value value of the function
65       * @param scale scaling factor to apply to all derivatives
66       * @param derivatives derivatives map, a deep copy will be performed,
67       * so the map given here will remain safe from changes in the new instance,
68       * may be null to create an empty derivatives map, i.e. a constant value
69       */
70      private SparseGradient(final double value, final double scale,
71                               final Map<Integer, Double> derivatives) {
72          this.value = value;
73          this.derivatives = new HashMap<>();
74          if (derivatives != null) {
75              for (final Map.Entry<Integer, Double> entry : derivatives.entrySet()) {
76                  this.derivatives.put(entry.getKey(), scale * entry.getValue());
77              }
78          }
79      }
80  
81      /** Factory method creating a constant.
82       * @param value value of the constant
83       * @return a new instance
84       */
85      public static SparseGradient createConstant(final double value) {
86          return new SparseGradient(value, Collections.<Integer, Double>emptyMap());
87      }
88  
89      /** Factory method creating an independent variable.
90       * @param idx index of the variable
91       * @param value value of the variable
92       * @return a new instance
93       */
94      public static SparseGradient createVariable(final int idx, final double value) {
95          return new SparseGradient(value, Collections.singletonMap(idx, 1.0));
96      }
97  
98      /**
99       * Find the number of variables.
100      * @return number of variables
101      */
102     public int numVars() {
103         return derivatives.size();
104     }
105 
106     /**
107      * Get the derivative with respect to a particular index variable.
108      *
109      * @param index index to differentiate with.
110      * @return derivative with respect to a particular index variable
111      */
112     public double getDerivative(final int index) {
113         final Double out = derivatives.get(index);
114         return (out == null) ? 0.0 : out;
115     }
116 
117     /**
118      * Get the value of the function.
119      * @return value of the function.
120      */
121     public double getValue() {
122         return value;
123     }
124 
125     /** {@inheritDoc} */
126     @Override
127     public double getReal() {
128         return value;
129     }
130 
131     /** {@inheritDoc} */
132     @Override
133     public SparseGradient add(final SparseGradient a) {
134         final SparseGradient out = new SparseGradient(value + a.value, derivatives);
135         for (Map.Entry<Integer, Double> entry : a.derivatives.entrySet()) {
136             final int id = entry.getKey();
137             final Double old = out.derivatives.get(id);
138             if (old == null) {
139                 out.derivatives.put(id, entry.getValue());
140             } else {
141                 out.derivatives.put(id, old + entry.getValue());
142             }
143         }
144 
145         return out;
146     }
147 
148     /**
149      * Add in place.
150      * <p>
151      * This method is designed to be faster when used multiple times in a loop.
152      * </p>
153      * <p>
154      * The instance is changed here, in order to not change the
155      * instance the {@link #add(SparseGradient)} method should
156      * be used.
157      * </p>
158      * @param a instance to add
159      */
160     public void addInPlace(final SparseGradient a) {
161         value += a.value;
162         for (final Map.Entry<Integer, Double> entry : a.derivatives.entrySet()) {
163             final int id = entry.getKey();
164             final Double old = derivatives.get(id);
165             if (old == null) {
166                 derivatives.put(id, entry.getValue());
167             } else {
168                 derivatives.put(id, old + entry.getValue());
169             }
170         }
171     }
172 
173     /** {@inheritDoc} */
174     @Override
175     public SparseGradient add(final double c) {
176         return new SparseGradient(value + c, derivatives);
177     }
178 
179     /** {@inheritDoc} */
180     @Override
181     public SparseGradient subtract(final SparseGradient a) {
182         final SparseGradient out = new SparseGradient(value - a.value, derivatives);
183         for (Map.Entry<Integer, Double> entry : a.derivatives.entrySet()) {
184             final int id = entry.getKey();
185             final Double old = out.derivatives.get(id);
186             if (old == null) {
187                 out.derivatives.put(id, -entry.getValue());
188             } else {
189                 out.derivatives.put(id, old - entry.getValue());
190             }
191         }
192         return out;
193     }
194 
195     /** {@inheritDoc} */
196     @Override
197     public SparseGradient subtract(double c) {
198         return new SparseGradient(value - c, derivatives);
199     }
200 
201     /** {@inheritDoc} */
202     @Override
203     public SparseGradient multiply(final SparseGradient a) {
204         final SparseGradient out =
205             new SparseGradient(value * a.value, Collections.<Integer, Double>emptyMap());
206 
207         // Derivatives.
208         for (Map.Entry<Integer, Double> entry : derivatives.entrySet()) {
209             out.derivatives.put(entry.getKey(), a.value * entry.getValue());
210         }
211         for (Map.Entry<Integer, Double> entry : a.derivatives.entrySet()) {
212             final int id = entry.getKey();
213             final Double old = out.derivatives.get(id);
214             if (old == null) {
215                 out.derivatives.put(id, value * entry.getValue());
216             } else {
217                 out.derivatives.put(id, old + value * entry.getValue());
218             }
219         }
220         return out;
221     }
222 
223     /**
224      * Multiply in place.
225      * <p>
226      * This method is designed to be faster when used multiple times in a loop.
227      * </p>
228      * <p>
229      * The instance is changed here, in order to not change the
230      * instance the {@link #add(SparseGradient)} method should
231      * be used.
232      * </p>
233      * @param a instance to multiply
234      */
235     public void multiplyInPlace(final SparseGradient a) {
236         // Derivatives.
237         for (Map.Entry<Integer, Double> entry : derivatives.entrySet()) {
238             derivatives.put(entry.getKey(), a.value * entry.getValue());
239         }
240         for (Map.Entry<Integer, Double> entry : a.derivatives.entrySet()) {
241             final int id = entry.getKey();
242             final Double old = derivatives.get(id);
243             if (old == null) {
244                 derivatives.put(id, value * entry.getValue());
245             } else {
246                 derivatives.put(id, old + value * entry.getValue());
247             }
248         }
249         value *= a.value;
250     }
251 
252     /** {@inheritDoc} */
253     @Override
254     public SparseGradient multiply(final double c) {
255         return new SparseGradient(value * c, c, derivatives);
256     }
257 
258     /** {@inheritDoc} */
259     @Override
260     public SparseGradient multiply(final int n) {
261         return new SparseGradient(value * n, n, derivatives);
262     }
263 
264     /** {@inheritDoc} */
265     @Override
266     public SparseGradient divide(final SparseGradient a) {
267         final SparseGradient out = new SparseGradient(value / a.value, Collections.<Integer, Double>emptyMap());
268 
269         // Derivatives.
270         for (Map.Entry<Integer, Double> entry : derivatives.entrySet()) {
271             out.derivatives.put(entry.getKey(), entry.getValue() / a.value);
272         }
273         for (Map.Entry<Integer, Double> entry : a.derivatives.entrySet()) {
274             final int id = entry.getKey();
275             final Double old = out.derivatives.get(id);
276             if (old == null) {
277                 out.derivatives.put(id, -out.value / a.value * entry.getValue());
278             } else {
279                 out.derivatives.put(id, old - out.value / a.value * entry.getValue());
280             }
281         }
282         return out;
283     }
284 
285     /** {@inheritDoc} */
286     @Override
287     public SparseGradient divide(final double c) {
288         return new SparseGradient(value / c, 1.0 / c, derivatives);
289     }
290 
291     /** {@inheritDoc} */
292     @Override
293     public SparseGradient negate() {
294         return new SparseGradient(-value, -1.0, derivatives);
295     }
296 
297     /** {@inheritDoc} */
298     @Override
299     public Field<SparseGradient> getField() {
300         return new Field<SparseGradient>() {
301 
302             /** {@inheritDoc} */
303             @Override
304             public SparseGradient getZero() {
305                 return createConstant(0);
306             }
307 
308             /** {@inheritDoc} */
309             @Override
310             public SparseGradient getOne() {
311                 return createConstant(1);
312             }
313 
314             /** {@inheritDoc} */
315             @Override
316             public Class<? extends FieldElement<SparseGradient>> getRuntimeClass() {
317                 return SparseGradient.class;
318             }
319         };
320     }
321 
322     /** {@inheritDoc} */
323     @Override
324     public SparseGradient remainder(final double a) {
325         return new SparseGradient(JdkMath.IEEEremainder(value, a), derivatives);
326     }
327 
328     /** {@inheritDoc} */
329     @Override
330     public SparseGradient remainder(final SparseGradient a) {
331 
332         // compute k such that lhs % rhs = lhs - k rhs
333         final double rem = JdkMath.IEEEremainder(value, a.value);
334         final double k   = JdkMath.rint((value - rem) / a.value);
335 
336         return subtract(a.multiply(k));
337     }
338 
339     /** {@inheritDoc} */
340     @Override
341     public SparseGradient abs() {
342         if (Double.doubleToLongBits(value) < 0) {
343             // we use the bits representation to also handle -0.0
344             return negate();
345         } else {
346             return this;
347         }
348     }
349 
350     /** {@inheritDoc} */
351     @Override
352     public SparseGradient ceil() {
353         return createConstant(JdkMath.ceil(value));
354     }
355 
356     /** {@inheritDoc} */
357     @Override
358     public SparseGradient floor() {
359         return createConstant(JdkMath.floor(value));
360     }
361 
362     /** {@inheritDoc} */
363     @Override
364     public SparseGradient rint() {
365         return createConstant(JdkMath.rint(value));
366     }
367 
368     /** {@inheritDoc} */
369     @Override
370     public long round() {
371         return JdkMath.round(value);
372     }
373 
374     /** {@inheritDoc} */
375     @Override
376     public SparseGradient signum() {
377         return createConstant(JdkMath.signum(value));
378     }
379 
380     /** {@inheritDoc} */
381     @Override
382     public SparseGradient copySign(final SparseGradient sign) {
383         final long m = Double.doubleToLongBits(value);
384         final long s = Double.doubleToLongBits(sign.value);
385         if ((m >= 0 && s >= 0) || (m < 0 && s < 0)) { // Sign is currently OK
386             return this;
387         }
388         return negate(); // flip sign
389     }
390 
391     /** {@inheritDoc} */
392     @Override
393     public SparseGradient copySign(final double sign) {
394         final long m = Double.doubleToLongBits(value);
395         final long s = Double.doubleToLongBits(sign);
396         if ((m >= 0 && s >= 0) || (m < 0 && s < 0)) { // Sign is currently OK
397             return this;
398         }
399         return negate(); // flip sign
400     }
401 
402     /** {@inheritDoc} */
403     @Override
404     public SparseGradient scalb(final int n) {
405         final SparseGradient out = new SparseGradient(JdkMath.scalb(value, n), Collections.<Integer, Double>emptyMap());
406         for (Map.Entry<Integer, Double> entry : derivatives.entrySet()) {
407             out.derivatives.put(entry.getKey(), JdkMath.scalb(entry.getValue(), n));
408         }
409         return out;
410     }
411 
412     /** {@inheritDoc} */
413     @Override
414     public SparseGradient hypot(final SparseGradient y) {
415         if (Double.isInfinite(value) || Double.isInfinite(y.value)) {
416             return createConstant(Double.POSITIVE_INFINITY);
417         } else if (Double.isNaN(value) || Double.isNaN(y.value)) {
418             return createConstant(Double.NaN);
419         } else {
420 
421             final int expX = JdkMath.getExponent(value);
422             final int expY = JdkMath.getExponent(y.value);
423             if (expX > expY + 27) {
424                 // y is negligible with respect to x
425                 return abs();
426             } else if (expY > expX + 27) {
427                 // x is negligible with respect to y
428                 return y.abs();
429             } else {
430 
431                 // find an intermediate scale to avoid both overflow and underflow
432                 final int middleExp = (expX + expY) / 2;
433 
434                 // scale parameters without losing precision
435                 final SparseGradient scaledX = scalb(-middleExp);
436                 final SparseGradient scaledY = y.scalb(-middleExp);
437 
438                 // compute scaled hypotenuse
439                 final SparseGradient scaledH =
440                         scaledX.multiply(scaledX).add(scaledY.multiply(scaledY)).sqrt();
441 
442                 // remove scaling
443                 return scaledH.scalb(middleExp);
444             }
445         }
446     }
447 
448     /**
449      * Returns the hypotenuse of a triangle with sides {@code x} and {@code y}
450      * - sqrt(<i>x</i><sup>2</sup>&nbsp;+<i>y</i><sup>2</sup>)
451      * avoiding intermediate overflow or underflow.
452      *
453      * <ul>
454      * <li> If either argument is infinite, then the result is positive infinity.</li>
455      * <li> else, if either argument is NaN then the result is NaN.</li>
456      * </ul>
457      *
458      * @param x a value
459      * @param y a value
460      * @return sqrt(<i>x</i><sup>2</sup>&nbsp;+<i>y</i><sup>2</sup>)
461      */
462     public static SparseGradient hypot(final SparseGradient x, final SparseGradient y) {
463         return x.hypot(y);
464     }
465 
466     /** {@inheritDoc} */
467     @Override
468     public SparseGradient reciprocal() {
469         return new SparseGradient(1.0 / value, -1.0 / (value * value), derivatives);
470     }
471 
472     /** {@inheritDoc} */
473     @Override
474     public SparseGradient sqrt() {
475         final double sqrt = JdkMath.sqrt(value);
476         return new SparseGradient(sqrt, 0.5 / sqrt, derivatives);
477     }
478 
479     /** {@inheritDoc} */
480     @Override
481     public SparseGradient cbrt() {
482         final double cbrt = JdkMath.cbrt(value);
483         return new SparseGradient(cbrt, 1.0 / (3 * cbrt * cbrt), derivatives);
484     }
485 
486     /** {@inheritDoc} */
487     @Override
488     public SparseGradient rootN(final int n) {
489         if (n == 2) {
490             return sqrt();
491         } else if (n == 3) {
492             return cbrt();
493         } else {
494             final double root = JdkMath.pow(value, 1.0 / n);
495             return new SparseGradient(root, 1.0 / (n * JdkMath.pow(root, n - 1)), derivatives);
496         }
497     }
498 
499     /** {@inheritDoc} */
500     @Override
501     public SparseGradient pow(final double p) {
502         return new SparseGradient(JdkMath.pow(value,  p), p * JdkMath.pow(value,  p - 1), derivatives);
503     }
504 
505     /** {@inheritDoc} */
506     @Override
507     public SparseGradient pow(final int n) {
508         if (n == 0) {
509             return getField().getOne();
510         } else {
511             final double valueNm1 = JdkMath.pow(value,  n - 1);
512             return new SparseGradient(value * valueNm1, n * valueNm1, derivatives);
513         }
514     }
515 
516     /** {@inheritDoc} */
517     @Override
518     public SparseGradient pow(final SparseGradient e) {
519         return log().multiply(e).exp();
520     }
521 
522     /** Compute a<sup>x</sup> where a is a double and x a {@link SparseGradient}.
523      * @param a number to exponentiate
524      * @param x power to apply
525      * @return a<sup>x</sup>
526      */
527     public static SparseGradient pow(final double a, final SparseGradient x) {
528         if (a == 0) {
529             if (x.value == 0) {
530                 return x.compose(1.0, Double.NEGATIVE_INFINITY);
531             } else if (x.value < 0) {
532                 return x.compose(Double.NaN, Double.NaN);
533             } else {
534                 return x.getField().getZero();
535             }
536         } else {
537             final double ax = JdkMath.pow(a, x.value);
538             return new SparseGradient(ax, ax * JdkMath.log(a), x.derivatives);
539         }
540     }
541 
542     /** {@inheritDoc} */
543     @Override
544     public SparseGradient exp() {
545         final double e = JdkMath.exp(value);
546         return new SparseGradient(e, e, derivatives);
547     }
548 
549     /** {@inheritDoc} */
550     @Override
551     public SparseGradient expm1() {
552         return new SparseGradient(JdkMath.expm1(value), JdkMath.exp(value), derivatives);
553     }
554 
555     /** {@inheritDoc} */
556     @Override
557     public SparseGradient log() {
558         return new SparseGradient(JdkMath.log(value), 1.0 / value, derivatives);
559     }
560 
561     /** Base 10 logarithm.
562      * @return base 10 logarithm of the instance
563      */
564     @Override
565     public SparseGradient log10() {
566         return new SparseGradient(JdkMath.log10(value), 1.0 / (JdkMath.log(10.0) * value), derivatives);
567     }
568 
569     /** {@inheritDoc} */
570     @Override
571     public SparseGradient log1p() {
572         return new SparseGradient(JdkMath.log1p(value), 1.0 / (1.0 + value), derivatives);
573     }
574 
575     /** {@inheritDoc} */
576     @Override
577     public SparseGradient cos() {
578         return new SparseGradient(JdkMath.cos(value), -JdkMath.sin(value), derivatives);
579     }
580 
581     /** {@inheritDoc} */
582     @Override
583     public SparseGradient sin() {
584         return new SparseGradient(JdkMath.sin(value), JdkMath.cos(value), derivatives);
585     }
586 
587     /** {@inheritDoc} */
588     @Override
589     public SparseGradient tan() {
590         final double t = JdkMath.tan(value);
591         return new SparseGradient(t, 1 + t * t, derivatives);
592     }
593 
594     /** {@inheritDoc} */
595     @Override
596     public SparseGradient acos() {
597         return new SparseGradient(JdkMath.acos(value), -1.0 / JdkMath.sqrt(1 - value * value), derivatives);
598     }
599 
600     /** {@inheritDoc} */
601     @Override
602     public SparseGradient asin() {
603         return new SparseGradient(JdkMath.asin(value), 1.0 / JdkMath.sqrt(1 - value * value), derivatives);
604     }
605 
606     /** {@inheritDoc} */
607     @Override
608     public SparseGradient atan() {
609         return new SparseGradient(JdkMath.atan(value), 1.0 / (1 + value * value), derivatives);
610     }
611 
612     /** {@inheritDoc} */
613     @Override
614     public SparseGradient atan2(final SparseGradient x) {
615 
616         // compute r = sqrt(x^2+y^2)
617         final SparseGradient r = multiply(this).add(x.multiply(x)).sqrt();
618 
619         final SparseGradient a;
620         if (x.value >= 0) {
621 
622             // compute atan2(y, x) = 2 atan(y / (r + x))
623             a = divide(r.add(x)).atan().multiply(2);
624         } else {
625 
626             // compute atan2(y, x) = +/- pi - 2 atan(y / (r - x))
627             final SparseGradient tmp = divide(r.subtract(x)).atan().multiply(-2);
628             a = tmp.add(tmp.value <= 0 ? -JdkMath.PI : JdkMath.PI);
629         }
630 
631         // fix value to take special cases (+0/+0, +0/-0, -0/+0, -0/-0, +/-infinity) correctly
632         a.value = JdkMath.atan2(value, x.value);
633 
634         return a;
635     }
636 
637     /** Two arguments arc tangent operation.
638      * @param y first argument of the arc tangent
639      * @param x second argument of the arc tangent
640      * @return atan2(y, x)
641      */
642     public static SparseGradient atan2(final SparseGradient y, final SparseGradient x) {
643         return y.atan2(x);
644     }
645 
646     /** {@inheritDoc} */
647     @Override
648     public SparseGradient cosh() {
649         return new SparseGradient(JdkMath.cosh(value), JdkMath.sinh(value), derivatives);
650     }
651 
652     /** {@inheritDoc} */
653     @Override
654     public SparseGradient sinh() {
655         return new SparseGradient(JdkMath.sinh(value), JdkMath.cosh(value), derivatives);
656     }
657 
658     /** {@inheritDoc} */
659     @Override
660     public SparseGradient tanh() {
661         final double t = JdkMath.tanh(value);
662         return new SparseGradient(t, 1 - t * t, derivatives);
663     }
664 
665     /** {@inheritDoc} */
666     @Override
667     public SparseGradient acosh() {
668         return new SparseGradient(JdkMath.acosh(value), 1.0 / JdkMath.sqrt(value * value - 1.0), derivatives);
669     }
670 
671     /** {@inheritDoc} */
672     @Override
673     public SparseGradient asinh() {
674         return new SparseGradient(JdkMath.asinh(value), 1.0 / JdkMath.sqrt(value * value + 1.0), derivatives);
675     }
676 
677     /** {@inheritDoc} */
678     @Override
679     public SparseGradient atanh() {
680         return new SparseGradient(JdkMath.atanh(value), 1.0 / (1.0 - value * value), derivatives);
681     }
682 
683     /** Convert radians to degrees, with error of less than 0.5 ULP.
684      *  @return instance converted into degrees
685      */
686     public SparseGradient toDegrees() {
687         return new SparseGradient(JdkMath.toDegrees(value), JdkMath.toDegrees(1.0), derivatives);
688     }
689 
690     /** Convert degrees to radians, with error of less than 0.5 ULP.
691      *  @return instance converted into radians
692      */
693     public SparseGradient toRadians() {
694         return new SparseGradient(JdkMath.toRadians(value), JdkMath.toRadians(1.0), derivatives);
695     }
696 
697     /** Evaluate Taylor expansion of a sparse gradient.
698      * @param delta parameters offsets (&Delta;x, &Delta;y, ...)
699      * @return value of the Taylor expansion at x + &Delta;x, y + &Delta;y, ...
700      */
701     public double taylor(final double ... delta) {
702         double y = value;
703         for (int i = 0; i < delta.length; ++i) {
704             y += delta[i] * getDerivative(i);
705         }
706         return y;
707     }
708 
709     /** Compute composition of the instance by a univariate function.
710      * @param f0 value of the function at (i.e. f({@link #getValue()}))
711      * @param f1 first derivative of the function at
712      * the current point (i.e. f'({@link #getValue()}))
713      * @return f(this)
714     */
715     public SparseGradient compose(final double f0, final double f1) {
716         return new SparseGradient(f0, f1, derivatives);
717     }
718 
719     /** {@inheritDoc} */
720     @Override
721     public SparseGradient linearCombination(final SparseGradient[] a,
722                                             final SparseGradient[] b) {
723         // compute a simple value, with all partial derivatives
724         SparseGradient out = a[0].getField().getZero();
725         for (int i = 0; i < a.length; ++i) {
726             out = out.add(a[i].multiply(b[i]));
727         }
728 
729         // recompute an accurate value, taking care of cancellations
730         final double[] aDouble = new double[a.length];
731         for (int i = 0; i < a.length; ++i) {
732             aDouble[i] = a[i].getValue();
733         }
734         final double[] bDouble = new double[b.length];
735         for (int i = 0; i < b.length; ++i) {
736             bDouble[i] = b[i].getValue();
737         }
738         out.value = Sum.ofProducts(aDouble, bDouble).getAsDouble();
739 
740         return out;
741     }
742 
743     /** {@inheritDoc} */
744     @Override
745     public SparseGradient linearCombination(final double[] a, final SparseGradient[] b) {
746 
747         // compute a simple value, with all partial derivatives
748         SparseGradient out = b[0].getField().getZero();
749         for (int i = 0; i < a.length; ++i) {
750             out = out.add(b[i].multiply(a[i]));
751         }
752 
753         // recompute an accurate value, taking care of cancellations
754         final double[] bDouble = new double[b.length];
755         for (int i = 0; i < b.length; ++i) {
756             bDouble[i] = b[i].getValue();
757         }
758         out.value = Sum.ofProducts(a, bDouble).getAsDouble();
759 
760         return out;
761     }
762 
763     /** {@inheritDoc} */
764     @Override
765     public SparseGradient linearCombination(final SparseGradient a1, final SparseGradient b1,
766                                             final SparseGradient a2, final SparseGradient b2) {
767 
768         // compute a simple value, with all partial derivatives
769         SparseGradient out = a1.multiply(b1).add(a2.multiply(b2));
770 
771         // recompute an accurate value, taking care of cancellations
772         out.value = Sum.create()
773             .addProduct(a1.value, b1.value)
774             .addProduct(a2.value, b2.value).getAsDouble();
775 
776         return out;
777     }
778 
779     /** {@inheritDoc} */
780     @Override
781     public SparseGradient linearCombination(final double a1, final SparseGradient b1,
782                                             final double a2, final SparseGradient b2) {
783 
784         // compute a simple value, with all partial derivatives
785         SparseGradient out = b1.multiply(a1).add(b2.multiply(a2));
786 
787         // recompute an accurate value, taking care of cancellations
788         out.value = Sum.create()
789             .addProduct(a1, b1.value)
790             .addProduct(a2, b2.value).getAsDouble();
791 
792         return out;
793     }
794 
795     /** {@inheritDoc} */
796     @Override
797     public SparseGradient linearCombination(final SparseGradient a1, final SparseGradient b1,
798                                             final SparseGradient a2, final SparseGradient b2,
799                                             final SparseGradient a3, final SparseGradient b3) {
800 
801         // compute a simple value, with all partial derivatives
802         SparseGradient out = a1.multiply(b1).add(a2.multiply(b2)).add(a3.multiply(b3));
803 
804         // recompute an accurate value, taking care of cancellations
805         out.value = Sum.create()
806             .addProduct(a1.value, b1.value)
807             .addProduct(a2.value, b2.value)
808             .addProduct(a3.value, b3.value).getAsDouble();
809 
810         return out;
811     }
812 
813     /** {@inheritDoc} */
814     @Override
815     public SparseGradient linearCombination(final double a1, final SparseGradient b1,
816                                             final double a2, final SparseGradient b2,
817                                             final double a3, final SparseGradient b3) {
818 
819         // compute a simple value, with all partial derivatives
820         SparseGradient out = b1.multiply(a1).add(b2.multiply(a2)).add(b3.multiply(a3));
821 
822         // recompute an accurate value, taking care of cancellations
823         out.value = Sum.create()
824             .addProduct(a1, b1.value)
825             .addProduct(a2, b2.value)
826             .addProduct(a3, b3.value).getAsDouble();
827 
828         return out;
829     }
830 
831     /** {@inheritDoc} */
832     @Override
833     public SparseGradient linearCombination(final SparseGradient a1, final SparseGradient b1,
834                                             final SparseGradient a2, final SparseGradient b2,
835                                             final SparseGradient a3, final SparseGradient b3,
836                                             final SparseGradient a4, final SparseGradient b4) {
837 
838         // compute a simple value, with all partial derivatives
839         SparseGradient out = a1.multiply(b1).add(a2.multiply(b2)).add(a3.multiply(b3)).add(a4.multiply(b4));
840 
841         // recompute an accurate value, taking care of cancellations
842         out.value = Sum.create()
843             .addProduct(a1.value, b1.value)
844             .addProduct(a2.value, b2.value)
845             .addProduct(a3.value, b3.value)
846             .addProduct(a4.value, b4.value).getAsDouble();
847 
848         return out;
849     }
850 
851     /** {@inheritDoc} */
852     @Override
853     public SparseGradient linearCombination(final double a1, final SparseGradient b1,
854                                               final double a2, final SparseGradient b2,
855                                               final double a3, final SparseGradient b3,
856                                               final double a4, final SparseGradient b4) {
857 
858         // compute a simple value, with all partial derivatives
859         SparseGradient out = b1.multiply(a1).add(b2.multiply(a2)).add(b3.multiply(a3)).add(b4.multiply(a4));
860 
861         // recompute an accurate value, taking care of cancellations
862         out.value = Sum.create()
863             .addProduct(a1, b1.value)
864             .addProduct(a2, b2.value)
865             .addProduct(a3, b3.value)
866             .addProduct(a4, b4.value).getAsDouble();
867 
868         return out;
869     }
870 
871     /**
872      * Test for the equality of two sparse gradients.
873      * <p>
874      * Sparse gradients are considered equal if they have the same value
875      * and the same derivatives.
876      * </p>
877      * @param other Object to test for equality to this
878      * @return true if two sparse gradients are equal
879      */
880     @Override
881     public boolean equals(Object other) {
882 
883         if (this == other) {
884             return true;
885         }
886 
887         if (other instanceof SparseGradient) {
888             final SparseGradient rhs = (SparseGradient)other;
889             if (!Precision.equals(value, rhs.value, 1)) {
890                 return false;
891             }
892             if (derivatives.size() != rhs.derivatives.size()) {
893                 return false;
894             }
895             for (final Map.Entry<Integer, Double> entry : derivatives.entrySet()) {
896                 if (!rhs.derivatives.containsKey(entry.getKey())) {
897                     return false;
898                 }
899                 if (!Precision.equals(entry.getValue(), rhs.derivatives.get(entry.getKey()), 1)) {
900                     return false;
901                 }
902             }
903             return true;
904         }
905 
906         return false;
907     }
908 
909     /**
910      * Get a hashCode for the derivative structure.
911      * @return a hash code value for this object
912      * @since 3.2
913      */
914     @Override
915     public int hashCode() {
916         return 743 + 809 * Double.hashCode(value) + 167 * derivatives.hashCode();
917     }
918 }