1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23 package org.hipparchus.optim.nonlinear.scalar.gradient;
24
25 import org.hipparchus.analysis.MultivariateFunction;
26 import org.hipparchus.analysis.MultivariateVectorFunction;
27 import org.hipparchus.exception.MathRuntimeException;
28 import org.hipparchus.geometry.euclidean.twod.Vector2D;
29 import org.hipparchus.linear.BlockRealMatrix;
30 import org.hipparchus.linear.RealMatrix;
31 import org.hipparchus.optim.InitialGuess;
32 import org.hipparchus.optim.MaxEval;
33 import org.hipparchus.optim.PointValuePair;
34 import org.hipparchus.optim.SimpleBounds;
35 import org.hipparchus.optim.SimpleValueChecker;
36 import org.hipparchus.optim.nonlinear.scalar.GoalType;
37 import org.hipparchus.optim.nonlinear.scalar.ObjectiveFunction;
38 import org.hipparchus.optim.nonlinear.scalar.ObjectiveFunctionGradient;
39 import org.junit.jupiter.api.Test;
40
41 import static org.junit.jupiter.api.Assertions.assertEquals;
42 import static org.junit.jupiter.api.Assertions.assertThrows;
43 import static org.junit.jupiter.api.Assertions.assertTrue;
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107 class NonLinearConjugateGradientOptimizerTest {
108 @Test
109 void testBoundsUnsupported() {
110 assertThrows(MathRuntimeException.class, () -> {
111 LinearProblem problem
112 = new LinearProblem(new double[][]{{2}}, new double[]{3});
113 NonLinearConjugateGradientOptimizer optimizer
114 = new NonLinearConjugateGradientOptimizer(NonLinearConjugateGradientOptimizer.Formula.POLAK_RIBIERE,
115 new SimpleValueChecker(1e-6, 1e-6),
116 1e-3, 1e-3, 1);
117 optimizer.optimize(new MaxEval(100),
118 problem.getObjectiveFunction(),
119 problem.getObjectiveFunctionGradient(),
120 GoalType.MINIMIZE,
121 new InitialGuess(new double[]{0}),
122 new SimpleBounds(new double[]{-1},
123 new double[]{1}));
124 });
125 }
126
127 @Test
128 void testTrivial() {
129 LinearProblem problem
130 = new LinearProblem(new double[][] { { 2 } }, new double[] { 3 });
131 NonLinearConjugateGradientOptimizer optimizer
132 = new NonLinearConjugateGradientOptimizer(NonLinearConjugateGradientOptimizer.Formula.POLAK_RIBIERE,
133 new SimpleValueChecker(1e-6, 1e-6),
134 1e-3, 1e-3, 1);
135 PointValuePair optimum
136 = optimizer.optimize(new MaxEval(100),
137 problem.getObjectiveFunction(),
138 problem.getObjectiveFunctionGradient(),
139 GoalType.MINIMIZE,
140 new InitialGuess(new double[] { 0 }));
141 assertEquals(1.5, optimum.getPoint()[0], 1.0e-10);
142 assertEquals(0.0, optimum.getValue(), 1.0e-10);
143
144
145 assertTrue(optimizer.getIterations() > 0);
146 }
147
148 @Test
149 void testColumnsPermutation() {
150 LinearProblem problem
151 = new LinearProblem(new double[][] { { 1.0, -1.0 }, { 0.0, 2.0 }, { 1.0, -2.0 } },
152 new double[] { 4.0, 6.0, 1.0 });
153
154 NonLinearConjugateGradientOptimizer optimizer
155 = new NonLinearConjugateGradientOptimizer(NonLinearConjugateGradientOptimizer.Formula.POLAK_RIBIERE,
156 new SimpleValueChecker(1e-6, 1e-6),
157 1e-3, 1e-3, 1);
158 PointValuePair optimum
159 = optimizer.optimize(new MaxEval(100),
160 problem.getObjectiveFunction(),
161 problem.getObjectiveFunctionGradient(),
162 GoalType.MINIMIZE,
163 new InitialGuess(new double[] { 0, 0 }));
164 assertEquals(7.0, optimum.getPoint()[0], 1.0e-10);
165 assertEquals(3.0, optimum.getPoint()[1], 1.0e-10);
166 assertEquals(0.0, optimum.getValue(), 1.0e-10);
167
168 }
169
170 @Test
171 void testNoDependency() {
172 LinearProblem problem = new LinearProblem(new double[][] {
173 { 2, 0, 0, 0, 0, 0 },
174 { 0, 2, 0, 0, 0, 0 },
175 { 0, 0, 2, 0, 0, 0 },
176 { 0, 0, 0, 2, 0, 0 },
177 { 0, 0, 0, 0, 2, 0 },
178 { 0, 0, 0, 0, 0, 2 }
179 }, new double[] { 0.0, 1.1, 2.2, 3.3, 4.4, 5.5 });
180 NonLinearConjugateGradientOptimizer optimizer
181 = new NonLinearConjugateGradientOptimizer(NonLinearConjugateGradientOptimizer.Formula.POLAK_RIBIERE,
182 new SimpleValueChecker(1e-6, 1e-6),
183 1e-3, 1e-3, 1);
184 PointValuePair optimum
185 = optimizer.optimize(new MaxEval(100),
186 problem.getObjectiveFunction(),
187 problem.getObjectiveFunctionGradient(),
188 GoalType.MINIMIZE,
189 new InitialGuess(new double[] { 0, 0, 0, 0, 0, 0 }));
190 for (int i = 0; i < problem.target.length; ++i) {
191 assertEquals(0.55 * i, optimum.getPoint()[i], 1.0e-10);
192 }
193 }
194
195 @Test
196 void testOneSet() {
197 LinearProblem problem = new LinearProblem(new double[][] {
198 { 1, 0, 0 },
199 { -1, 1, 0 },
200 { 0, -1, 1 }
201 }, new double[] { 1, 1, 1});
202 NonLinearConjugateGradientOptimizer optimizer
203 = new NonLinearConjugateGradientOptimizer(NonLinearConjugateGradientOptimizer.Formula.POLAK_RIBIERE,
204 new SimpleValueChecker(1e-6, 1e-6),
205 1e-3, 1e-3, 1);
206 PointValuePair optimum
207 = optimizer.optimize(new MaxEval(100),
208 problem.getObjectiveFunction(),
209 problem.getObjectiveFunctionGradient(),
210 GoalType.MINIMIZE,
211 new InitialGuess(new double[] { 0, 0, 0 }));
212 assertEquals(1.0, optimum.getPoint()[0], 1.0e-10);
213 assertEquals(2.0, optimum.getPoint()[1], 1.0e-10);
214 assertEquals(3.0, optimum.getPoint()[2], 1.0e-10);
215
216 }
217
218 @Test
219 void testTwoSets() {
220 final double epsilon = 1.0e-7;
221 LinearProblem problem = new LinearProblem(new double[][] {
222 { 2, 1, 0, 4, 0, 0 },
223 { -4, -2, 3, -7, 0, 0 },
224 { 4, 1, -2, 8, 0, 0 },
225 { 0, -3, -12, -1, 0, 0 },
226 { 0, 0, 0, 0, epsilon, 1 },
227 { 0, 0, 0, 0, 1, 1 }
228 }, new double[] { 2, -9, 2, 2, 1 + epsilon * epsilon, 2});
229
230 final Preconditioner preconditioner
231 = new Preconditioner() {
232 public double[] precondition(double[] point, double[] r) {
233 double[] d = r.clone();
234 d[0] /= 72.0;
235 d[1] /= 30.0;
236 d[2] /= 314.0;
237 d[3] /= 260.0;
238 d[4] /= 2 * (1 + epsilon * epsilon);
239 d[5] /= 4.0;
240 return d;
241 }
242 };
243
244 NonLinearConjugateGradientOptimizer optimizer
245 = new NonLinearConjugateGradientOptimizer(NonLinearConjugateGradientOptimizer.Formula.POLAK_RIBIERE,
246 new SimpleValueChecker(1e-13, 1e-13),
247 1e-7, 1e-7, 1,
248 preconditioner);
249
250 PointValuePair optimum
251 = optimizer.optimize(new MaxEval(100),
252 problem.getObjectiveFunction(),
253 problem.getObjectiveFunctionGradient(),
254 GoalType.MINIMIZE,
255 new InitialGuess(new double[] { 0, 0, 0, 0, 0, 0 }));
256
257 final double[] result = optimum.getPoint();
258 final double[] expected = {3, 4, -1, -2, 1 + epsilon, 1 - epsilon};
259
260 assertEquals(expected[0], result[0], 1.0e-7);
261 assertEquals(expected[1], result[1], 1.0e-7);
262 assertEquals(expected[2], result[2], 1.0e-9);
263 assertEquals(expected[3], result[3], 1.0e-8);
264 assertEquals(expected[4] + epsilon, result[4], 1.0e-6);
265 assertEquals(expected[5] - epsilon, result[5], 1.0e-6);
266
267 }
268
269 @Test
270 void testNonInversible() {
271 LinearProblem problem = new LinearProblem(new double[][] {
272 { 1, 2, -3 },
273 { 2, 1, 3 },
274 { -3, 0, -9 }
275 }, new double[] { 1, 1, 1 });
276 NonLinearConjugateGradientOptimizer optimizer
277 = new NonLinearConjugateGradientOptimizer(NonLinearConjugateGradientOptimizer.Formula.POLAK_RIBIERE,
278 new SimpleValueChecker(1e-6, 1e-6),
279 1e-3, 1e-3, 1);
280 PointValuePair optimum
281 = optimizer.optimize(new MaxEval(100),
282 problem.getObjectiveFunction(),
283 problem.getObjectiveFunctionGradient(),
284 GoalType.MINIMIZE,
285 new InitialGuess(new double[] { 0, 0, 0 }));
286 assertTrue(optimum.getValue() > 0.5);
287 }
288
289 @Test
290 void testIllConditioned() {
291 LinearProblem problem1 = new LinearProblem(new double[][] {
292 { 10.0, 7.0, 8.0, 7.0 },
293 { 7.0, 5.0, 6.0, 5.0 },
294 { 8.0, 6.0, 10.0, 9.0 },
295 { 7.0, 5.0, 9.0, 10.0 }
296 }, new double[] { 32, 23, 33, 31 });
297 NonLinearConjugateGradientOptimizer optimizer
298 = new NonLinearConjugateGradientOptimizer(NonLinearConjugateGradientOptimizer.Formula.POLAK_RIBIERE,
299 new SimpleValueChecker(1e-13, 1e-13),
300 1e-15, 1e-15, 1);
301 PointValuePair optimum1
302 = optimizer.optimize(new MaxEval(200),
303 problem1.getObjectiveFunction(),
304 problem1.getObjectiveFunctionGradient(),
305 GoalType.MINIMIZE,
306 new InitialGuess(new double[] { 0, 1, 2, 3 }));
307 assertEquals(1.0, optimum1.getPoint()[0], 1.0e-4);
308 assertEquals(1.0, optimum1.getPoint()[1], 1.0e-3);
309 assertEquals(1.0, optimum1.getPoint()[2], 1.0e-4);
310 assertEquals(1.0, optimum1.getPoint()[3], 1.0e-4);
311
312 LinearProblem problem2 = new LinearProblem(new double[][] {
313 { 10.00, 7.00, 8.10, 7.20 },
314 { 7.08, 5.04, 6.00, 5.00 },
315 { 8.00, 5.98, 9.89, 9.00 },
316 { 6.99, 4.99, 9.00, 9.98 }
317 }, new double[] { 32, 23, 33, 31 });
318 PointValuePair optimum2
319 = optimizer.optimize(new MaxEval(200),
320 problem2.getObjectiveFunction(),
321 problem2.getObjectiveFunctionGradient(),
322 GoalType.MINIMIZE,
323 new InitialGuess(new double[] { 0, 1, 2, 3 }));
324
325 final double[] result2 = optimum2.getPoint();
326 final double[] expected2 = {-81, 137, -34, 22};
327
328 assertEquals(expected2[0], result2[0], 2);
329 assertEquals(expected2[1], result2[1], 4);
330 assertEquals(expected2[2], result2[2], 1);
331 assertEquals(expected2[3], result2[3], 1);
332 }
333
334 @Test
335 void testMoreEstimatedParametersSimple() {
336 LinearProblem problem = new LinearProblem(new double[][] {
337 { 3.0, 2.0, 0.0, 0.0 },
338 { 0.0, 1.0, -1.0, 1.0 },
339 { 2.0, 0.0, 1.0, 0.0 }
340 }, new double[] { 7.0, 3.0, 5.0 });
341
342 NonLinearConjugateGradientOptimizer optimizer
343 = new NonLinearConjugateGradientOptimizer(NonLinearConjugateGradientOptimizer.Formula.POLAK_RIBIERE,
344 new SimpleValueChecker(1e-6, 1e-6),
345 1e-3, 1e-3, 1);
346 PointValuePair optimum
347 = optimizer.optimize(new MaxEval(100),
348 problem.getObjectiveFunction(),
349 problem.getObjectiveFunctionGradient(),
350 GoalType.MINIMIZE,
351 new InitialGuess(new double[] { 7, 6, 5, 4 }));
352 assertEquals(0, optimum.getValue(), 1.0e-10);
353
354 }
355
356 @Test
357 void testMoreEstimatedParametersUnsorted() {
358 LinearProblem problem = new LinearProblem(new double[][] {
359 { 1.0, 1.0, 0.0, 0.0, 0.0, 0.0 },
360 { 0.0, 0.0, 1.0, 1.0, 1.0, 0.0 },
361 { 0.0, 0.0, 0.0, 0.0, 1.0, -1.0 },
362 { 0.0, 0.0, -1.0, 1.0, 0.0, 1.0 },
363 { 0.0, 0.0, 0.0, -1.0, 1.0, 0.0 }
364 }, new double[] { 3.0, 12.0, -1.0, 7.0, 1.0 });
365 NonLinearConjugateGradientOptimizer optimizer
366 = new NonLinearConjugateGradientOptimizer(NonLinearConjugateGradientOptimizer.Formula.POLAK_RIBIERE,
367 new SimpleValueChecker(1e-6, 1e-6),
368 1e-3, 1e-3, 1);
369 PointValuePair optimum
370 = optimizer.optimize(new MaxEval(100),
371 problem.getObjectiveFunction(),
372 problem.getObjectiveFunctionGradient(),
373 GoalType.MINIMIZE,
374 new InitialGuess(new double[] { 2, 2, 2, 2, 2, 2 }));
375 assertEquals(0, optimum.getValue(), 1.0e-10);
376 }
377
378 @Test
379 void testRedundantEquations() {
380 LinearProblem problem = new LinearProblem(new double[][] {
381 { 1.0, 1.0 },
382 { 1.0, -1.0 },
383 { 1.0, 3.0 }
384 }, new double[] { 3.0, 1.0, 5.0 });
385
386 NonLinearConjugateGradientOptimizer optimizer
387 = new NonLinearConjugateGradientOptimizer(NonLinearConjugateGradientOptimizer.Formula.POLAK_RIBIERE,
388 new SimpleValueChecker(1e-6, 1e-6),
389 1e-3, 1e-3, 1);
390 PointValuePair optimum
391 = optimizer.optimize(new MaxEval(100),
392 problem.getObjectiveFunction(),
393 problem.getObjectiveFunctionGradient(),
394 GoalType.MINIMIZE,
395 new InitialGuess(new double[] { 1, 1 }));
396 assertEquals(2.0, optimum.getPoint()[0], 1.0e-8);
397 assertEquals(1.0, optimum.getPoint()[1], 1.0e-8);
398
399 }
400
401 @Test
402 void testInconsistentEquations() {
403 LinearProblem problem = new LinearProblem(new double[][] {
404 { 1.0, 1.0 },
405 { 1.0, -1.0 },
406 { 1.0, 3.0 }
407 }, new double[] { 3.0, 1.0, 4.0 });
408
409 NonLinearConjugateGradientOptimizer optimizer
410 = new NonLinearConjugateGradientOptimizer(NonLinearConjugateGradientOptimizer.Formula.POLAK_RIBIERE,
411 new SimpleValueChecker(1e-6, 1e-6),
412 1e-3, 1e-3, 1);
413 PointValuePair optimum
414 = optimizer.optimize(new MaxEval(100),
415 problem.getObjectiveFunction(),
416 problem.getObjectiveFunctionGradient(),
417 GoalType.MINIMIZE,
418 new InitialGuess(new double[] { 1, 1 }));
419 assertTrue(optimum.getValue() > 0.1);
420
421 }
422
423 @Test
424 void testCircleFitting() {
425 CircleScalar problem = new CircleScalar();
426 problem.addPoint( 30.0, 68.0);
427 problem.addPoint( 50.0, -6.0);
428 problem.addPoint(110.0, -20.0);
429 problem.addPoint( 35.0, 15.0);
430 problem.addPoint( 45.0, 97.0);
431 NonLinearConjugateGradientOptimizer optimizer
432 = new NonLinearConjugateGradientOptimizer(NonLinearConjugateGradientOptimizer.Formula.POLAK_RIBIERE,
433 new SimpleValueChecker(1e-30, 1e-30),
434 1e-15, 1e-13, 1);
435 PointValuePair optimum
436 = optimizer.optimize(new MaxEval(100),
437 problem.getObjectiveFunction(),
438 problem.getObjectiveFunctionGradient(),
439 GoalType.MINIMIZE,
440 new InitialGuess(new double[] { 98.680, 47.345 }));
441 Vector2D center = new Vector2D(optimum.getPointRef()[0], optimum.getPointRef()[1]);
442 assertEquals(69.960161753, problem.getRadius(center), 1.0e-8);
443 assertEquals(96.075902096, center.getX(), 1.0e-7);
444 assertEquals(48.135167894, center.getY(), 1.0e-6);
445 }
446
447 private static class LinearProblem {
448 final RealMatrix factors;
449 final double[] target;
450
451 public LinearProblem(double[][] factors,
452 double[] target) {
453 this.factors = new BlockRealMatrix(factors);
454 this.target = target;
455 }
456
457 public ObjectiveFunction getObjectiveFunction() {
458 return new ObjectiveFunction(new MultivariateFunction() {
459 public double value(double[] point) {
460 double[] y = factors.operate(point);
461 double sum = 0;
462 for (int i = 0; i < y.length; ++i) {
463 double ri = y[i] - target[i];
464 sum += ri * ri;
465 }
466 return sum;
467 }
468 });
469 }
470
471 public ObjectiveFunctionGradient getObjectiveFunctionGradient() {
472 return new ObjectiveFunctionGradient(new MultivariateVectorFunction() {
473 public double[] value(double[] point) {
474 double[] r = factors.operate(point);
475 for (int i = 0; i < r.length; ++i) {
476 r[i] -= target[i];
477 }
478 double[] p = factors.transpose().operate(r);
479 for (int i = 0; i < p.length; ++i) {
480 p[i] *= 2;
481 }
482 return p;
483 }
484 });
485 }
486 }
487 }