1 #ifndef SimTK_SIMMATH_OPTIMIZER_SYSTEMS_H_
2 #define SimTK_SIMMATH_OPTIMIZER_SYSTEMS_H_
3 
4 /* -------------------------------------------------------------------------- *
5  *                        Simbody(tm): SimTKmath                              *
6  * -------------------------------------------------------------------------- *
7  * This is part of the SimTK biosimulation toolkit originating from           *
8  * Simbios, the NIH National Center for Physics-Based Simulation of           *
9  * Biological Structures at Stanford, funded under the NIH Roadmap for        *
10  * Medical Research, grant U54 GM072970. See https://simtk.org/home/simbody.  *
11  *                                                                            *
12  * Portions copyright (c) 2006-14 Stanford University and the Authors.        *
13  * Authors: Chris Dembia                                                      *
14  * Contributors:                                                              *
15  *                                                                            *
16  * Licensed under the Apache License, Version 2.0 (the "License"); you may    *
17  * not use this file except in compliance with the License. You may obtain a  *
18  * copy of the License at http://www.apache.org/licenses/LICENSE-2.0.         *
19  *                                                                            *
20  * Unless required by applicable law or agreed to in writing, software        *
21  * distributed under the License is distributed on an "AS IS" BASIS,          *
22  * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.   *
23  * See the License for the specific language governing permissions and        *
24  * limitations under the License.                                             *
25  * -------------------------------------------------------------------------- */
26 
27 // These websites list test functions for optimization.
28 // https://en.wikipedia.org/wiki/Test_functions_for_optimization
29 // http://www.sfu.ca/~ssurjano/optimization.html
30 
31 #include "SimTKmath.h"
32 
33 using SimTK::Vector;
34 using SimTK::Real;
35 using SimTK::OptimizerSystem;
36 using SimTK::Pi;
37 using SimTK::square;
38 using SimTK::sqrt;
39 
40 class TestOptimizerSystem : public OptimizerSystem {
41 public:
TestOptimizerSystem(int nParameters)42     TestOptimizerSystem(int nParameters)
43         : OptimizerSystem(nParameters) {}
optimalValue()44     virtual Real optimalValue() const {
45         Real f;
46         objectiveFunc(optimalParameters(), true, f);
47         return f;
48     }
49     virtual Vector optimalParameters() const = 0;
50 };
51 
52 // This function come from Nikolaus Hansen's source code
53 // (https://github.com/cma-es). I think it is supposed to be cigar-shaped.
54 class Cigtab : public TestOptimizerSystem {
55 public:
Cigtab(int nParameters)56     Cigtab(int nParameters) : TestOptimizerSystem(nParameters) {}
objectiveFunc(const Vector & x,bool new_parameters,Real & f)57     int objectiveFunc(const Vector& x, bool new_parameters,
58             Real& f) const override {
59         f = 1e4 * x[0] * x[0] + 1e-4 * x[1] * x[1];
60         for (int i = 0; i < getNumParameters(); ++i) {
61             f += x[i] * x[i];
62         }
63         return 0;
64     }
optimalParameters()65     Vector optimalParameters() const override {
66         Vector x(getNumParameters());
67         x.setToZero();
68         return x;
69     }
70 };
71 
72 // A function with many local minima. http://www.sfu.ca/~ssurjano/ackley.html
73 class Ackley : public TestOptimizerSystem {
74 public:
Ackley(int nParameters)75     Ackley(int nParameters) : TestOptimizerSystem(nParameters),
76             a(20), b(0.2), c(2 * Pi) {
77         // The website above says this function usually has the following
78         // bounds:
79         Vector limits(nParameters);
80         limits.setTo(32.768);
81         setParameterLimits(-limits, limits);
82     }
objectiveFunc(const Vector & x,bool new_parameters,Real & f)83     int objectiveFunc(const Vector& x, bool new_parameters, Real& f) const override {
84         const Real n = getNumParameters();
85         Real sumcos = 0;
86         for (int i = 0; i < n; ++i) {
87             sumcos += cos(c * x[i]);
88         }
89         f = -a * exp(-b * x.normRMS()) - exp(sumcos / n) + a + SimTK::E;
90         return 0;
91     }
optimalParameters()92     Vector optimalParameters() const override {
93         Vector x(getNumParameters());
94         x.setToZero();
95         return x;
96     }
97 private:
98     const Real a;
99     const Real b;
100     const Real c;
101 };
102 
103 // A very complex 2D function http://www.sfu.ca/~ssurjano/drop.html
104 class DropWave : public TestOptimizerSystem {
105 public:
DropWave()106     DropWave() : TestOptimizerSystem(2) {
107         // The website above says this function usually has the following
108         // bounds:
109         Vector limits(getNumParameters());
110         limits.setTo(5.12);
111         setParameterLimits(-limits, limits);
112     }
objectiveFunc(const Vector & x,bool new_parameters,Real & f)113     int objectiveFunc(const Vector& x, bool new_parameters, Real& f) const override {
114         const Real dotprod = x[0] * x[0] + x[1] * x[1];
115         f = -(1 + cos(12 * sqrt(dotprod))) / (0.5 * dotprod + 2);
116         return 0;
117     }
optimalParameters()118     Vector optimalParameters() const override {
119         Vector x(getNumParameters());
120         x.setToZero();
121         return x;
122     }
123 };
124 
125 // Looks like a curved valley.
126 // https://en.wikipedia.org/wiki/Test_functions_for_optimization
127 class Rosenbrock : public TestOptimizerSystem {
128 public:
Rosenbrock(int nParameters)129     Rosenbrock(int nParameters) : TestOptimizerSystem(nParameters) {}
objectiveFunc(const Vector & x,bool new_parameters,Real & f)130     int objectiveFunc(const Vector& x, bool new_parameters, Real& f) const override {
131         f = 0;
132         for (int i = 0; i < getNumParameters() - 1; ++i) {
133             f += 100.0 * square(x[i+1] - square(x[i])) + square(x[i] - 1);
134         }
135         return 0;
136     }
optimalParameters()137     Vector optimalParameters() const override {
138         Vector x(getNumParameters());
139         x.setTo(1);
140         return x;
141     }
142 };
143 
144 // http://www.sfu.ca/~ssurjano/schwef.html
145 class Schwefel : public TestOptimizerSystem {
146 public:
Schwefel(int nParameters)147     Schwefel(int nParameters) : TestOptimizerSystem(nParameters) {
148         // The website above says this function usually has the following
149         // bounds:
150         Vector limits(nParameters);
151         limits.setTo(500);
152         setParameterLimits(-limits, limits);
153     }
objectiveFunc(const Vector & x,bool new_parameters,Real & f)154     int objectiveFunc(const Vector& x, bool new_parameters, Real& f) const override {
155         Real sum = 0;
156         for (int i = 0; i < getNumParameters(); ++i) {
157             sum += x[i] * sin(sqrt(std::abs(x[i])));
158         }
159         f = 418.9829 * getNumParameters() -sum;
160         return 0;
161     }
optimalParameters()162     Vector optimalParameters() const override {
163         Vector x(getNumParameters());
164         x.setTo(420.9687);
165         return x;
166     }
167 };
168 
169 class Easom : public TestOptimizerSystem {
170 public:
Easom()171     Easom() : TestOptimizerSystem(2) {
172         // The website above says this function usually has the following
173         // bounds:
174         Vector limits(getNumParameters());
175         limits.setTo(100);
176         setParameterLimits(-limits, limits);
177     }
objectiveFunc(const Vector & x,bool new_parameters,Real & f)178     int objectiveFunc(const Vector& x, bool new_parameters, Real& f) const override {
179         f = -cos(x[0]) * cos(x[1]) *
180             exp(-pow(x[0] - Pi, 2) - pow(x[1] - Pi, 2));
181         return 0;
182     }
optimalParameters()183     Vector optimalParameters() const override {
184         Vector x(getNumParameters());
185         x.setTo(Pi);
186         return x;
187     }
188 };
189 
190 
191 #endif // SimTK_SIMMATH_OPTIMIZER_SYSTEMS_H_
192