    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
    018    package org.apache.commons.math3.optimization.direct;
    019
    020    import java.util.Comparator;
    021
    022    import org.apache.commons.math3.optimization.PointValuePair;
    023    import org.apache.commons.math3.analysis.MultivariateFunction;
    024
    025    /**
    026     * This class implements the Nelder-Mead simplex algorithm.
    027     *
    028     * @version $Id: NelderMeadSimplex.java 1422230 2012-12-15 12:11:13Z erans $
    029     * @deprecated As of 3.1 (to be removed in 4.0).
    030     * @since 3.0
    031     */
    032    @Deprecated
    033    public class NelderMeadSimplex extends AbstractSimplex {
    034        /** Default value for {@link #rho}: {@value}. */
    035        private static final double DEFAULT_RHO = 1;
    036        /** Default value for {@link #khi}: {@value}. */
    037        private static final double DEFAULT_KHI = 2;
    038        /** Default value for {@link #gamma}: {@value}. */
    039        private static final double DEFAULT_GAMMA = 0.5;
    040        /** Default value for {@link #sigma}: {@value}. */
    041        private static final double DEFAULT_SIGMA = 0.5;
    042        /** Reflection coefficient. */
    043        private final double rho;
    044        /** Expansion coefficient. */
    045        private final double khi;
    046        /** Contraction coefficient. */
    047        private final double gamma;
    048        /** Shrinkage coefficient. */
    049        private final double sigma;
    050
    051        /**
    052         * Build a Nelder-Mead simplex with default coefficients.
    053         * The default coefficients are 1.0 for rho, 2.0 for khi and 0.5
    054         * for both gamma and sigma.
    055         *
    056         * @param n Dimension of the simplex.
    057         */
    058        public NelderMeadSimplex(final int n) {
    059            this(n, 1d);
    060        }
    061
    062        /**
    063         * Build a Nelder-Mead simplex with default coefficients.
    064         * The default coefficients are 1.0 for rho, 2.0 for khi and 0.5
    065         * for both gamma and sigma.
    066         *
    067         * @param n Dimension of the simplex.
    068         * @param sideLength Length of the sides of the default (hypercube)
    069         * simplex. See {@link AbstractSimplex#AbstractSimplex(int,double)}.
    070         */
    071        public NelderMeadSimplex(final int n, double sideLength) {
    072            this(n, sideLength,
    073                 DEFAULT_RHO, DEFAULT_KHI, DEFAULT_GAMMA, DEFAULT_SIGMA);
    074        }
    075
    076        /**
    077         * Build a Nelder-Mead simplex with specified coefficients.
    078         *
    079         * @param n Dimension of the simplex. See
    080         * {@link AbstractSimplex#AbstractSimplex(int,double)}.
    081         * @param sideLength Length of the sides of the default (hypercube)
    082         * simplex. See {@link AbstractSimplex#AbstractSimplex(int,double)}.
    083         * @param rho Reflection coefficient.
    084         * @param khi Expansion coefficient.
    085         * @param gamma Contraction coefficient.
    086         * @param sigma Shrinkage coefficient.
    087         */
    088        public NelderMeadSimplex(final int n, double sideLength,
    089                                 final double rho, final double khi,
    090                                 final double gamma, final double sigma) {
    091            super(n, sideLength);
    092
    093            this.rho = rho;
    094            this.khi = khi;
    095            this.gamma = gamma;
    096            this.sigma = sigma;
    097        }
    098
    099        /**
    100         * Build a Nelder-Mead simplex with specified coefficients.
    101         *
    102         * @param n Dimension of the simplex. See
    103         * {@link AbstractSimplex#AbstractSimplex(int)}.
    104         * @param rho Reflection coefficient.
    105         * @param khi Expansion coefficient.
    106         * @param gamma Contraction coefficient.
    107         * @param sigma Shrinkage coefficient.
    108         */
    109        public NelderMeadSimplex(final int n,
    110                                 final double rho, final double khi,
    111                                 final double gamma, final double sigma) {
    112            this(n, 1d, rho, khi, gamma, sigma);
    113        }
    114
    115        /**
    116         * Build a Nelder-Mead simplex with default coefficients.
    117         * The default coefficients are 1.0 for rho, 2.0 for khi and 0.5
    118         * for both gamma and sigma.
    119         *
    120         * @param steps Steps along the canonical axes representing box edges.
    121         * They may be negative but not zero. See
    122         */
    123        public NelderMeadSimplex(final double[] steps) {
    124            this(steps, DEFAULT_RHO, DEFAULT_KHI, DEFAULT_GAMMA, DEFAULT_SIGMA);
    125        }
    126
    127        /**
    128         * Build a Nelder-Mead simplex with specified coefficients.
    129         *
    130         * @param steps Steps along the canonical axes representing box edges.
    131         * They may be negative but not zero. See
    132         * {@link AbstractSimplex#AbstractSimplex(double[])}.
    133         * @param rho Reflection coefficient.
    134         * @param khi Expansion coefficient.
    135         * @param gamma Contraction coefficient.
    136         * @param sigma Shrinkage coefficient.
    137         * @throws IllegalArgumentException if one of the steps is zero.
    138         */
    139        public NelderMeadSimplex(final double[] steps,
    140                                 final double rho, final double khi,
    141                                 final double gamma, final double sigma) {
    142            super(steps);
    143
    144            this.rho = rho;
    145            this.khi = khi;
    146            this.gamma = gamma;
    147            this.sigma = sigma;
    148        }
    149
    150        /**
    151         * Build a Nelder-Mead simplex with default coefficients.
    152         * The default coefficients are 1.0 for rho, 2.0 for khi and 0.5
    153         * for both gamma and sigma.
    154         *
    155         * @param referenceSimplex Reference simplex. See
    156         * {@link AbstractSimplex#AbstractSimplex(double[][])}.
    157         */
    158        public NelderMeadSimplex(final double[][] referenceSimplex) {
    159            this(referenceSimplex, DEFAULT_RHO, DEFAULT_KHI, DEFAULT_GAMMA, DEFAULT_SIGMA);
    160        }
    161
    162        /**
    163         * Build a Nelder-Mead simplex with specified coefficients.
    164         *
    165         * @param referenceSimplex Reference simplex. See
    166         * {@link AbstractSimplex#AbstractSimplex(double[][])}.
    167         * @param rho Reflection coefficient.
    168         * @param khi Expansion coefficient.
    169         * @param gamma Contraction coefficient.
    170         * @param sigma Shrinkage coefficient.
    171         * @throws org.apache.commons.math3.exception.NotStrictlyPositiveException
    172         * if the reference simplex does not contain at least one point.
    173         * @throws org.apache.commons.math3.exception.DimensionMismatchException
    174         * if there is a dimension mismatch in the reference simplex.
    175         */
    176        public NelderMeadSimplex(final double[][] referenceSimplex,
    177                                 final double rho, final double khi,
    178                                 final double gamma, final double sigma) {
    179            super(referenceSimplex);
    180
    181            this.rho = rho;
    182            this.khi = khi;
    183            this.gamma = gamma;
    184            this.sigma = sigma;
    185        }
    186
    187        /** {@inheritDoc} */
    188        @Override
    189        public void iterate(final MultivariateFunction evaluationFunction,
    190                            final Comparator<PointValuePair> comparator) {
    191            // The simplex has n + 1 points if dimension is n.
    192            final int n = getDimension();
    193
    194            // Interesting values.
    195            final PointValuePair best = getPoint(0);
    196            final PointValuePair secondBest = getPoint(n - 1);
    197            final PointValuePair worst = getPoint(n);
    198            final double[] xWorst = worst.getPointRef();
    199
    200            // Compute the centroid of the best vertices (dismissing the worst
    201            // point at index n).
    202            final double[] centroid = new double[n];
    203            for (int i = 0; i < n; i++) {
    204                final double[] x = getPoint(i).getPointRef();
    205                for (int j = 0; j < n; j++) {
    206                    centroid[j] += x[j];
    207                }
    208            }
    209            final double scaling = 1.0 / n;
    210            for (int j = 0; j < n; j++) {
    211                centroid[j] *= scaling;
    212            }
    213
    214            // compute the reflection point
    215            final double[] xR = new double[n];
    216            for (int j = 0; j < n; j++) {
    217                xR[j] = centroid[j] + rho * (centroid[j] - xWorst[j]);
    218            }
    219            final PointValuePair reflected
    220                = new PointValuePair(xR, evaluationFunction.value(xR), false);
    221
    222            if (comparator.compare(best, reflected) <= 0 &&
    223                comparator.compare(reflected, secondBest) < 0) {
    224                // Accept the reflected point.
    225                replaceWorstPoint(reflected, comparator);
    226            } else if (comparator.compare(reflected, best) < 0) {
    227                // Compute the expansion point.
    228                final double[] xE = new double[n];
    229                for (int j = 0; j < n; j++) {
    230                    xE[j] = centroid[j] + khi * (xR[j] - centroid[j]);
    231                }
    232                final PointValuePair expanded
    233                    = new PointValuePair(xE, evaluationFunction.value(xE), false);
    234
    235                if (comparator.compare(expanded, reflected) < 0) {
    236                    // Accept the expansion point.
    237                    replaceWorstPoint(expanded, comparator);
    238                } else {
    239                    // Accept the reflected point.
    240                    replaceWorstPoint(reflected, comparator);
    241                }
    242            } else {
    243                if (comparator.compare(reflected, worst) < 0) {
    244                    // Perform an outside contraction.
    245                    final double[] xC = new double[n];
    246                    for (int j = 0; j < n; j++) {
    247                        xC[j] = centroid[j] + gamma * (xR[j] - centroid[j]);
    248                    }
    249                    final PointValuePair outContracted
    250                        = new PointValuePair(xC, evaluationFunction.value(xC), false);
    251                    if (comparator.compare(outContracted, reflected) <= 0) {
    252                        // Accept the contraction point.
    253                        replaceWorstPoint(outContracted, comparator);
    254                        return;
    255                    }
    256                } else {
    257                    // Perform an inside contraction.
    258                    final double[] xC = new double[n];
    259                    for (int j = 0; j < n; j++) {
    260                        xC[j] = centroid[j] - gamma * (centroid[j] - xWorst[j]);
    261                    }
    262                    final PointValuePair inContracted
    263                        = new PointValuePair(xC, evaluationFunction.value(xC), false);
    264
    265                    if (comparator.compare(inContracted, worst) < 0) {
    266                        // Accept the contraction point.
    267                        replaceWorstPoint(inContracted, comparator);
    268                        return;
    269                    }
    270                }
    271
    272                // Perform a shrink.
    273                final double[] xSmallest = getPoint(0).getPointRef();
    274                for (int i = 1; i <= n; i++) {
    275                    final double[] x = getPoint(i).getPoint();
    276                    for (int j = 0; j < n; j++) {
    277                        x[j] = xSmallest[j] + sigma * (x[j] - xSmallest[j]);
    278                    }
    279                    setPoint(i, new PointValuePair(x, Double.NaN, false));
    280                }
    281                evaluate(evaluationFunction, comparator);
    282            }
    283        }
    284    }