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 */ 017package org.apache.commons.math4.legacy.optim.nonlinear.scalar.noderiv; 018 019import java.util.Comparator; 020import java.util.List; 021import java.util.ArrayList; 022import java.util.Collections; 023import java.util.function.UnaryOperator; 024import java.util.function.DoublePredicate; 025import org.apache.commons.rng.UniformRandomProvider; 026import org.apache.commons.rng.simple.RandomSource; 027import org.apache.commons.rng.sampling.distribution.ContinuousSampler; 028import org.apache.commons.rng.sampling.distribution.ContinuousUniformSampler; 029import org.apache.commons.math4.legacy.analysis.MultivariateFunction; 030import org.apache.commons.math4.legacy.optim.PointValuePair; 031 032/** 033 * DSSA algorithm. 034 * 035 * Described in 036 * <blockquote> 037 * <em>Abdel-Rahman Hedar and Masao Fukushima (2002)</em>, 038 * <b> 039 * Hybrid simulated annealing and direct search method 040 * for nonlinear unconstrained global optimization 041 * </b>, 042 * Optimization Methods and Software, 17:5, 891-912, 043 * DOI: 10.1080/1055678021000030084 044 * </blockquote> 045 * 046 * <p> 047 * A note about the {@link #HedarFukushimaTransform(double) "shrink" factor}: 048 * Per DSSA's description, the simplex must keep its size during the simulated 049 * annealing (SA) phase to avoid premature convergence. This assumes that the 050 * best candidates from the SA phase will each subsequently serve as starting 051 * point for another optimization to hone in on the local optimum. 052 * Values lower than 1 and no subsequent "best list" search correspond to the 053 * "SSA" algorithm in the above paper. 054 */ 055public class HedarFukushimaTransform 056 implements Simplex.TransformFactory { 057 /** Shrinkage coefficient. */ 058 private final double sigma; 059 /** Sampler for reflection coefficient. */ 060 private final ContinuousSampler alphaSampler; 061 /** No shrink indicator. */ 062 private final boolean noShrink; 063 064 /** 065 * @param sigma Shrink factor. 066 * @param rng Random generator. 067 * @throws IllegalArgumentException if {@code sigma <= 0} or 068 * {@code sigma > 1}. 069 */ 070 public HedarFukushimaTransform(double sigma, 071 UniformRandomProvider rng) { 072 if (sigma <= 0 || 073 sigma > 1) { 074 throw new IllegalArgumentException("Shrink factor out of range: " + 075 sigma); 076 } 077 078 this.sigma = sigma; 079 alphaSampler = ContinuousUniformSampler.of(rng, 0.9, 1.1); 080 noShrink = sigma == 1d; 081 } 082 083 /** 084 * @param sigma Shrink factor. 085 * @throws IllegalArgumentException if {@code sigma <= 0} or 086 * {@code sigma > 1}. 087 */ 088 public HedarFukushimaTransform(double sigma) { 089 this(sigma, RandomSource.KISS.create()); 090 } 091 092 /** 093 * Disable shrinking of the simplex (as mandated by DSSA). 094 */ 095 public HedarFukushimaTransform() { 096 this(1d); 097 } 098 099 /** {@inheritDoc} */ 100 @Override 101 public UnaryOperator<Simplex> create(final MultivariateFunction evaluationFunction, 102 final Comparator<PointValuePair> comparator, 103 final DoublePredicate saAcceptance) { 104 if (saAcceptance == null) { 105 throw new IllegalArgumentException("Missing SA acceptance test"); 106 } 107 108 return original -> transform(original, 109 saAcceptance, 110 evaluationFunction, 111 comparator); 112 } 113 114 /** 115 * Simulated annealing step (at fixed temperature). 116 * 117 * @param original Current simplex. Points must be sorted from best to worst. 118 * @param sa Simulated annealing acceptance test. 119 * @param eval Evaluation function. 120 * @param comp Objective function comparator. 121 * @return a new instance. 122 */ 123 private Simplex transform(Simplex original, 124 DoublePredicate sa, 125 MultivariateFunction eval, 126 Comparator<PointValuePair> comp) { 127 final int size = original.getSize(); 128 // Current best point. 129 final PointValuePair best = original.get(0); 130 final double bestValue = best.getValue(); 131 132 for (int k = 1; k < size; k++) { 133 // Perform reflections of the "k" worst points. 134 final List<PointValuePair> reflected = reflectPoints(original, k, eval); 135 Collections.sort(reflected, comp); 136 137 // Check whether the best of the reflected points is better than the 138 // current overall best. 139 final PointValuePair candidate = reflected.get(0); 140 final boolean candidateIsBetter = comp.compare(candidate, best) < 0; 141 final boolean candidateIsAccepted = candidateIsBetter || 142 sa.test(candidate.getValue() - bestValue); 143 144 if (candidateIsAccepted) { 145 // Replace worst points with the reflected points. 146 return original.replaceLast(reflected); 147 } 148 } 149 150 // No direction provided a better point. 151 return noShrink ? 152 original : 153 original.shrink(sigma, eval); 154 } 155 156 /** 157 * @param simplex Current simplex (whose points must be sorted, from best 158 * to worst). 159 * @param nPoints Number of points to reflect. 160 * The {@code nPoints} worst points will be reflected through the centroid 161 * of the {@code n + 1 - nPoints} best points. 162 * @param eval Evaluation function. 163 * @return the (unsorted) list of reflected points. 164 * @throws IllegalArgumentException if {@code nPoints < 1} or 165 * {@code nPoints > n}. 166 */ 167 private List<PointValuePair> reflectPoints(Simplex simplex, 168 int nPoints, 169 MultivariateFunction eval) { 170 final int size = simplex.getSize(); 171 if (nPoints < 1 || 172 nPoints >= size) { 173 throw new IllegalArgumentException("Out of range: " + nPoints); 174 } 175 176 final int nCentroid = size - nPoints; 177 final List<PointValuePair> centroidList = simplex.asList(0, nCentroid); 178 final List<PointValuePair> reflectList = simplex.asList(nCentroid, size); 179 180 final double[] centroid = Simplex.centroid(centroidList); 181 182 final List<PointValuePair> reflected = new ArrayList<>(nPoints); 183 for (int i = 0; i < reflectList.size(); i++) { 184 reflected.add(newReflectedPoint(reflectList.get(i), 185 centroid, 186 eval)); 187 } 188 189 return reflected; 190 } 191 192 /** 193 * @param point Current point. 194 * @param centroid Coordinates through which reflection must be performed. 195 * @param eval Evaluation function. 196 * @return a new point with Cartesian coordinates set to the reflection 197 * of {@code point} through {@code centroid}. 198 */ 199 private PointValuePair newReflectedPoint(PointValuePair point, 200 double[] centroid, 201 MultivariateFunction eval) { 202 final double alpha = alphaSampler.sample(); 203 return Simplex.newPoint(centroid, 204 -alpha, 205 point.getPoint(), 206 eval); 207 } 208 209 /** {@inheritDoc} */ 210 @Override 211 public String toString() { 212 return "Hedar-Fukushima [s=" + sigma + "]"; 213 } 214}