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 package org.apache.commons.math3.optimization.univariate;
018
019 import org.apache.commons.math3.util.Incrementor;
020 import org.apache.commons.math3.exception.NotStrictlyPositiveException;
021 import org.apache.commons.math3.exception.TooManyEvaluationsException;
022 import org.apache.commons.math3.exception.MaxCountExceededException;
023 import org.apache.commons.math3.analysis.UnivariateFunction;
024 import org.apache.commons.math3.optimization.GoalType;
025
026 /**
027 * Provide an interval that brackets a local optimum of a function.
028 * This code is based on a Python implementation (from <em>SciPy</em>,
029 * module {@code optimize.py} v0.5).
030 *
031 * @version $Id: BracketFinder.java 1422230 2012-12-15 12:11:13Z erans $
032 * @deprecated As of 3.1 (to be removed in 4.0).
033 * @since 2.2
034 */
035 @Deprecated
036 public class BracketFinder {
037 /** Tolerance to avoid division by zero. */
038 private static final double EPS_MIN = 1e-21;
039 /**
040 * Golden section.
041 */
042 private static final double GOLD = 1.618034;
043 /**
044 * Factor for expanding the interval.
045 */
046 private final double growLimit;
047 /**
048 * Counter for function evaluations.
049 */
050 private final Incrementor evaluations = new Incrementor();
051 /**
052 * Lower bound of the bracket.
053 */
054 private double lo;
055 /**
056 * Higher bound of the bracket.
057 */
058 private double hi;
059 /**
060 * Point inside the bracket.
061 */
062 private double mid;
063 /**
064 * Function value at {@link #lo}.
065 */
066 private double fLo;
067 /**
068 * Function value at {@link #hi}.
069 */
070 private double fHi;
071 /**
072 * Function value at {@link #mid}.
073 */
074 private double fMid;
075
076 /**
077 * Constructor with default values {@code 100, 50} (see the
078 * {@link #BracketFinder(double,int) other constructor}).
079 */
080 public BracketFinder() {
081 this(100, 50);
082 }
083
084 /**
085 * Create a bracketing interval finder.
086 *
087 * @param growLimit Expanding factor.
088 * @param maxEvaluations Maximum number of evaluations allowed for finding
089 * a bracketing interval.
090 */
091 public BracketFinder(double growLimit,
092 int maxEvaluations) {
093 if (growLimit <= 0) {
094 throw new NotStrictlyPositiveException(growLimit);
095 }
096 if (maxEvaluations <= 0) {
097 throw new NotStrictlyPositiveException(maxEvaluations);
098 }
099
100 this.growLimit = growLimit;
101 evaluations.setMaximalCount(maxEvaluations);
102 }
103
104 /**
105 * Search new points that bracket a local optimum of the function.
106 *
107 * @param func Function whose optimum should be bracketed.
108 * @param goal {@link GoalType Goal type}.
109 * @param xA Initial point.
110 * @param xB Initial point.
111 * @throws TooManyEvaluationsException if the maximum number of evaluations
112 * is exceeded.
113 */
114 public void search(UnivariateFunction func, GoalType goal, double xA, double xB) {
115 evaluations.resetCount();
116 final boolean isMinim = goal == GoalType.MINIMIZE;
117
118 double fA = eval(func, xA);
119 double fB = eval(func, xB);
120 if (isMinim ?
121 fA < fB :
122 fA > fB) {
123
124 double tmp = xA;
125 xA = xB;
126 xB = tmp;
127
128 tmp = fA;
129 fA = fB;
130 fB = tmp;
131 }
132
133 double xC = xB + GOLD * (xB - xA);
134 double fC = eval(func, xC);
135
136 while (isMinim ? fC < fB : fC > fB) {
137 double tmp1 = (xB - xA) * (fB - fC);
138 double tmp2 = (xB - xC) * (fB - fA);
139
140 double val = tmp2 - tmp1;
141 double denom = Math.abs(val) < EPS_MIN ? 2 * EPS_MIN : 2 * val;
142
143 double w = xB - ((xB - xC) * tmp2 - (xB - xA) * tmp1) / denom;
144 double wLim = xB + growLimit * (xC - xB);
145
146 double fW;
147 if ((w - xC) * (xB - w) > 0) {
148 fW = eval(func, w);
149 if (isMinim ?
150 fW < fC :
151 fW > fC) {
152 xA = xB;
153 xB = w;
154 fA = fB;
155 fB = fW;
156 break;
157 } else if (isMinim ?
158 fW > fB :
159 fW < fB) {
160 xC = w;
161 fC = fW;
162 break;
163 }
164 w = xC + GOLD * (xC - xB);
165 fW = eval(func, w);
166 } else if ((w - wLim) * (wLim - xC) >= 0) {
167 w = wLim;
168 fW = eval(func, w);
169 } else if ((w - wLim) * (xC - w) > 0) {
170 fW = eval(func, w);
171 if (isMinim ?
172 fW < fC :
173 fW > fC) {
174 xB = xC;
175 xC = w;
176 w = xC + GOLD * (xC - xB);
177 fB = fC;
178 fC =fW;
179 fW = eval(func, w);
180 }
181 } else {
182 w = xC + GOLD * (xC - xB);
183 fW = eval(func, w);
184 }
185
186 xA = xB;
187 fA = fB;
188 xB = xC;
189 fB = fC;
190 xC = w;
191 fC = fW;
192 }
193
194 lo = xA;
195 fLo = fA;
196 mid = xB;
197 fMid = fB;
198 hi = xC;
199 fHi = fC;
200
201 if (lo > hi) {
202 double tmp = lo;
203 lo = hi;
204 hi = tmp;
205
206 tmp = fLo;
207 fLo = fHi;
208 fHi = tmp;
209 }
210 }
211
212 /**
213 * @return the number of evalutations.
214 */
215 public int getMaxEvaluations() {
216 return evaluations.getMaximalCount();
217 }
218
219 /**
220 * @return the number of evalutations.
221 */
222 public int getEvaluations() {
223 return evaluations.getCount();
224 }
225
226 /**
227 * @return the lower bound of the bracket.
228 * @see #getFLo()
229 */
230 public double getLo() {
231 return lo;
232 }
233
234 /**
235 * Get function value at {@link #getLo()}.
236 * @return function value at {@link #getLo()}
237 */
238 public double getFLo() {
239 return fLo;
240 }
241
242 /**
243 * @return the higher bound of the bracket.
244 * @see #getFHi()
245 */
246 public double getHi() {
247 return hi;
248 }
249
250 /**
251 * Get function value at {@link #getHi()}.
252 * @return function value at {@link #getHi()}
253 */
254 public double getFHi() {
255 return fHi;
256 }
257
258 /**
259 * @return a point in the middle of the bracket.
260 * @see #getFMid()
261 */
262 public double getMid() {
263 return mid;
264 }
265
266 /**
267 * Get function value at {@link #getMid()}.
268 * @return function value at {@link #getMid()}
269 */
270 public double getFMid() {
271 return fMid;
272 }
273
274 /**
275 * @param f Function.
276 * @param x Argument.
277 * @return {@code f(x)}
278 * @throws TooManyEvaluationsException if the maximal number of evaluations is
279 * exceeded.
280 */
281 private double eval(UnivariateFunction f, double x) {
282 try {
283 evaluations.incrementCount();
284 } catch (MaxCountExceededException e) {
285 throw new TooManyEvaluationsException(e.getMax());
286 }
287 return f.value(x);
288 }
289 }