001/*
002 * Licensed to the Apache Software Foundation (ASF) under one or more
003 * contributor license agreements.  See the NOTICE file distributed with
004 * this work for additional information regarding copyright ownership.
005 * The ASF licenses this file to You under the Apache License, Version 2.0
006 * (the "License"); you may not use this file except in compliance with
007 * the License.  You may obtain a copy of the License at
008 *
009 *      http://www.apache.org/licenses/LICENSE-2.0
010 *
011 * Unless required by applicable law or agreed to in writing, software
012 * distributed under the License is distributed on an "AS IS" BASIS,
013 * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
014 * See the License for the specific language governing permissions and
015 * limitations under the License.
016 */
017
018package org.apache.commons.math4.util;
019
020import java.lang.reflect.Array;
021import java.util.ArrayList;
022import java.util.Arrays;
023import java.util.Collections;
024import java.util.Comparator;
025import java.util.Iterator;
026import java.util.List;
027import java.util.TreeSet;
028
029import org.apache.commons.numbers.core.Precision;
030import org.apache.commons.math4.Field;
031import org.apache.commons.math4.exception.DimensionMismatchException;
032import org.apache.commons.math4.exception.MathArithmeticException;
033import org.apache.commons.math4.exception.MathIllegalArgumentException;
034import org.apache.commons.math4.exception.MathInternalError;
035import org.apache.commons.math4.exception.NoDataException;
036import org.apache.commons.math4.exception.NonMonotonicSequenceException;
037import org.apache.commons.math4.exception.NotANumberException;
038import org.apache.commons.math4.exception.NotPositiveException;
039import org.apache.commons.math4.exception.NotStrictlyPositiveException;
040import org.apache.commons.math4.exception.NullArgumentException;
041import org.apache.commons.math4.exception.NumberIsTooLargeException;
042import org.apache.commons.math4.exception.util.LocalizedFormats;
043
044/**
045 * Arrays utilities.
046 *
047 * @since 3.0
048 */
049public class MathArrays {
050
051    /**
052     * Private constructor.
053     */
054    private MathArrays() {}
055
056    /**
057     * Real-valued function that operate on an array or a part of it.
058     * @since 3.1
059     */
060    public interface Function {
061        /**
062         * Operates on an entire array.
063         *
064         * @param array Array to operate on.
065         * @return the result of the operation.
066         */
067        double evaluate(double[] array);
068        /**
069         * @param array Array to operate on.
070         * @param startIndex Index of the first element to take into account.
071         * @param numElements Number of elements to take into account.
072         * @return the result of the operation.
073         */
074        double evaluate(double[] array,
075                        int startIndex,
076                        int numElements);
077    }
078
079    /**
080     * Create a copy of an array scaled by a value.
081     *
082     * @param arr Array to scale.
083     * @param val Scalar.
084     * @return scaled copy of array with each entry multiplied by val.
085     * @since 3.2
086     */
087    public static double[] scale(double val, final double[] arr) {
088        double[] newArr = new double[arr.length];
089        for (int i = 0; i < arr.length; i++) {
090            newArr[i] = arr[i] * val;
091        }
092        return newArr;
093    }
094
095    /**
096     * <p>Multiply each element of an array by a value.</p>
097     *
098     * <p>The array is modified in place (no copy is created).</p>
099     *
100     * @param arr Array to scale
101     * @param val Scalar
102     * @since 3.2
103     */
104    public static void scaleInPlace(double val, final double[] arr) {
105        for (int i = 0; i < arr.length; i++) {
106            arr[i] *= val;
107        }
108    }
109
110    /**
111     * Creates an array whose contents will be the element-by-element
112     * addition of the arguments.
113     *
114     * @param a First term of the addition.
115     * @param b Second term of the addition.
116     * @return a new array {@code r} where {@code r[i] = a[i] + b[i]}.
117     * @throws DimensionMismatchException if the array lengths differ.
118     * @since 3.1
119     */
120    public static double[] ebeAdd(double[] a, double[] b)
121        throws DimensionMismatchException {
122        checkEqualLength(a, b);
123
124        final double[] result = a.clone();
125        for (int i = 0; i < a.length; i++) {
126            result[i] += b[i];
127        }
128        return result;
129    }
130    /**
131     * Creates an array whose contents will be the element-by-element
132     * subtraction of the second argument from the first.
133     *
134     * @param a First term.
135     * @param b Element to be subtracted.
136     * @return a new array {@code r} where {@code r[i] = a[i] - b[i]}.
137     * @throws DimensionMismatchException if the array lengths differ.
138     * @since 3.1
139     */
140    public static double[] ebeSubtract(double[] a, double[] b)
141        throws DimensionMismatchException {
142        checkEqualLength(a, b);
143
144        final double[] result = a.clone();
145        for (int i = 0; i < a.length; i++) {
146            result[i] -= b[i];
147        }
148        return result;
149    }
150    /**
151     * Creates an array whose contents will be the element-by-element
152     * multiplication of the arguments.
153     *
154     * @param a First factor of the multiplication.
155     * @param b Second factor of the multiplication.
156     * @return a new array {@code r} where {@code r[i] = a[i] * b[i]}.
157     * @throws DimensionMismatchException if the array lengths differ.
158     * @since 3.1
159     */
160    public static double[] ebeMultiply(double[] a, double[] b)
161        throws DimensionMismatchException {
162        checkEqualLength(a, b);
163
164        final double[] result = a.clone();
165        for (int i = 0; i < a.length; i++) {
166            result[i] *= b[i];
167        }
168        return result;
169    }
170    /**
171     * Creates an array whose contents will be the element-by-element
172     * division of the first argument by the second.
173     *
174     * @param a Numerator of the division.
175     * @param b Denominator of the division.
176     * @return a new array {@code r} where {@code r[i] = a[i] / b[i]}.
177     * @throws DimensionMismatchException if the array lengths differ.
178     * @since 3.1
179     */
180    public static double[] ebeDivide(double[] a, double[] b)
181        throws DimensionMismatchException {
182        checkEqualLength(a, b);
183
184        final double[] result = a.clone();
185        for (int i = 0; i < a.length; i++) {
186            result[i] /= b[i];
187        }
188        return result;
189    }
190
191    /**
192     * Calculates the L<sub>1</sub> (sum of abs) distance between two points.
193     *
194     * @param p1 the first point
195     * @param p2 the second point
196     * @return the L<sub>1</sub> distance between the two points
197     * @throws DimensionMismatchException if the array lengths differ.
198     */
199    public static double distance1(double[] p1, double[] p2)
200    throws DimensionMismatchException {
201        checkEqualLength(p1, p2);
202        double sum = 0;
203        for (int i = 0; i < p1.length; i++) {
204            sum += FastMath.abs(p1[i] - p2[i]);
205        }
206        return sum;
207    }
208
209    /**
210     * Calculates the L<sub>1</sub> (sum of abs) distance between two points.
211     *
212     * @param p1 the first point
213     * @param p2 the second point
214     * @return the L<sub>1</sub> distance between the two points
215     * @throws DimensionMismatchException if the array lengths differ.
216     */
217    public static int distance1(int[] p1, int[] p2)
218    throws DimensionMismatchException {
219        checkEqualLength(p1, p2);
220        int sum = 0;
221        for (int i = 0; i < p1.length; i++) {
222            sum += FastMath.abs(p1[i] - p2[i]);
223        }
224        return sum;
225    }
226
227    /**
228     * Calculates the L<sub>2</sub> (Euclidean) distance between two points.
229     *
230     * @param p1 the first point
231     * @param p2 the second point
232     * @return the L<sub>2</sub> distance between the two points
233     * @throws DimensionMismatchException if the array lengths differ.
234     */
235    public static double distance(double[] p1, double[] p2)
236    throws DimensionMismatchException {
237        checkEqualLength(p1, p2);
238        double sum = 0;
239        for (int i = 0; i < p1.length; i++) {
240            final double dp = p1[i] - p2[i];
241            sum += dp * dp;
242        }
243        return FastMath.sqrt(sum);
244    }
245
246    /**
247     * Calculates the L<sub>2</sub> (Euclidean) distance between two points.
248     *
249     * @param p1 the first point
250     * @param p2 the second point
251     * @return the L<sub>2</sub> distance between the two points
252     * @throws DimensionMismatchException if the array lengths differ.
253     */
254    public static double distance(int[] p1, int[] p2)
255    throws DimensionMismatchException {
256      checkEqualLength(p1, p2);
257      double sum = 0;
258      for (int i = 0; i < p1.length; i++) {
259          final double dp = p1[i] - p2[i];
260          sum += dp * dp;
261      }
262      return FastMath.sqrt(sum);
263    }
264
265    /**
266     * Calculates the L<sub>&infin;</sub> (max of abs) distance between two points.
267     *
268     * @param p1 the first point
269     * @param p2 the second point
270     * @return the L<sub>&infin;</sub> distance between the two points
271     * @throws DimensionMismatchException if the array lengths differ.
272     */
273    public static double distanceInf(double[] p1, double[] p2)
274    throws DimensionMismatchException {
275        checkEqualLength(p1, p2);
276        double max = 0;
277        for (int i = 0; i < p1.length; i++) {
278            max = FastMath.max(max, FastMath.abs(p1[i] - p2[i]));
279        }
280        return max;
281    }
282
283    /**
284     * Calculates the L<sub>&infin;</sub> (max of abs) distance between two points.
285     *
286     * @param p1 the first point
287     * @param p2 the second point
288     * @return the L<sub>&infin;</sub> distance between the two points
289     * @throws DimensionMismatchException if the array lengths differ.
290     */
291    public static int distanceInf(int[] p1, int[] p2)
292    throws DimensionMismatchException {
293        checkEqualLength(p1, p2);
294        int max = 0;
295        for (int i = 0; i < p1.length; i++) {
296            max = FastMath.max(max, FastMath.abs(p1[i] - p2[i]));
297        }
298        return max;
299    }
300
301    /**
302     * Specification of ordering direction.
303     */
304    public enum OrderDirection {
305        /** Constant for increasing direction. */
306        INCREASING,
307        /** Constant for decreasing direction. */
308        DECREASING
309    }
310
311    /**
312     * Check that an array is monotonically increasing or decreasing.
313     *
314     * @param <T> the type of the elements in the specified array
315     * @param val Values.
316     * @param dir Ordering direction.
317     * @param strict Whether the order should be strict.
318     * @return {@code true} if sorted, {@code false} otherwise.
319     */
320    public static  <T extends Comparable<? super T>> boolean isMonotonic(T[] val,
321                                      OrderDirection dir,
322                                      boolean strict) {
323        T previous = val[0];
324        final int max = val.length;
325        for (int i = 1; i < max; i++) {
326            final int comp;
327            switch (dir) {
328            case INCREASING:
329                comp = previous.compareTo(val[i]);
330                if (strict) {
331                    if (comp >= 0) {
332                        return false;
333                    }
334                } else {
335                    if (comp > 0) {
336                        return false;
337                    }
338                }
339                break;
340            case DECREASING:
341                comp = val[i].compareTo(previous);
342                if (strict) {
343                    if (comp >= 0) {
344                        return false;
345                    }
346                } else {
347                    if (comp > 0) {
348                       return false;
349                    }
350                }
351                break;
352            default:
353                // Should never happen.
354                throw new MathInternalError();
355            }
356
357            previous = val[i];
358        }
359        return true;
360    }
361
362    /**
363     * Check that an array is monotonically increasing or decreasing.
364     *
365     * @param val Values.
366     * @param dir Ordering direction.
367     * @param strict Whether the order should be strict.
368     * @return {@code true} if sorted, {@code false} otherwise.
369     */
370    public static boolean isMonotonic(double[] val, OrderDirection dir, boolean strict) {
371        return checkOrder(val, dir, strict, false);
372    }
373
374    /**
375     * Check that both arrays have the same length.
376     *
377     * @param a Array.
378     * @param b Array.
379     * @param abort Whether to throw an exception if the check fails.
380     * @return {@code true} if the arrays have the same length.
381     * @throws DimensionMismatchException if the lengths differ and
382     * {@code abort} is {@code true}.
383     * @since 3.6
384     */
385    public static boolean checkEqualLength(double[] a,
386                                           double[] b,
387                                           boolean abort) {
388        if (a.length == b.length) {
389            return true;
390        } else {
391            if (abort) {
392                throw new DimensionMismatchException(a.length, b.length);
393            }
394            return false;
395        }
396    }
397
398    /**
399     * Check that both arrays have the same length.
400     *
401     * @param a Array.
402     * @param b Array.
403     * @throws DimensionMismatchException if the lengths differ.
404     * @since 3.6
405     */
406    public static void checkEqualLength(double[] a,
407                                        double[] b) {
408        checkEqualLength(a, b, true);
409    }
410
411
412    /**
413     * Check that both arrays have the same length.
414     *
415     * @param a Array.
416     * @param b Array.
417     * @param abort Whether to throw an exception if the check fails.
418     * @return {@code true} if the arrays have the same length.
419     * @throws DimensionMismatchException if the lengths differ and
420     * {@code abort} is {@code true}.
421     * @since 3.6
422     */
423    public static boolean checkEqualLength(int[] a,
424                                           int[] b,
425                                           boolean abort) {
426        if (a.length == b.length) {
427            return true;
428        } else {
429            if (abort) {
430                throw new DimensionMismatchException(a.length, b.length);
431            }
432            return false;
433        }
434    }
435
436    /**
437     * Check that both arrays have the same length.
438     *
439     * @param a Array.
440     * @param b Array.
441     * @throws DimensionMismatchException if the lengths differ.
442     * @since 3.6
443     */
444    public static void checkEqualLength(int[] a,
445                                        int[] b) {
446        checkEqualLength(a, b, true);
447    }
448
449    /**
450     * Check that the given array is sorted.
451     *
452     * @param val Values.
453     * @param dir Ordering direction.
454     * @param strict Whether the order should be strict.
455     * @param abort Whether to throw an exception if the check fails.
456     * @return {@code true} if the array is sorted.
457     * @throws NonMonotonicSequenceException if the array is not sorted
458     * and {@code abort} is {@code true}.
459     */
460    public static boolean checkOrder(double[] val, OrderDirection dir,
461                                     boolean strict, boolean abort)
462        throws NonMonotonicSequenceException {
463        double previous = val[0];
464        final int max = val.length;
465
466        int index;
467        ITEM:
468        for (index = 1; index < max; index++) {
469            switch (dir) {
470            case INCREASING:
471                if (strict) {
472                    if (val[index] <= previous) {
473                        break ITEM;
474                    }
475                } else {
476                    if (val[index] < previous) {
477                        break ITEM;
478                    }
479                }
480                break;
481            case DECREASING:
482                if (strict) {
483                    if (val[index] >= previous) {
484                        break ITEM;
485                    }
486                } else {
487                    if (val[index] > previous) {
488                        break ITEM;
489                    }
490                }
491                break;
492            default:
493                // Should never happen.
494                throw new MathInternalError();
495            }
496
497            previous = val[index];
498        }
499
500        if (index == max) {
501            // Loop completed.
502            return true;
503        }
504
505        // Loop early exit means wrong ordering.
506        if (abort) {
507            throw new NonMonotonicSequenceException(val[index], previous, index, dir, strict);
508        } else {
509            return false;
510        }
511    }
512
513    /**
514     * Check that the given array is sorted.
515     *
516     * @param val Values.
517     * @param dir Ordering direction.
518     * @param strict Whether the order should be strict.
519     * @throws NonMonotonicSequenceException if the array is not sorted.
520     * @since 2.2
521     */
522    public static void checkOrder(double[] val, OrderDirection dir,
523                                  boolean strict) throws NonMonotonicSequenceException {
524        checkOrder(val, dir, strict, true);
525    }
526
527    /**
528     * Check that the given array is sorted in strictly increasing order.
529     *
530     * @param val Values.
531     * @throws NonMonotonicSequenceException if the array is not sorted.
532     * @since 2.2
533     */
534    public static void checkOrder(double[] val) throws NonMonotonicSequenceException {
535        checkOrder(val, OrderDirection.INCREASING, true);
536    }
537
538    /**
539     * Throws DimensionMismatchException if the input array is not rectangular.
540     *
541     * @param in array to be tested
542     * @throws NullArgumentException if input array is null
543     * @throws DimensionMismatchException if input array is not rectangular
544     * @since 3.1
545     */
546    public static void checkRectangular(final long[][] in)
547        throws NullArgumentException, DimensionMismatchException {
548        MathUtils.checkNotNull(in);
549        for (int i = 1; i < in.length; i++) {
550            if (in[i].length != in[0].length) {
551                throw new DimensionMismatchException(
552                        LocalizedFormats.DIFFERENT_ROWS_LENGTHS,
553                        in[i].length, in[0].length);
554            }
555        }
556    }
557
558    /**
559     * Check that all entries of the input array are strictly positive.
560     *
561     * @param in Array to be tested
562     * @throws NotStrictlyPositiveException if any entries of the array are not
563     * strictly positive.
564     * @since 3.1
565     */
566    public static void checkPositive(final double[] in)
567        throws NotStrictlyPositiveException {
568        for (int i = 0; i < in.length; i++) {
569            if (in[i] <= 0) {
570                throw new NotStrictlyPositiveException(in[i]);
571            }
572        }
573    }
574
575    /**
576     * Check that no entry of the input array is {@code NaN}.
577     *
578     * @param in Array to be tested.
579     * @throws NotANumberException if an entry is {@code NaN}.
580     * @since 3.4
581     */
582    public static void checkNotNaN(final double[] in)
583        throws NotANumberException {
584        for(int i = 0; i < in.length; i++) {
585            if (Double.isNaN(in[i])) {
586                throw new NotANumberException();
587            }
588        }
589    }
590
591    /**
592     * Check that all entries of the input array are &gt;= 0.
593     *
594     * @param in Array to be tested
595     * @throws NotPositiveException if any array entries are less than 0.
596     * @since 3.1
597     */
598    public static void checkNonNegative(final long[] in)
599        throws NotPositiveException {
600        for (int i = 0; i < in.length; i++) {
601            if (in[i] < 0) {
602                throw new NotPositiveException(in[i]);
603            }
604        }
605    }
606
607    /**
608     * Check all entries of the input array are &gt;= 0.
609     *
610     * @param in Array to be tested
611     * @throws NotPositiveException if any array entries are less than 0.
612     * @since 3.1
613     */
614    public static void checkNonNegative(final long[][] in)
615        throws NotPositiveException {
616        for (int i = 0; i < in.length; i ++) {
617            for (int j = 0; j < in[i].length; j++) {
618                if (in[i][j] < 0) {
619                    throw new NotPositiveException(in[i][j]);
620                }
621            }
622        }
623    }
624
625    /**
626     * Sort an array in ascending order in place and perform the same reordering
627     * of entries on other arrays. For example, if
628     * {@code x = [3, 1, 2], y = [1, 2, 3]} and {@code z = [0, 5, 7]}, then
629     * {@code sortInPlace(x, y, z)} will update {@code x} to {@code [1, 2, 3]},
630     * {@code y} to {@code [2, 3, 1]} and {@code z} to {@code [5, 7, 0]}.
631     *
632     * @param x Array to be sorted and used as a pattern for permutation
633     * of the other arrays.
634     * @param yList Set of arrays whose permutations of entries will follow
635     * those performed on {@code x}.
636     * @throws DimensionMismatchException if any {@code y} is not the same
637     * size as {@code x}.
638     * @throws NullArgumentException if {@code x} or any {@code y} is null.
639     * @since 3.0
640     */
641    public static void sortInPlace(double[] x, double[] ... yList)
642        throws DimensionMismatchException, NullArgumentException {
643        sortInPlace(x, OrderDirection.INCREASING, yList);
644    }
645
646    /**
647     * Helper data structure holding a (double, integer) pair.
648     */
649    private static class PairDoubleInteger {
650        /** Key */
651        private final double key;
652        /** Value */
653        private final int value;
654
655        /**
656         * @param key Key.
657         * @param value Value.
658         */
659        PairDoubleInteger(double key, int value) {
660            this.key = key;
661            this.value = value;
662        }
663
664        /** @return the key. */
665        public double getKey() {
666            return key;
667        }
668
669        /** @return the value. */
670        public int getValue() {
671            return value;
672        }
673    }
674
675    /**
676     * Sort an array in place and perform the same reordering of entries on
677     * other arrays.  This method works the same as the other
678     * {@link #sortInPlace(double[], double[][]) sortInPlace} method, but
679     * allows the order of the sort to be provided in the {@code dir}
680     * parameter.
681     *
682     * @param x Array to be sorted and used as a pattern for permutation
683     * of the other arrays.
684     * @param dir Order direction.
685     * @param yList Set of arrays whose permutations of entries will follow
686     * those performed on {@code x}.
687     * @throws DimensionMismatchException if any {@code y} is not the same
688     * size as {@code x}.
689     * @throws NullArgumentException if {@code x} or any {@code y} is null
690     * @since 3.0
691     */
692    public static void sortInPlace(double[] x,
693                                   final OrderDirection dir,
694                                   double[] ... yList)
695        throws NullArgumentException,
696               DimensionMismatchException {
697
698        // Consistency checks.
699        if (x == null) {
700            throw new NullArgumentException();
701        }
702
703        final int yListLen = yList.length;
704        final int len = x.length;
705
706        for (int j = 0; j < yListLen; j++) {
707            final double[] y = yList[j];
708            if (y == null) {
709                throw new NullArgumentException();
710            }
711            if (y.length != len) {
712                throw new DimensionMismatchException(y.length, len);
713            }
714        }
715
716        // Associate each abscissa "x[i]" with its index "i".
717        final List<PairDoubleInteger> list
718            = new ArrayList<>(len);
719        for (int i = 0; i < len; i++) {
720            list.add(new PairDoubleInteger(x[i], i));
721        }
722
723        // Create comparators for increasing and decreasing orders.
724        final Comparator<PairDoubleInteger> comp
725            = dir == MathArrays.OrderDirection.INCREASING ?
726            new Comparator<PairDoubleInteger>() {
727            /** {@inheritDoc} */
728            @Override
729            public int compare(PairDoubleInteger o1,
730                               PairDoubleInteger o2) {
731                return Double.compare(o1.getKey(), o2.getKey());
732            }
733        } : new Comparator<PairDoubleInteger>() {
734            /** {@inheritDoc} */
735            @Override
736            public int compare(PairDoubleInteger o1,
737                               PairDoubleInteger o2) {
738                return Double.compare(o2.getKey(), o1.getKey());
739            }
740        };
741
742        // Sort.
743        Collections.sort(list, comp);
744
745        // Modify the original array so that its elements are in
746        // the prescribed order.
747        // Retrieve indices of original locations.
748        final int[] indices = new int[len];
749        for (int i = 0; i < len; i++) {
750            final PairDoubleInteger e = list.get(i);
751            x[i] = e.getKey();
752            indices[i] = e.getValue();
753        }
754
755        // In each of the associated arrays, move the
756        // elements to their new location.
757        for (int j = 0; j < yListLen; j++) {
758            // Input array will be modified in place.
759            final double[] yInPlace = yList[j];
760            final double[] yOrig = yInPlace.clone();
761
762            for (int i = 0; i < len; i++) {
763                yInPlace[i] = yOrig[indices[i]];
764            }
765        }
766    }
767
768    /**
769     * Returns true iff both arguments are null or have same dimensions and all
770     * their elements are equal as defined by
771     * {@link Precision#equals(float,float)}.
772     *
773     * @param x first array
774     * @param y second array
775     * @return true if the values are both null or have same dimension
776     * and equal elements.
777     */
778    public static boolean equals(float[] x, float[] y) {
779        if ((x == null) || (y == null)) {
780            return !((x == null) ^ (y == null));
781        }
782        if (x.length != y.length) {
783            return false;
784        }
785        for (int i = 0; i < x.length; ++i) {
786            if (!Precision.equals(x[i], y[i])) {
787                return false;
788            }
789        }
790        return true;
791    }
792
793    /**
794     * Returns true iff both arguments are null or have same dimensions and all
795     * their elements are equal as defined by
796     * {@link Precision#equalsIncludingNaN(double,double) this method}.
797     *
798     * @param x first array
799     * @param y second array
800     * @return true if the values are both null or have same dimension and
801     * equal elements
802     * @since 2.2
803     */
804    public static boolean equalsIncludingNaN(float[] x, float[] y) {
805        if ((x == null) || (y == null)) {
806            return !((x == null) ^ (y == null));
807        }
808        if (x.length != y.length) {
809            return false;
810        }
811        for (int i = 0; i < x.length; ++i) {
812            if (!Precision.equalsIncludingNaN(x[i], y[i])) {
813                return false;
814            }
815        }
816        return true;
817    }
818
819    /**
820     * Returns {@code true} iff both arguments are {@code null} or have same
821     * dimensions and all their elements are equal as defined by
822     * {@link Precision#equals(double,double)}.
823     *
824     * @param x First array.
825     * @param y Second array.
826     * @return {@code true} if the values are both {@code null} or have same
827     * dimension and equal elements.
828     */
829    public static boolean equals(double[] x, double[] y) {
830        if ((x == null) || (y == null)) {
831            return !((x == null) ^ (y == null));
832        }
833        if (x.length != y.length) {
834            return false;
835        }
836        for (int i = 0; i < x.length; ++i) {
837            if (!Precision.equals(x[i], y[i])) {
838                return false;
839            }
840        }
841        return true;
842    }
843
844    /**
845     * Returns {@code true} iff both arguments are {@code null} or have same
846     * dimensions and all their elements are equal as defined by
847     * {@link Precision#equalsIncludingNaN(double,double) this method}.
848     *
849     * @param x First array.
850     * @param y Second array.
851     * @return {@code true} if the values are both {@code null} or have same
852     * dimension and equal elements.
853     * @since 2.2
854     */
855    public static boolean equalsIncludingNaN(double[] x, double[] y) {
856        if ((x == null) || (y == null)) {
857            return !((x == null) ^ (y == null));
858        }
859        if (x.length != y.length) {
860            return false;
861        }
862        for (int i = 0; i < x.length; ++i) {
863            if (!Precision.equalsIncludingNaN(x[i], y[i])) {
864                return false;
865            }
866        }
867        return true;
868    }
869
870    /**
871     * Normalizes an array to make it sum to a specified value.
872     * Returns the result of the transformation
873     * <pre>
874     *    x |-&gt; x * normalizedSum / sum
875     * </pre>
876     * applied to each non-NaN element x of the input array, where sum is the
877     * sum of the non-NaN entries in the input array.
878     * <p>
879     * Throws IllegalArgumentException if {@code normalizedSum} is infinite
880     * or NaN and ArithmeticException if the input array contains any infinite elements
881     * or sums to 0.
882     * <p>
883     * Ignores (i.e., copies unchanged to the output array) NaNs in the input array.
884     *
885     * @param values Input array to be normalized
886     * @param normalizedSum Target sum for the normalized array
887     * @return the normalized array.
888     * @throws MathArithmeticException if the input array contains infinite
889     * elements or sums to zero.
890     * @throws MathIllegalArgumentException if the target sum is infinite or {@code NaN}.
891     * @since 2.1
892     */
893    public static double[] normalizeArray(double[] values, double normalizedSum)
894        throws MathIllegalArgumentException, MathArithmeticException {
895        if (Double.isInfinite(normalizedSum)) {
896            throw new MathIllegalArgumentException(LocalizedFormats.NORMALIZE_INFINITE);
897        }
898        if (Double.isNaN(normalizedSum)) {
899            throw new MathIllegalArgumentException(LocalizedFormats.NORMALIZE_NAN);
900        }
901        double sum = 0d;
902        final int len = values.length;
903        double[] out = new double[len];
904        for (int i = 0; i < len; i++) {
905            if (Double.isInfinite(values[i])) {
906                throw new MathIllegalArgumentException(LocalizedFormats.INFINITE_ARRAY_ELEMENT, values[i], i);
907            }
908            if (!Double.isNaN(values[i])) {
909                sum += values[i];
910            }
911        }
912        if (sum == 0) {
913            throw new MathArithmeticException(LocalizedFormats.ARRAY_SUMS_TO_ZERO);
914        }
915        final double scale = normalizedSum / sum;
916        for (int i = 0; i < len; i++) {
917            if (Double.isNaN(values[i])) {
918                out[i] = Double.NaN;
919            } else {
920                out[i] = values[i] * scale;
921            }
922        }
923        return out;
924    }
925
926    /** Build an array of elements.
927     * <p>
928     * Arrays are filled with field.getZero()
929     *
930     * @param <T> the type of the field elements
931     * @param field field to which array elements belong
932     * @param length of the array
933     * @return a new array
934     * @since 3.2
935     */
936    public static <T> T[] buildArray(final Field<T> field, final int length) {
937        @SuppressWarnings("unchecked") // OK because field must be correct class
938        T[] array = (T[]) Array.newInstance(field.getRuntimeClass(), length);
939        Arrays.fill(array, field.getZero());
940        return array;
941    }
942
943    /** Build a double dimension  array of elements.
944     * <p>
945     * Arrays are filled with field.getZero()
946     *
947     * @param <T> the type of the field elements
948     * @param field field to which array elements belong
949     * @param rows number of rows in the array
950     * @param columns number of columns (may be negative to build partial
951     * arrays in the same way <code>new Field[rows][]</code> works)
952     * @return a new array
953     * @since 3.2
954     */
955    @SuppressWarnings("unchecked")
956    public static <T> T[][] buildArray(final Field<T> field, final int rows, final int columns) {
957        final T[][] array;
958        if (columns < 0) {
959            T[] dummyRow = buildArray(field, 0);
960            array = (T[][]) Array.newInstance(dummyRow.getClass(), rows);
961        } else {
962            array = (T[][]) Array.newInstance(field.getRuntimeClass(),
963                                              new int[] {
964                                                  rows, columns
965                                              });
966            for (int i = 0; i < rows; ++i) {
967                Arrays.fill(array[i], field.getZero());
968            }
969        }
970        return array;
971    }
972
973    /**
974     * Calculates the <a href="http://en.wikipedia.org/wiki/Convolution">
975     * convolution</a> between two sequences.
976     * <p>
977     * The solution is obtained via straightforward computation of the
978     * convolution sum (and not via FFT). Whenever the computation needs
979     * an element that would be located at an index outside the input arrays,
980     * the value is assumed to be zero.
981     *
982     * @param x First sequence.
983     * Typically, this sequence will represent an input signal to a system.
984     * @param h Second sequence.
985     * Typically, this sequence will represent the impulse response of the system.
986     * @return the convolution of {@code x} and {@code h}.
987     * This array's length will be {@code x.length + h.length - 1}.
988     * @throws NullArgumentException if either {@code x} or {@code h} is {@code null}.
989     * @throws NoDataException if either {@code x} or {@code h} is empty.
990     *
991     * @since 3.3
992     */
993    public static double[] convolve(double[] x, double[] h)
994        throws NullArgumentException,
995               NoDataException {
996        MathUtils.checkNotNull(x);
997        MathUtils.checkNotNull(h);
998
999        final int xLen = x.length;
1000        final int hLen = h.length;
1001
1002        if (xLen == 0 || hLen == 0) {
1003            throw new NoDataException();
1004        }
1005
1006        // initialize the output array
1007        final int totalLength = xLen + hLen - 1;
1008        final double[] y = new double[totalLength];
1009
1010        // straightforward implementation of the convolution sum
1011        for (int n = 0; n < totalLength; n++) {
1012            double yn = 0;
1013            int k = FastMath.max(0, n + 1 - xLen);
1014            int j = n - k;
1015            while (k < hLen && j >= 0) {
1016                yn += x[j--] * h[k++];
1017            }
1018            y[n] = yn;
1019        }
1020
1021        return y;
1022    }
1023
1024    /**
1025     * Returns an array representing the natural number {@code n}.
1026     *
1027     * @param n Natural number.
1028     * @return an array whose entries are the numbers 0, 1, ..., {@code n}-1.
1029     * If {@code n == 0}, the returned array is empty.
1030     */
1031    public static int[] natural(int n) {
1032        return sequence(n, 0, 1);
1033    }
1034    /**
1035     * Returns an array of {@code size} integers starting at {@code start},
1036     * skipping {@code stride} numbers.
1037     *
1038     * @param size Natural number.
1039     * @param start Natural number.
1040     * @param stride Natural number.
1041     * @return an array whose entries are the numbers
1042     * {@code start, start + stride, ..., start + (size - 1) * stride}.
1043     * If {@code size == 0}, the returned array is empty.
1044     *
1045     * @since 3.4
1046     */
1047    public static int[] sequence(int size,
1048                                 int start,
1049                                 int stride) {
1050        final int[] a = new int[size];
1051        for (int i = 0; i < size; i++) {
1052            a[i] = start + i * stride;
1053        }
1054        return a;
1055    }
1056    /**
1057     * This method is used
1058     * to verify that the input parameters designate a subarray of positive length.
1059     * <ul>
1060     * <li>returns <code>true</code> iff the parameters designate a subarray of
1061     * positive length</li>
1062     * <li>throws <code>MathIllegalArgumentException</code> if the array is null or
1063     * or the indices are invalid</li>
1064     * <li>returns <code>false</code> if the array is non-null, but
1065     * <code>length</code> is 0.</li>
1066     * </ul>
1067     *
1068     * @param values the input array
1069     * @param begin index of the first array element to include
1070     * @param length the number of elements to include
1071     * @return true if the parameters are valid and designate a subarray of positive length
1072     * @throws MathIllegalArgumentException if the indices are invalid or the array is null
1073     * @since 3.3
1074     */
1075    public static boolean verifyValues(final double[] values, final int begin, final int length)
1076            throws MathIllegalArgumentException {
1077        return verifyValues(values, begin, length, false);
1078    }
1079
1080    /**
1081     * This method is used
1082     * to verify that the input parameters designate a subarray of positive length.
1083     * <ul>
1084     * <li>returns <code>true</code> iff the parameters designate a subarray of
1085     * non-negative length</li>
1086     * <li>throws <code>IllegalArgumentException</code> if the array is null or
1087     * or the indices are invalid</li>
1088     * <li>returns <code>false</code> if the array is non-null, but
1089     * <code>length</code> is 0 unless <code>allowEmpty</code> is <code>true</code></li>
1090     * </ul>
1091     *
1092     * @param values the input array
1093     * @param begin index of the first array element to include
1094     * @param length the number of elements to include
1095     * @param allowEmpty if <code>true</code> then zero length arrays are allowed
1096     * @return true if the parameters are valid
1097     * @throws MathIllegalArgumentException if the indices are invalid or the array is null
1098     * @since 3.3
1099     */
1100    public static boolean verifyValues(final double[] values, final int begin,
1101            final int length, final boolean allowEmpty) throws MathIllegalArgumentException {
1102
1103        if (values == null) {
1104            throw new NullArgumentException(LocalizedFormats.INPUT_ARRAY);
1105        }
1106
1107        if (begin < 0) {
1108            throw new NotPositiveException(LocalizedFormats.START_POSITION, Integer.valueOf(begin));
1109        }
1110
1111        if (length < 0) {
1112            throw new NotPositiveException(LocalizedFormats.LENGTH, Integer.valueOf(length));
1113        }
1114
1115        if (begin + length > values.length) {
1116            throw new NumberIsTooLargeException(LocalizedFormats.SUBARRAY_ENDS_AFTER_ARRAY_END,
1117                    Integer.valueOf(begin + length), Integer.valueOf(values.length), true);
1118        }
1119
1120        if (length == 0 && !allowEmpty) {
1121            return false;
1122        }
1123
1124        return true;
1125
1126    }
1127
1128    /**
1129     * This method is used
1130     * to verify that the begin and length parameters designate a subarray of positive length
1131     * and the weights are all non-negative, non-NaN, finite, and not all zero.
1132     * <ul>
1133     * <li>returns <code>true</code> iff the parameters designate a subarray of
1134     * positive length and the weights array contains legitimate values.</li>
1135     * <li>throws <code>IllegalArgumentException</code> if any of the following are true:
1136     * <ul><li>the values array is null</li>
1137     *     <li>the weights array is null</li>
1138     *     <li>the weights array does not have the same length as the values array</li>
1139     *     <li>the weights array contains one or more infinite values</li>
1140     *     <li>the weights array contains one or more NaN values</li>
1141     *     <li>the weights array contains negative values</li>
1142     *     <li>the start and length arguments do not determine a valid array</li></ul>
1143     * </li>
1144     * <li>returns <code>false</code> if the array is non-null, but
1145     * <code>length</code> is 0.</li>
1146     * </ul>
1147     *
1148     * @param values the input array
1149     * @param weights the weights array
1150     * @param begin index of the first array element to include
1151     * @param length the number of elements to include
1152     * @return true if the parameters are valid and designate a subarray of positive length
1153     * @throws MathIllegalArgumentException if the indices are invalid or the array is null
1154     * @since 3.3
1155     */
1156    public static boolean verifyValues(
1157        final double[] values,
1158        final double[] weights,
1159        final int begin,
1160        final int length) throws MathIllegalArgumentException {
1161        return verifyValues(values, weights, begin, length, false);
1162    }
1163
1164    /**
1165     * This method is used
1166     * to verify that the begin and length parameters designate a subarray of positive length
1167     * and the weights are all non-negative, non-NaN, finite, and not all zero.
1168     * <ul>
1169     * <li>returns <code>true</code> iff the parameters designate a subarray of
1170     * non-negative length and the weights array contains legitimate values.</li>
1171     * <li>throws <code>MathIllegalArgumentException</code> if any of the following are true:
1172     * <ul><li>the values array is null</li>
1173     *     <li>the weights array is null</li>
1174     *     <li>the weights array does not have the same length as the values array</li>
1175     *     <li>the weights array contains one or more infinite values</li>
1176     *     <li>the weights array contains one or more NaN values</li>
1177     *     <li>the weights array contains negative values</li>
1178     *     <li>the start and length arguments do not determine a valid array</li></ul>
1179     * </li>
1180     * <li>returns <code>false</code> if the array is non-null, but
1181     * <code>length</code> is 0 unless <code>allowEmpty</code> is <code>true</code>.</li>
1182     * </ul>
1183     *
1184     * @param values the input array.
1185     * @param weights the weights array.
1186     * @param begin index of the first array element to include.
1187     * @param length the number of elements to include.
1188     * @param allowEmpty if {@code true} than allow zero length arrays to pass.
1189     * @return {@code true} if the parameters are valid.
1190     * @throws NullArgumentException if either of the arrays are null
1191     * @throws MathIllegalArgumentException if the array indices are not valid,
1192     * the weights array contains NaN, infinite or negative elements, or there
1193     * are no positive weights.
1194     * @since 3.3
1195     */
1196    public static boolean verifyValues(final double[] values, final double[] weights,
1197            final int begin, final int length, final boolean allowEmpty) throws MathIllegalArgumentException {
1198
1199        if (weights == null || values == null) {
1200            throw new NullArgumentException(LocalizedFormats.INPUT_ARRAY);
1201        }
1202
1203        checkEqualLength(weights, values);
1204
1205        boolean containsPositiveWeight = false;
1206        for (int i = begin; i < begin + length; i++) {
1207            final double weight = weights[i];
1208            if (Double.isNaN(weight)) {
1209                throw new MathIllegalArgumentException(LocalizedFormats.NAN_ELEMENT_AT_INDEX, Integer.valueOf(i));
1210            }
1211            if (Double.isInfinite(weight)) {
1212                throw new MathIllegalArgumentException(LocalizedFormats.INFINITE_ARRAY_ELEMENT, Double.valueOf(weight), Integer.valueOf(i));
1213            }
1214            if (weight < 0) {
1215                throw new MathIllegalArgumentException(LocalizedFormats.NEGATIVE_ELEMENT_AT_INDEX, Integer.valueOf(i), Double.valueOf(weight));
1216            }
1217            if (!containsPositiveWeight && weight > 0.0) {
1218                containsPositiveWeight = true;
1219            }
1220        }
1221
1222        if (!containsPositiveWeight) {
1223            throw new MathIllegalArgumentException(LocalizedFormats.WEIGHT_AT_LEAST_ONE_NON_ZERO);
1224        }
1225
1226        return verifyValues(values, begin, length, allowEmpty);
1227    }
1228
1229    /**
1230     * Concatenates a sequence of arrays. The return array consists of the
1231     * entries of the input arrays concatenated in the order they appear in
1232     * the argument list.  Null arrays cause NullPointerExceptions; zero
1233     * length arrays are allowed (contributing nothing to the output array).
1234     *
1235     * @param x list of double[] arrays to concatenate
1236     * @return a new array consisting of the entries of the argument arrays
1237     * @throws NullPointerException if any of the arrays are null
1238     * @since 3.6
1239     */
1240    public static double[] concatenate(double[] ...x) {
1241        int combinedLength = 0;
1242        for (double[] a : x) {
1243            combinedLength += a.length;
1244        }
1245        int offset = 0;
1246        int curLength = 0;
1247        final double[] combined = new double[combinedLength];
1248        for (int i = 0; i < x.length; i++) {
1249            curLength = x[i].length;
1250            System.arraycopy(x[i], 0, combined, offset, curLength);
1251            offset += curLength;
1252        }
1253        return combined;
1254    }
1255
1256    /**
1257     * Returns an array consisting of the unique values in {@code data}.
1258     * The return array is sorted in descending order.  Empty arrays
1259     * are allowed, but null arrays result in NullPointerException.
1260     * Infinities are allowed.  NaN values are allowed with maximum
1261     * sort order - i.e., if there are NaN values in {@code data},
1262     * {@code Double.NaN} will be the first element of the output array,
1263     * even if the array also contains {@code Double.POSITIVE_INFINITY}.
1264     *
1265     * @param data array to scan
1266     * @return descending list of values included in the input array
1267     * @throws NullPointerException if data is null
1268     * @since 3.6
1269     */
1270    public static double[] unique(double[] data) {
1271        TreeSet<Double> values = new TreeSet<>();
1272        for (int i = 0; i < data.length; i++) {
1273            values.add(data[i]);
1274        }
1275        final int count = values.size();
1276        final double[] out = new double[count];
1277        Iterator<Double> iterator = values.descendingIterator();
1278        int i = 0;
1279        while (iterator.hasNext()) {
1280            out[i++] = iterator.next();
1281        }
1282        return out;
1283    }
1284}