source: src/main/java/agents/anac/y2019/harddealer/math3/optim/univariate/BracketFinder.java

Last change on this file was 204, checked in by Katsuhide Fujita, 5 years ago

Fixed errors of ANAC2019 agents
  • Property svn:executable set to *
File size: 7.9 KB
Line 
1/*
2 * Licensed to the Apache Software Foundation (ASF) under one or more
3 * contributor license agreements. See the NOTICE file distributed with
4 * this work for additional information regarding copyright ownership.
5 * The ASF licenses this file to You under the Apache License, Version 2.0
6 * (the "License"); you may not use this file except in compliance with
7 * the License. You may obtain a copy of the License at
8 *
9 * http://www.apache.org/licenses/LICENSE-2.0
10 *
11 * Unless required by applicable law or agreed to in writing, software
12 * distributed under the License is distributed on an "AS IS" BASIS,
13 * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
14 * See the License for the specific language governing permissions and
15 * limitations under the License.
16 */
17package agents.anac.y2019.harddealer.math3.optim.univariate;
18
19import agents.anac.y2019.harddealer.math3.util.FastMath;
20import agents.anac.y2019.harddealer.math3.util.IntegerSequence;
21import agents.anac.y2019.harddealer.math3.exception.NotStrictlyPositiveException;
22import agents.anac.y2019.harddealer.math3.exception.TooManyEvaluationsException;
23import agents.anac.y2019.harddealer.math3.exception.MaxCountExceededException;
24import agents.anac.y2019.harddealer.math3.analysis.UnivariateFunction;
25import agents.anac.y2019.harddealer.math3.optim.nonlinear.scalar.GoalType;
26
27/**
28 * Provide an interval that brackets a local optimum of a function.
29 * This code is based on a Python implementation (from <em>SciPy</em>,
30 * module {@code optimize.py} v0.5).
31 *
32 * @since 2.2
33 */
34public class BracketFinder {
35 /** Tolerance to avoid division by zero. */
36 private static final double EPS_MIN = 1e-21;
37 /**
38 * Golden section.
39 */
40 private static final double GOLD = 1.618034;
41 /**
42 * Factor for expanding the interval.
43 */
44 private final double growLimit;
45 /**
46 * Counter for function evaluations.
47 */
48 private IntegerSequence.Incrementor evaluations;
49 /**
50 * Lower bound of the bracket.
51 */
52 private double lo;
53 /**
54 * Higher bound of the bracket.
55 */
56 private double hi;
57 /**
58 * Point inside the bracket.
59 */
60 private double mid;
61 /**
62 * Function value at {@link #lo}.
63 */
64 private double fLo;
65 /**
66 * Function value at {@link #hi}.
67 */
68 private double fHi;
69 /**
70 * Function value at {@link #mid}.
71 */
72 private double fMid;
73
74 /**
75 * Constructor with default values {@code 100, 500} (see the
76 * {@link #BracketFinder(double,int) other constructor}).
77 */
78 public BracketFinder() {
79 this(100, 500);
80 }
81
82 /**
83 * Create a bracketing interval finder.
84 *
85 * @param growLimit Expanding factor.
86 * @param maxEvaluations Maximum number of evaluations allowed for finding
87 * a bracketing interval.
88 */
89 public BracketFinder(double growLimit,
90 int maxEvaluations) {
91 if (growLimit <= 0) {
92 throw new NotStrictlyPositiveException(growLimit);
93 }
94 if (maxEvaluations <= 0) {
95 throw new NotStrictlyPositiveException(maxEvaluations);
96 }
97
98 this.growLimit = growLimit;
99 evaluations = IntegerSequence.Incrementor.create().withMaximalCount(maxEvaluations);
100 }
101
102 /**
103 * Search new points that bracket a local optimum of the function.
104 *
105 * @param func Function whose optimum should be bracketed.
106 * @param goal {@link GoalType Goal type}.
107 * @param xA Initial point.
108 * @param xB Initial point.
109 * @throws TooManyEvaluationsException if the maximum number of evaluations
110 * is exceeded.
111 */
112 public void search(UnivariateFunction func,
113 GoalType goal,
114 double xA,
115 double xB) {
116 evaluations = evaluations.withStart(0);
117 final boolean isMinim = goal == GoalType.MINIMIZE;
118
119 double fA = eval(func, xA);
120 double fB = eval(func, xB);
121 if (isMinim ?
122 fA < fB :
123 fA > fB) {
124
125 double tmp = xA;
126 xA = xB;
127 xB = tmp;
128
129 tmp = fA;
130 fA = fB;
131 fB = tmp;
132 }
133
134 double xC = xB + GOLD * (xB - xA);
135 double fC = eval(func, xC);
136
137 while (isMinim ? fC < fB : fC > fB) {
138 double tmp1 = (xB - xA) * (fB - fC);
139 double tmp2 = (xB - xC) * (fB - fA);
140
141 double val = tmp2 - tmp1;
142 double denom = FastMath.abs(val) < EPS_MIN ? 2 * EPS_MIN : 2 * val;
143
144 double w = xB - ((xB - xC) * tmp2 - (xB - xA) * tmp1) / denom;
145 double wLim = xB + growLimit * (xC - xB);
146
147 double fW;
148 if ((w - xC) * (xB - w) > 0) {
149 fW = eval(func, w);
150 if (isMinim ?
151 fW < fC :
152 fW > fC) {
153 xA = xB;
154 xB = w;
155 fA = fB;
156 fB = fW;
157 break;
158 } else if (isMinim ?
159 fW > fB :
160 fW < fB) {
161 xC = w;
162 fC = fW;
163 break;
164 }
165 w = xC + GOLD * (xC - xB);
166 fW = eval(func, w);
167 } else if ((w - wLim) * (wLim - xC) >= 0) {
168 w = wLim;
169 fW = eval(func, w);
170 } else if ((w - wLim) * (xC - w) > 0) {
171 fW = eval(func, w);
172 if (isMinim ?
173 fW < fC :
174 fW > fC) {
175 xB = xC;
176 xC = w;
177 w = xC + GOLD * (xC - xB);
178 fB = fC;
179 fC =fW;
180 fW = eval(func, w);
181 }
182 } else {
183 w = xC + GOLD * (xC - xB);
184 fW = eval(func, w);
185 }
186
187 xA = xB;
188 fA = fB;
189 xB = xC;
190 fB = fC;
191 xC = w;
192 fC = fW;
193 }
194
195 lo = xA;
196 fLo = fA;
197 mid = xB;
198 fMid = fB;
199 hi = xC;
200 fHi = fC;
201
202 if (lo > hi) {
203 double tmp = lo;
204 lo = hi;
205 hi = tmp;
206
207 tmp = fLo;
208 fLo = fHi;
209 fHi = tmp;
210 }
211 }
212
213 /**
214 * @return the number of evalutations.
215 */
216 public int getMaxEvaluations() {
217 return evaluations.getMaximalCount();
218 }
219
220 /**
221 * @return the number of evalutations.
222 */
223 public int getEvaluations() {
224 return evaluations.getCount();
225 }
226
227 /**
228 * @return the lower bound of the bracket.
229 * @see #getFLo()
230 */
231 public double getLo() {
232 return lo;
233 }
234
235 /**
236 * Get function value at {@link #getLo()}.
237 * @return function value at {@link #getLo()}
238 */
239 public double getFLo() {
240 return fLo;
241 }
242
243 /**
244 * @return the higher bound of the bracket.
245 * @see #getFHi()
246 */
247 public double getHi() {
248 return hi;
249 }
250
251 /**
252 * Get function value at {@link #getHi()}.
253 * @return function value at {@link #getHi()}
254 */
255 public double getFHi() {
256 return fHi;
257 }
258
259 /**
260 * @return a point in the middle of the bracket.
261 * @see #getFMid()
262 */
263 public double getMid() {
264 return mid;
265 }
266
267 /**
268 * Get function value at {@link #getMid()}.
269 * @return function value at {@link #getMid()}
270 */
271 public double getFMid() {
272 return fMid;
273 }
274
275 /**
276 * @param f Function.
277 * @param x Argument.
278 * @return {@code f(x)}
279 * @throws TooManyEvaluationsException if the maximal number of evaluations is
280 * exceeded.
281 */
282 private double eval(UnivariateFunction f, double x) {
283 try {
284 evaluations.increment();
285 } catch (MaxCountExceededException e) {
286 throw new TooManyEvaluationsException(e.getMax());
287 }
288 return f.value(x);
289 }
290}
Note: See TracBrowser for help on using the repository browser.