chiark / gitweb /
added simple convex test function
[nlopt.git] / test / testfuncs.c
1 #include <stdio.h>
2 #include <math.h>
3
4 #include "testfuncs.h"
5 #include "config.h"
6
7 #define UNUSED(x) (void) x
8
9 static double sqr(double x) { return x * x; }
10
11 int testfuncs_verbose = 0;
12 int testfuncs_counter = 0;
13
14 static double testfuncs_status(int n, const double *x, double f)
15 {
16      ++testfuncs_counter;
17      if (testfuncs_verbose) {
18           int i;
19           printf("f_%d (%g", testfuncs_counter, x[0]);
20           for (i = 1; i < n; ++i) printf(", %g", x[i]);
21           printf(") = %g\n", f);
22      }
23      return f;
24 }
25
26 #define RETURN(f) return testfuncs_status(n, x, f);
27
28 /****************************************************************************/
29 static double rosenbrock_f(int n, const double *x, double *grad, void *data)
30 {
31      double a = x[1] - x[0] * x[0], b = 1 - x[0];
32      UNUSED(data);
33      if (grad) {
34           grad[0] = -400 * a * x[0] - 2*b;
35           grad[1] = -200 * a;
36      }
37      RETURN(100 * sqr(a) + sqr(b));
38 }
39
40 static const double rosenbrock_lb[2] = {-2, -2};
41 static const double rosenbrock_ub[2] = {2, 2};
42 static const double rosenbrock_xmin[2] = {1, 1};
43
44 /****************************************************************************/
45 static double mccormic_f(int n, const double *x, double *grad, void *data)
46 {
47      double a = x[0] + x[1], b = x[0] - x[1];
48      UNUSED(data);
49      if (grad) {
50           grad[0] = cos(a) + 2*b - 1.5;
51           grad[1] = cos(a) - 2*b + 2.5;
52      }
53      RETURN(sin(a) + sqr(b) - 1.5*x[0] + 2.5*x[1] + 1);
54 }
55
56 static const double mccormic_lb[2] = {-1.5, -3};
57 static const double mccormic_ub[2] = {4, 4};
58 static const double mccormic_xmin[2] = {-0.54719, 1.54719};
59
60 /****************************************************************************/
61 static double boxbetts_f(int n, const double *x, double *grad, void *data)
62 {
63      int i;
64      double f = 0;
65      UNUSED(data);
66      if (grad)
67           grad[0] = grad[1] = grad[2] = 0;
68      for (i = 1; i <= 10; ++i) {
69           double e0 = exp(-0.1*i*x[0]);
70           double e1 = exp(-0.1*i*x[1]);
71           double e2 = exp(-0.1*i) - exp((double) -i);
72           double g = e0 - e1 - e2 * x[2];
73           f += sqr(g);
74           if (grad) {
75                grad[0] += (2 * g) * (-0.1*i*e0);
76                grad[1] += (2 * g) * (0.1*i*e1);
77                grad[2] += -(2 * g) * e2;
78           }
79      }
80      RETURN(f);
81 }
82
83 static const double boxbetts_lb[3] = {0.9,9,0.9};
84 static const double boxbetts_ub[3] = {1.2,11.2,1.2};
85 static const double boxbetts_xmin[3] = {1,10,1};
86
87 /****************************************************************************/
88 static double paviani_f(int n, const double *x, double *grad, void *data)
89 {
90      int i;
91      double f = 0, prod = 1;
92      UNUSED(data);
93      if (grad) for (i = 0; i < 10; ++i) grad[i] = 0;
94      for (i = 0; i < 10; ++i) {
95           double ln1 = log(x[i] - 2);
96           double ln2 = log(10 - x[i]);
97           f += sqr(ln1) + sqr(ln2);
98           if (grad)
99                grad[i] += 2 * ln1 / (x[i] - 2) - 2 * ln2 / (10 - x[i]);
100           prod *= x[i];
101      }
102      f -= (prod = pow(prod, 0.2));
103      if (grad)
104           for (i = 0; i < 10; ++i)
105                grad[i] -= 0.2 * prod / x[i];
106      RETURN(f);
107 }
108
109 static const double paviani_lb[10] = {2.001,2.001,2.001,2.001,2.001,2.001,2.001,2.001,2.001,2.001};
110 static const double paviani_ub[10] = {9.999,9.999,9.999,9.999,9.999,9.999,9.999,9.999,9.999,9.999};
111 static const double paviani_xmin[10] = {9.350266,9.350266,9.350266,9.350266,9.350266,9.350266,9.350266,9.350266,9.350266,9.350266};
112
113 /****************************************************************************/
114 static double grosenbrock_f(int n, const double *x, double *grad, void *data)
115 {
116      int i;
117      double f = 0;
118      UNUSED(data);
119      if (grad) grad[0] = 0;
120      for (i = 0; i < 29; ++i) {
121           double a = x[i+1] - x[i] * x[i], b = 1 - x[i];
122           if (grad) {
123                grad[i] += -400 * a * x[i] - 2*b;
124                grad[i+1] = -200 * a;
125           }
126           f += 100 * sqr(a) + sqr(b);
127      }
128      RETURN(f);
129 }
130
131 static const double grosenbrock_lb[30] = {-30,-30,-30,-30,-30,-30,-30,-30,-30,-30,-30,-30,-30,-30,-30,-30,-30,-30,-30,-30,-30,-30,-30,-30,-30,-30,-30,-30,-30,-30};
132 static const double grosenbrock_ub[30] = {30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30,30};
133 static const double grosenbrock_xmin[30] = {1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1};
134
135 /****************************************************************************/
136 static double goldsteinprice_f(int n, const double *x, double *grad, void *data)
137 {
138      double x0, x1, a1, a12, a2, b1, b12, b2;
139      UNUSED(data);
140      x0 = x[0]; x1 = x[1];
141      a1 = x0+x1+1; a12 = sqr(a1);
142      a2 = 19 - 14*x0 + 3*x0*x0 - 14*x1 + 6*x0*x1 + 3*x1*x1;
143      b1 = 2*x0-3*x1; b12 = sqr(b1);
144      b2 = 18 - 32*x0 + 12*x0*x0 + 48*x1 - 36*x0*x1 + 27*x1*x1;
145      if (grad) {
146           grad[0] = (1 + a12 * a2) * (2 * b1 * 2 * b2 
147                                       + b12 * (-32 + 24*x0 - 36*x1))
148                + (2 * a1 * a2 + a12 * (-14 + 6*x0 + 6*x1)) * (30 + b12 * b2);
149           grad[1] = (1 + a12 * a2) * (2 * b1 * (-3) * b2
150                                       + b12 * (48 - 36*x0 + 54 * x1))
151                + (2 * a1 * a2 + a12 * (-14 + 6*x0 + 6*x1)) * (30 + b12 * b2);
152      }
153      RETURN((1 + a12 * a2) * (30 + b12 * b2));
154 }
155
156 static const double goldsteinprice_lb[2] = {-2, -2};
157 static const double goldsteinprice_ub[2] = {2, 2};
158 static const double goldsteinprice_xmin[2] = {0, -1};
159
160 /****************************************************************************/
161 static double shekel_f(int n, const double *x, double *grad, void *data)
162 {
163      static const double A[10][4] = { {4,4,4,4},
164                                       {1,1,1,1},
165                                       {8,8,8,8},
166                                       {6,6,6,6},
167                                       {3,7,3,7},
168                                       {2,9,2,9},
169                                       {5,5,3,3},
170                                       {8,1,8,1},
171                                       {6,2,6,2},
172                                       {7,3.6,7,3.6} };
173      static const double c[10] = {.1,.2,.2,.4,.4,.6,.3,.7,.5,.5};
174      int i;
175      double f = 0;
176      if (grad) for (i = 0; i < n; ++i) grad[i] = 0;
177      int m = *((int *) data);
178      for (i = 0; i < m; ++i) {
179           double fi = 1.0 / (c[i] 
180                              + sqr(x[0]-A[i][0])
181                              + sqr(x[1]-A[i][1])
182                              + sqr(x[2]-A[i][2])
183                              + sqr(x[3]-A[i][3]));
184           f -= fi;
185           if (grad) {
186                grad[0] += (2*fi*fi) * (x[0]-A[i][0]);
187                grad[1] += (2*fi*fi) * (x[1]-A[i][1]);
188                grad[2] += (2*fi*fi) * (x[2]-A[i][2]);
189                grad[3] += (2*fi*fi) * (x[3]-A[i][3]);
190           }
191      }
192      RETURN(f);
193 }
194
195 static int shekel_m[3] = {5,7,10};
196 static const double shekel_lb[4] = {0,0,0,0};
197 static const double shekel_ub[4] = {10,10,10,10};
198 static const double shekel0_xmin[4] = {4.00004,4.00013,4.00004,4.00013};
199 static const double shekel1_xmin[4] = {4.00057,4.00069,3.99949,3.99961};
200 static const double shekel2_xmin[4] = {4.00075,4.00059,3.99966,3.99951};
201
202 /****************************************************************************/
203 #define PI3 9.424777960769379 /* 3*pi */
204 #define PI2 6.283185307179586 /* 2*pi */
205 static double levy_f(int n, const double *x, double *grad, void *data)
206 {
207      UNUSED(data);
208      int i;
209      double a = x[n-1] - 1, b = 1 + sqr(sin(PI2*x[n-1]));
210      double f = sqr(sin(PI3*x[0])) + a * b;
211      if (grad) {
212           for (i = 0; i < n; ++i) grad[i] = 0;
213           grad[0] = 2 * PI3 * sin(PI3*x[0]) * cos(PI3*x[0]);
214           grad[n-1] += b + a * 2 * PI2 * sin(PI2*x[n-1]) * cos(PI2*x[n-1]);
215      }
216      for (i = 0; i < n-1; ++i) {
217           a = x[i] - 1;
218           b = 1 + sqr(sin(PI3*x[i+1]));
219           f += sqr(a) * b;
220           if (grad) {
221                grad[i] += 2 * a * b;
222                grad[i+1] += 2*PI3 * sqr(a) * sin(PI3*x[i+1])*cos(PI3*x[i+1]);
223           }
224      }
225      RETURN(f);
226 }
227
228 static const double levy_lb[7] = {-5,-5,-5,-5,-5,-5,-5};
229 static const double levy_ub[7] = {5,5,5,5,5,5,5};
230 static const double levy_xmin[7] = {1,1,1,1,1,1,-4.754402};
231 static const double levy4_lb[4] = {-10,-10,-10,-10};
232 static const double levy4_ub[4] = {10,10,10,10};
233 static const double levy4_xmin[4] = {1,1,1,-9.752356};
234
235 /****************************************************************************/
236 static double griewank_f(int n, const double *x, double *grad, void *data)
237 {
238      int i;
239      double f = 1, p = 1;
240      UNUSED(data);
241      for (i = 0; i < n; ++i) {
242           f += sqr(x[i]) * 0.00025;
243           p *= cos(x[i] / sqrt(i + 1.));
244           if (grad) grad[i] = x[i] * 0.0005;
245      }
246      f -= p;
247      if (grad)
248           for (i = 0; i < n; ++i)
249                grad[i] += p * tan(x[i] / sqrt(i + 1.)) / sqrt(i + 1.);
250      RETURN(f);
251 }
252
253 static const double griewank_lb[10] = {-500,-500,-500,-500,-500,-500,-500,-500,-500,-500};
254 static const double griewank_ub[10] = {600,600,600,600,600,600,600,600,600,600};
255 static const double griewank_xmin[10] = {0,0,0,0,0,0,0,0,0,0};
256
257 /****************************************************************************/
258 static double sixhumpcamel_f(int n, const double *x, double *grad, void *data)
259 {
260      UNUSED(data);
261      if (grad) {
262           grad[0] = 8*x[0] - 2.1*4*pow(x[0],3.) + 2*pow(x[0],5.) + x[1];
263           grad[1] = x[0] - 8*x[1] + 16*pow(x[1],3.);
264      }
265      RETURN(4*sqr(x[0]) - 2.1 * pow(x[0],4.) + pow(x[0],6.)/3. 
266             + x[0]*x[1] - 4*sqr(x[1]) + 4*pow(x[1],4.));
267 }
268
269 static const double sixhumpcamel_lb[2] = {-5,-5};
270 static const double sixhumpcamel_ub[2] = {5,5};
271 static const double sixhumpcamel_xmin[2] = {0.08984, -0.71266};
272
273 /****************************************************************************/
274 static double convexcosh_f(int n, const double *x, double *grad, void *data)
275 {
276      int i;
277      double f = 1;
278      UNUSED(data);
279      for (i = 0; i < n; ++i)
280           f *= cosh((x[i] - i) * (i+1));
281      if (grad)
282           for (i = 0; i < n; ++i)
283                grad[i] = f * tanh((x[i] - i) * (i+1)) * (i+1);
284      RETURN(f);
285 }
286
287 static const double convexcosh_lb[10] = {-1,0,0,0,0,0,0,0,0,0};
288 static const double convexcosh_ub[10] = {2,3,6,7,8,10,11,13,14,16};
289 static const double convexcosh_xmin[10] = {0,1,2,3,4,5,6,7,8,9};
290
291 /****************************************************************************/
292 /****************************************************************************/
293
294 const testfunc testfuncs[NTESTFUNCS] = {
295      { rosenbrock_f, NULL, 1, 2,
296        rosenbrock_lb, rosenbrock_ub, rosenbrock_xmin,
297        0.0, "Rosenbrock function" },
298      { mccormic_f, NULL, 1, 2,
299        mccormic_lb, mccormic_ub, mccormic_xmin,
300        -1.9133, "McCormic function" },
301      { boxbetts_f, NULL, 1, 3,
302        boxbetts_lb, boxbetts_ub, boxbetts_xmin,
303        0.0, "Box and Betts exponential quadratic sum" },
304      { paviani_f, NULL, 1, 10,
305        paviani_lb, paviani_ub, paviani_xmin,
306        -45.778470, "Paviani function" },
307      { grosenbrock_f, NULL, 1, 30,
308        grosenbrock_lb, grosenbrock_ub, grosenbrock_xmin,
309        0.0, "Generalized Rosenbrock function" },
310      { goldsteinprice_f, NULL, 1, 2,
311        goldsteinprice_lb, goldsteinprice_ub, goldsteinprice_xmin,
312        3.0, "Goldstein and Price function" },
313      { shekel_f, shekel_m + 0, 1, 4,
314        shekel_lb, shekel_ub, shekel0_xmin,
315        -10.1532, "Shekel m=5 function" },
316      { shekel_f, shekel_m + 1, 1, 4,
317        shekel_lb, shekel_ub, shekel1_xmin,
318        -10.4029, "Shekel m=7 function" },
319      { shekel_f, shekel_m + 2, 1, 4,
320        shekel_lb, shekel_ub, shekel2_xmin,
321        -10.5364, "Shekel m=10 function" },
322      { levy_f, NULL, 1, 4,
323        levy4_lb, levy4_ub, levy4_xmin,
324        -21.502356, "Levy n=4 function" },
325      { levy_f, NULL, 1, 5,
326        levy_lb, levy_ub, levy_xmin+2,
327        -11.504403, "Levy n=5 function" },
328      { levy_f, NULL, 1, 6,
329        levy_lb, levy_ub, levy_xmin+1,
330        -11.504403, "Levy n=6 function" },
331      { levy_f, NULL, 1, 7,
332        levy_lb, levy_ub, levy_xmin,
333        -11.504403, "Levy n=7 function" },
334      { griewank_f, NULL, 1, 10,
335        griewank_lb, griewank_ub, griewank_xmin,
336        0.0, "Griewank function" },
337      { sixhumpcamel_f, NULL, 1, 2,
338        sixhumpcamel_lb, sixhumpcamel_ub, sixhumpcamel_xmin,
339        -1.03163, "Six-hump camel back function" },
340      { convexcosh_f, NULL, 1, 10,
341        convexcosh_lb, convexcosh_ub, convexcosh_xmin,
342        1.0, "Convex product of cosh functions" }
343 };