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.legacy.optim.nonlinear.scalar.noderiv;
019
020import java.util.List;
021import java.util.ArrayList;
022import java.util.Arrays;
023import java.util.Comparator;
024import java.util.Collections;
025import java.util.function.UnaryOperator;
026import java.util.function.DoublePredicate;
027
028import org.apache.commons.math4.legacy.analysis.MultivariateFunction;
029import org.apache.commons.math4.legacy.exception.DimensionMismatchException;
030import org.apache.commons.math4.legacy.exception.MathIllegalArgumentException;
031import org.apache.commons.math4.legacy.exception.NotStrictlyPositiveException;
032import org.apache.commons.math4.legacy.exception.ZeroException;
033import org.apache.commons.math4.legacy.exception.util.LocalizedFormats;
034import org.apache.commons.math4.legacy.optim.OptimizationData;
035import org.apache.commons.math4.legacy.optim.PointValuePair;
036
037/**
038 * Represents a <a href="https://en.wikipedia.org/wiki/Simplex">simplex</a>.
039 *
040 * @see SimplexOptimizer
041 */
042public final class Simplex implements OptimizationData {
043    /** Coordinates. */
044    private final List<PointValuePair> points;
045
046    /**
047     * Builds from a given set of coordinates.
048     *
049     * @param referenceSimplex Reference simplex.
050     * @throws NotStrictlyPositiveException if the reference simplex does not
051     * contain at least one point.
052     * @throws DimensionMismatchException if there is a dimension mismatch
053     * in the reference simplex.
054     * @throws IllegalArgumentException if one of its vertices is duplicated.
055     */
056    private Simplex(double[][] referenceSimplex) {
057        if (referenceSimplex.length <= 0) {
058            throw new NotStrictlyPositiveException(LocalizedFormats.SIMPLEX_NEED_ONE_POINT,
059                                                   referenceSimplex.length);
060        }
061        final int len = referenceSimplex.length;
062        points = new ArrayList<>(len);
063
064        final int dim = len - 1;
065
066        // Loop over vertices.
067        for (int i = 0; i < len; i++) {
068            final double[] refI = referenceSimplex[i];
069
070            // Safety checks.
071            if (refI.length != dim) {
072                throw new DimensionMismatchException(refI.length, dim);
073            }
074            for (int j = 1; j < i; j++) {
075                final double[] refJ = referenceSimplex[j];
076                boolean allEquals = true;
077                for (int k = 0; k < dim; k++) {
078                    if (refI[k] != refJ[k]) {
079                        allEquals = false;
080                        break;
081                    }
082                }
083                if (allEquals) {
084                    throw new MathIllegalArgumentException(LocalizedFormats.EQUAL_VERTICES_IN_SIMPLEX,
085                                                           i, j);
086                }
087            }
088
089            points.add(new PointValuePair(refI, Double.NaN));
090        }
091    }
092
093    /**
094     * Builds from an existing simplex.
095     *
096     * @param points Simplex data.  Reference will be stored in the newly
097     * constructed instance.
098     */
099    private Simplex(List<PointValuePair> points) {
100        this.points = points;
101    }
102
103    /**
104     * Builds from a given set of coordinates.
105     *
106     * @param simplex Simplex coordinates.
107     * @return a new instance.
108     * @throws NotStrictlyPositiveException if the reference simplex does not
109     * contain at least one point.
110     * @throws DimensionMismatchException if there is a dimension mismatch
111     * in the reference simplex.
112     * @throws IllegalArgumentException if one of its vertices is duplicated.
113     */
114    public static Simplex of(double[][] simplex) {
115        return new Simplex(simplex);
116    }
117
118    /**
119     * Builds simplex with the given side length.
120     *
121     * @param dim Space dimensions.
122     * @param sideLength Length of the sides of the hypercube.
123     * @return a new instance.
124     */
125    public static Simplex equalSidesAlongAxes(int dim,
126                                              double sideLength) {
127        final double[] steps = new double[dim];
128        Arrays.fill(steps, sideLength);
129        return alongAxes(steps);
130    }
131
132    /**
133     * The start configuration for simplex is built from a box parallel to
134     * the canonical axes of the space. The simplex is the subset of vertices
135     * of a box parallel to the canonical axes. It is built as the path followed
136     * while traveling from one vertex of the box to the diagonally opposite
137     * vertex moving only along the box edges. The first vertex of the box will
138     * be located at the origin of the coordinate system.
139     *
140     * To be used for simplex-based optimization, the simplex must be
141     * {@link #translate(double[]) translated} so that its first vertex will be
142     * the {@link org.apache.commons.math4.legacy.optim.InitialGuess initial guess}.
143     *
144     * For example, in dimension 3 a simplex has 4 vertices. Setting the
145     * steps to (1, 10, 2) and the start point to (1, 1, 1) would imply the
146     * initial simplex would be:
147     * <ol>
148     *  <li>(1, 1, 1),</li>
149     *  <li>(2, 1, 1),</li>
150     *  <li>(2, 11, 1),</li>
151     *  <li>(2, 11, 3).</li>
152     * </ol>
153     *
154     * @param steps Steps along the canonical axes representing box edges.
155     * They may be negative but not zero.
156     * @throws ZeroException if one of the steps is zero.
157     * @return a new instance.
158     */
159    public static Simplex alongAxes(double[] steps) {
160        if (steps.length == 0) {
161            throw new ZeroException();
162        }
163        final int dim = steps.length;
164        final int len = dim + 1;
165
166        // Only the relative position of the n final vertices with respect
167        // to the first one are stored.
168        final double[][] simplex = new double[len][dim];
169        for (int i = 1; i < len; i++) { // First point is the origin (zero).
170            final double[] vertexI = simplex[i];
171            for (int j = 0; j < i; j++) {
172                if (steps[j] == 0) {
173                    throw new ZeroException();
174                }
175                System.arraycopy(steps, 0, vertexI, 0, j + 1);
176            }
177        }
178
179        return new Simplex(simplex);
180    }
181
182    /**
183     * Returns the space dimension.
184     *
185     * @return the dimension of the simplex.
186     */
187    public int getDimension() {
188        return points.size() - 1;
189    }
190
191    /**
192     * Returns the number of vertices.
193     *
194     * @return the size of the simplex.
195     */
196    public int getSize() {
197        return points.size();
198    }
199
200    /**
201     * Evaluates the (non-evaluated) simplex points and returns a new instance
202     * with vertices sorted from best to worst.
203     *
204     * @param function Evaluation function.
205     * @param comparator Comparator for sorting vertices, from best to worst.
206     * @return a new instance in which the vertices are sorted according to
207     * the given {@code comparator}.
208     */
209    public Simplex evaluate(MultivariateFunction function,
210                            Comparator<PointValuePair> comparator) {
211        final List<PointValuePair> newPoints = new ArrayList<>(points.size());
212        for (PointValuePair pv : points) {
213            final double[] coord = pv.getPoint();
214            final double value = Double.isNaN(pv.getValue()) ?
215                function.value(coord) :
216                pv.getValue();
217
218            newPoints.add(new PointValuePair(coord, value, false));
219        }
220
221        Collections.sort(newPoints, comparator);
222        return new Simplex(newPoints);
223    }
224
225    /**
226     * Retrieves a copy of the simplex point stored at {@code index}.
227     *
228     * @param index Location.
229     * @return the point at location {@code index}.
230     */
231    public PointValuePair get(int index) {
232        final PointValuePair p = points.get(index);
233        return new PointValuePair(p.getPoint(), p.getValue());
234    }
235
236    /**
237     * Creates a (deep) copy of the simplex points.
238     *
239     * @return the points.
240     */
241    public List<PointValuePair> asList() {
242        return asList(0, points.size());
243    }
244
245    /**
246     * Generator of simplex transform.
247     *
248     * @see MultiDirectionalTransform
249     * @see NelderMeadTransform
250     * @see HedarFukushimaTransform
251     */
252    public interface TransformFactory extends OptimizationData {
253        /**
254         * Creates a simplex transformation.
255         *
256         * @param evaluationFunction Evaluation function.
257         * @param comparator Vertex fitness comparator.
258         * @param saAcceptance Simulated annealing acceptance test.
259         * @return the simplex transform operator.
260         */
261        UnaryOperator<Simplex> create(MultivariateFunction evaluationFunction,
262                                      Comparator<PointValuePair> comparator,
263                                      DoublePredicate saAcceptance);
264    }
265
266    /**
267     * Creates a (deep) copy of the simplex points within slots
268     * {@code from} (included) and {@code to} (excluded).
269     *
270     * @param from Index of the first point to retrieve.
271     * @param to One past the index of the last point to retrieve.
272     * @return the points.
273     * @throws IllegalArgumentException if {@code from} and {@code to} are
274     * not within the {@code [0, n + 1]} interval (where {@code n} is the
275     * space dimension) or {@code from > to}.
276     */
277    /* package private */ List<PointValuePair> asList(int from,
278                                                      int to) {
279        if (from < 0 ||
280            to > points.size() ||
281            from > to) {
282            throw new IllegalArgumentException("Index");
283        }
284
285        final int len = to - from;
286        final List<PointValuePair> copy = new ArrayList<>(len);
287        for (int i = from; i < to; i++) {
288            copy.add(get(i));
289        }
290
291        return copy;
292    }
293
294    /**
295     * Utility for evaluating a point with coordinates \( a_i + s (b_i - a_i) \).
296     *
297     * @param a Cartesian coordinates.
298     * @param s Scaling factor.
299     * @param b Cartesian coordinates.
300     * @param function Evaluation function.
301     * @return a new point.
302     */
303    /* package private */ static PointValuePair newPoint(double[] a,
304                                                         double s,
305                                                         double[] b,
306                                                         MultivariateFunction function) {
307        final int dim = a.length;
308        final double[] r = new double[dim];
309        for (int i = 0; i < dim; i++) {
310            final double m = a[i];
311            r[i] = m + s * (b[i] - m);
312        }
313
314        return new PointValuePair(r, function.value(r), false);
315    }
316
317    /**
318     * Utility for the "shrinking" a simplex: All the points will be
319     * transformed except the one at index 0.
320     *
321     * @param sigma Shrink factor.
322     * @param function Evaluation function.
323     * @return a new instance.
324     */
325    /* package private */ Simplex shrink(double sigma,
326                                         MultivariateFunction function) {
327        final int replSize = getSize() - 1;
328        final List<PointValuePair> replacement = new ArrayList<>();
329        final double[] bestPoint = get(0).getPoint();
330        for (int i = 0; i < replSize; i++) {
331            replacement.add(Simplex.newPoint(bestPoint,
332                                             sigma,
333                                             get(i + 1).getPoint(),
334                                             function));
335        }
336
337        return replaceLast(replacement);
338    }
339
340    /**
341     * Translates the simplex such that the first point's new coordinates
342     * will be at the given {@code point}.
343     *
344     * @param point Coordinates of the new simplex's first point.
345     * @return the translated points.
346     * @throws DimensionMismatchException if the dimensions do not match.
347     */
348    /* package private */ Simplex translate(double[] point) {
349        final int dim = point.length;
350        if (getDimension() != dim) {
351            throw new DimensionMismatchException(getDimension(), dim);
352        }
353        final int len = points.size();
354        final double[][] coordinates = new double[len][dim];
355        final double[] current0 = points.get(0).getPoint(); // Current first point.
356
357        // Set new vertices.
358        for (int i = 0; i < len; i++) {
359            final double[] currentI = points.get(i).getPoint();
360
361            final double[] newI = coordinates[i];
362            for (int k = 0; k < dim; k++) {
363                newI[k] = point[k] + currentI[k] - current0[k];
364            }
365        }
366
367        return new Simplex(coordinates);
368    }
369
370    /**
371     * Creates a new simplex where the given {@code point} replaces the one at the
372     * last position.
373     * Caveat: No check is done that the resulting set of points forms is a simplex.
374     *
375     * @param point Point.
376     * @return a new instance.
377     */
378    /* package private */ Simplex replaceLast(PointValuePair point) {
379        final List<PointValuePair> newPoints = asList(0, getDimension()); // Deep copy.
380        newPoints.add(new PointValuePair(point.getPoint(), // Deep copy.
381                                         point.getValue(),
382                                         false));
383
384        return new Simplex(newPoints);
385    }
386
387    /**
388     * Replace the last points of the simplex with the points from the given
389     * {@code replacement} list.
390     * Caveat: No check is done that the resulting set of points is a simplex.
391     *
392     * @param replacement List of points that will replace the last points of
393     * the {@code simplex}.
394     * @return a new instance.
395     */
396    /* package private */ Simplex replaceLast(List<PointValuePair> replacement) {
397        final int nPoints = replacement.size();
398        final int from = points.size() - nPoints;
399        final List<PointValuePair> newPoints = asList(0, from); // Deep copy.
400
401
402        for (int i = 0; i < nPoints; i++) {
403            final PointValuePair p = replacement.get(i);
404            newPoints.add(new PointValuePair(p.getPoint(), // Deep copy.
405                                             p.getValue(),
406                                             false));
407        }
408
409        return new Simplex(newPoints);
410    }
411
412    /**
413     * @param list List of simplex points.
414     * @return the centroid of the points in the given {@code list}.
415     */
416    /* package private */ static double[] centroid(List<PointValuePair> list) {
417        final double[] centroid = list.get(0).getPoint();
418
419        final int nPoints = list.size();
420        final int dim = centroid.length;
421        for (int i = 1; i < nPoints; i++) {
422            final double[] p = list.get(i).getPoint();
423            for (int k = 0; k < dim; k++) {
424                centroid[k] += p[k];
425            }
426        }
427
428        for (int k = 0; k < dim; k++) {
429            centroid[k] /= nPoints;
430        }
431
432        return centroid;
433    }
434}