View Javadoc
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    *      https://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   */
17  
18  /*
19   * This is not the original file distributed by the Apache Software Foundation
20   * It has been modified by the Hipparchus project
21   */
22  
23  package org.hipparchus.optim.nonlinear.vector.leastsquares;
24  
25  import org.hipparchus.analysis.MultivariateMatrixFunction;
26  import org.hipparchus.analysis.MultivariateVectorFunction;
27  import org.hipparchus.exception.MathIllegalStateException;
28  import org.hipparchus.linear.DiagonalMatrix;
29  import org.hipparchus.optim.nonlinear.vector.leastsquares.LeastSquaresOptimizer.Optimum;
30  import org.hipparchus.util.FastMath;
31  import org.junit.jupiter.api.Test;
32  
33  import java.util.Arrays;
34  
35  import static org.junit.jupiter.api.Assertions.assertEquals;
36  import static org.junit.jupiter.api.Assertions.assertFalse;
37  import static org.junit.jupiter.api.Assertions.assertTrue;
38  
39  /**
40   * <p>Some of the unit tests are re-implementations of the MINPACK <a
41   * href="http://www.netlib.org/minpack/ex/file17">file17</a> and <a
42   * href="http://www.netlib.org/minpack/ex/file22">file22</a> test files.
43   * The redistribution policy for MINPACK is available <a
44   * href="http://www.netlib.org/minpack/disclaimer">here</a>, for
45   * convenience, it is reproduced below.</p>
46   
47   * <table border="0" width="80%" cellpadding="10" align="center" bgcolor="#E0E0E0">
48   * <tr><td>
49   *    Minpack Copyright Notice (1999) University of Chicago.
50   *    All rights reserved
51   * </td></tr>
52   * <tr><td>
53   * Redistribution and use in source and binary forms, with or without
54   * modification, are permitted provided that the following conditions
55   * are met:
56   * <ol>
57   *  <li>Redistributions of source code must retain the above copyright
58   *      notice, this list of conditions and the following disclaimer.</li>
59   * <li>Redistributions in binary form must reproduce the above
60   *     copyright notice, this list of conditions and the following
61   *     disclaimer in the documentation and/or other materials provided
62   *     with the distribution.</li>
63   * <li>The end-user documentation included with the redistribution, if any,
64   *     must include the following acknowledgment:
65   *     <code>This product includes software developed by the University of
66   *           Chicago, as Operator of Argonne National Laboratory.</code>
67   *     Alternately, this acknowledgment may appear in the software itself,
68   *     if and wherever such third-party acknowledgments normally appear.</li>
69   * <li><strong>WARRANTY DISCLAIMER. THE SOFTWARE IS SUPPLIED "AS IS"
70   *     WITHOUT WARRANTY OF ANY KIND. THE COPYRIGHT HOLDER, THE
71   *     UNITED STATES, THE UNITED STATES DEPARTMENT OF ENERGY, AND
72   *     THEIR EMPLOYEES: (1) DISCLAIM ANY WARRANTIES, EXPRESS OR
73   *     IMPLIED, INCLUDING BUT NOT LIMITED TO ANY IMPLIED WARRANTIES
74   *     OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE, TITLE
75   *     OR NON-INFRINGEMENT, (2) DO NOT ASSUME ANY LEGAL LIABILITY
76   *     OR RESPONSIBILITY FOR THE ACCURACY, COMPLETENESS, OR
77   *     USEFULNESS OF THE SOFTWARE, (3) DO NOT REPRESENT THAT USE OF
78   *     THE SOFTWARE WOULD NOT INFRINGE PRIVATELY OWNED RIGHTS, (4)
79   *     DO NOT WARRANT THAT THE SOFTWARE WILL FUNCTION
80   *     UNINTERRUPTED, THAT IT IS ERROR-FREE OR THAT ANY ERRORS WILL
81   *     BE CORRECTED.</strong></li>
82   * <li><strong>LIMITATION OF LIABILITY. IN NO EVENT WILL THE COPYRIGHT
83   *     HOLDER, THE UNITED STATES, THE UNITED STATES DEPARTMENT OF
84   *     ENERGY, OR THEIR EMPLOYEES: BE LIABLE FOR ANY INDIRECT,
85   *     INCIDENTAL, CONSEQUENTIAL, SPECIAL OR PUNITIVE DAMAGES OF
86   *     ANY KIND OR NATURE, INCLUDING BUT NOT LIMITED TO LOSS OF
87   *     PROFITS OR LOSS OF DATA, FOR ANY REASON WHATSOEVER, WHETHER
88   *     SUCH LIABILITY IS ASSERTED ON THE BASIS OF CONTRACT, TORT
89   *     (INCLUDING NEGLIGENCE OR STRICT LIABILITY), OR OTHERWISE,
90   *     EVEN IF ANY OF SAID PARTIES HAS BEEN WARNED OF THE
91   *     POSSIBILITY OF SUCH LOSS OR DAMAGES.</strong></li>
92   * <ol></td></tr>
93   * </table>
94   
95   * @author Argonne National Laboratory. MINPACK project. March 1980 (original fortran minpack tests)
96   * @author Burton S. Garbow (original fortran minpack tests)
97   * @author Kenneth E. Hillstrom (original fortran minpack tests)
98   * @author Jorge J. More (original fortran minpack tests)
99   * @author Luc Maisonobe (non-minpack tests and minpack tests Java translation)
100  */
101 class MinpackTest {
102 
103     @Test
104     void testMinpackLinearFullRank() {
105         minpackTest(new LinearFullRankFunction(10, 5, 1.0,
106                                                5.0, 2.23606797749979), false);
107         minpackTest(new LinearFullRankFunction(50, 5, 1.0,
108                                                8.06225774829855, 6.70820393249937), false);
109     }
110 
111     @Test
112     void testMinpackLinearRank1() {
113         minpackTest(new LinearRank1Function(10, 5, 1.0,
114                                             291.521868819476, 1.4638501094228), false);
115         minpackTest(new LinearRank1Function(50, 5, 1.0,
116                                             3101.60039334535, 3.48263016573496), false);
117     }
118 
119     @Test
120     void testMinpackLinearRank1ZeroColsAndRows() {
121         minpackTest(new LinearRank1ZeroColsAndRowsFunction(10, 5, 1.0), false);
122         minpackTest(new LinearRank1ZeroColsAndRowsFunction(50, 5, 1.0), false);
123     }
124 
125     @Test
126     void testMinpackRosenbrok() {
127         minpackTest(new RosenbrockFunction(new double[] { -1.2, 1.0 },
128                                            FastMath.sqrt(24.2)), false);
129         minpackTest(new RosenbrockFunction(new double[] { -12.0, 10.0 },
130                                            FastMath.sqrt(1795769.0)), false);
131         minpackTest(new RosenbrockFunction(new double[] { -120.0, 100.0 },
132                                            11.0 * FastMath.sqrt(169000121.0)), false);
133     }
134 
135     @Test
136     void testMinpackHelicalValley() {
137         minpackTest(new HelicalValleyFunction(new double[] { -1.0, 0.0, 0.0 },
138                                               50.0), false);
139         minpackTest(new HelicalValleyFunction(new double[] { -10.0, 0.0, 0.0 },
140                                               102.95630140987), false);
141         minpackTest(new HelicalValleyFunction(new double[] { -100.0, 0.0, 0.0},
142                                               991.261822123701), false);
143     }
144 
145     @Test
146     void testMinpackPowellSingular() {
147         minpackTest(new PowellSingularFunction(new double[] { 3.0, -1.0, 0.0, 1.0 },
148                                                14.6628782986152), false);
149         minpackTest(new PowellSingularFunction(new double[] { 30.0, -10.0, 0.0, 10.0 },
150                                                1270.9838708654), false);
151         minpackTest(new PowellSingularFunction(new double[] { 300.0, -100.0, 0.0, 100.0 },
152                                                126887.903284750), false);
153     }
154 
155     @Test
156     void testMinpackFreudensteinRoth() {
157         minpackTest(new FreudensteinRothFunction(new double[] { 0.5, -2.0 },
158                                                  20.0124960961895, 6.99887517584575,
159                                                  new double[] {
160                                                      11.4124844654993,
161                                                      -0.896827913731509
162                                                  }), false);
163         minpackTest(new FreudensteinRothFunction(new double[] { 5.0, -20.0 },
164                                                  12432.833948863, 6.9988751744895,
165                                                  new double[] {
166                                                      11.41300466147456,
167                                                      -0.896796038685959
168                                                  }), false);
169         minpackTest(new FreudensteinRothFunction(new double[] { 50.0, -200.0 },
170                                                  11426454.595762, 6.99887517242903,
171                                                  new double[] {
172                                                      11.412781785788564,
173                                                      -0.8968051074920405
174                                                  }), false);
175     }
176 
177     @Test
178     void testMinpackBard() {
179         minpackTest(new BardFunction(1.0, 6.45613629515967, 0.0906359603390466,
180                                      new double[] {
181                                          0.0824105765758334,
182                                          1.1330366534715,
183                                          2.34369463894115
184                                      }), false);
185         minpackTest(new BardFunction(10.0, 36.1418531596785, 4.17476870138539,
186                                      new double[] {
187                                          0.840666673818329,
188                                          -158848033.259565,
189                                          -164378671.653535
190                                      }), false);
191         minpackTest(new BardFunction(100.0, 384.114678637399, 4.17476870135969,
192                                      new double[] {
193                                          0.840666673867645,
194                                          -158946167.205518,
195                                          -164464906.857771
196                                      }), false);
197     }
198 
199     @Test
200     void testMinpackKowalikOsborne() {
201         minpackTest(new KowalikOsborneFunction(new double[] { 0.25, 0.39, 0.415, 0.39 },
202                                                0.0728915102882945,
203                                                0.017535837721129,
204                                                new double[] {
205                                                    0.192807810476249,
206                                                    0.191262653354071,
207                                                    0.123052801046931,
208                                                    0.136053221150517
209                                                }), false);
210         minpackTest(new KowalikOsborneFunction(new double[] { 2.5, 3.9, 4.15, 3.9 },
211                                                2.97937007555202,
212                                                0.032052192917937,
213                                                new double[] {
214                                                    728675.473768287,
215                                                    -14.0758803129393,
216                                                    -32977797.7841797,
217                                                    -20571594.1977912
218                                                }), false);
219         minpackTest(new KowalikOsborneFunction(new double[] { 25.0, 39.0, 41.5, 39.0 },
220                                                29.9590617016037,
221                                                0.0175364017658228,
222                                                new double[] {
223                                                    0.192948328597594,
224                                                    0.188053165007911,
225                                                    0.122430604321144,
226                                                    0.134575665392506
227                                                }), false);
228     }
229 
230     @Test
231     void testMinpackMeyer() {
232         minpackTest(new MeyerFunction(new double[] { 0.02, 4000.0, 250.0 },
233                                       41153.4665543031, 9.37794514651874,
234                                       new double[] {
235                                           0.00560963647102661,
236                                           6181.34634628659,
237                                           345.223634624144
238                                       }), false);
239         minpackTest(new MeyerFunction(new double[] { 0.2, 40000.0, 2500.0 },
240                                       4168216.89130846, 792.917871779501,
241                                       new double[] {
242                                           1.42367074157994e-11,
243                                           33695.7133432541,
244                                           901.268527953801
245                                       }), true);
246     }
247 
248     @Test
249     void testMinpackWatson() {
250         minpackTest(new WatsonFunction(6, 0.0,
251                                        5.47722557505166, 0.0478295939097601,
252                                        new double[] {
253                                            -0.0157249615083782, 1.01243488232965,
254                                            -0.232991722387673,  1.26043101102818,
255                                            -1.51373031394421,   0.99299727291842
256                                        }), false);
257         minpackTest(new WatsonFunction(6, 10.0,
258                                        6433.12578950026, 0.0478295939096951,
259                                        new double[] {
260                                            -0.0157251901386677, 1.01243485860105,
261                                            -0.232991545843829,  1.26042932089163,
262                                            -1.51372776706575,   0.99299573426328
263                                        }), false);
264         minpackTest(new WatsonFunction(6, 100.0,
265                                        674256.040605213, 0.047829593911544,
266                                        new double[] {
267                                            -0.0157247019712586, 1.01243490925658,
268                                            -0.232991922761641,  1.26043292929555,
269                                            -1.51373320452707,   0.99299901922322
270                                        }), false);
271         minpackTest(new WatsonFunction(9, 0.0,
272                                        5.47722557505166, 0.00118311459212420,
273                                        new double[] {
274                                            -0.153070644166722e-4, 0.999789703934597,
275                                            0.0147639634910978,   0.146342330145992,
276                                            1.00082109454817,    -2.61773112070507,
277                                            4.10440313943354,    -3.14361226236241,
278                                            1.05262640378759
279                                        }), false);
280         minpackTest(new WatsonFunction(9, 10.0,
281                                        12088.127069307, 0.00118311459212513,
282                                        new double[] {
283                                            -0.153071334849279e-4, 0.999789703941234,
284                                            0.0147639629786217,   0.146342334818836,
285                                            1.00082107321386,    -2.61773107084722,
286                                            4.10440307655564,    -3.14361222178686,
287                                            1.05262639322589
288                                        }), false);
289         minpackTest(new WatsonFunction(9, 100.0,
290                                        1269109.29043834, 0.00118311459212384,
291                                        new double[] {
292                                            -0.153069523352176e-4, 0.999789703958371,
293                                            0.0147639625185392,   0.146342341096326,
294                                            1.00082104729164,    -2.61773101573645,
295                                            4.10440301427286,    -3.14361218602503,
296                                            1.05262638516774
297                                        }), false);
298         minpackTest(new WatsonFunction(12, 0.0,
299                                        5.47722557505166, 0.217310402535861e-4,
300                                        new double[] {
301                                            -0.660266001396382e-8, 1.00000164411833,
302                                            -0.000563932146980154, 0.347820540050756,
303                                            -0.156731500244233,    1.05281515825593,
304                                            -3.24727109519451,     7.2884347837505,
305                                            -10.271848098614,       9.07411353715783,
306                                            -4.54137541918194,     1.01201187975044
307                                        }), false);
308         minpackTest(new WatsonFunction(12, 10.0,
309                                        19220.7589790951, 0.217310402518509e-4,
310                                        new double[] {
311                                            -0.663710223017410e-8, 1.00000164411787,
312                                            -0.000563932208347327, 0.347820540486998,
313                                            -0.156731503955652,    1.05281517654573,
314                                            -3.2472711515214,      7.28843489430665,
315                                            -10.2718482369638,      9.07411364383733,
316                                            -4.54137546533666,     1.01201188830857
317                                        }), false);
318         minpackTest(new WatsonFunction(12, 100.0,
319                                        2018918.04462367, 0.217310402539845e-4,
320                                        new double[] {
321                                            -0.663806046485249e-8, 1.00000164411786,
322                                            -0.000563932210324959, 0.347820540503588,
323                                            -0.156731504091375,    1.05281517718031,
324                                            -3.24727115337025,     7.28843489775302,
325                                            -10.2718482410813,      9.07411364688464,
326                                            -4.54137546660822,     1.0120118885369
327                                        }), false);
328     }
329 
330     @Test
331     void testMinpackBox3Dimensional() {
332         minpackTest(new Box3DimensionalFunction(10, new double[] { 0.0, 10.0, 20.0 },
333                                                 32.1115837449572), false);
334     }
335 
336     @Test
337     void testMinpackJennrichSampson() {
338         minpackTest(new JennrichSampsonFunction(10, new double[] { 0.3, 0.4 },
339                                                 64.5856498144943, 11.1517793413499,
340                                                 new double[] {
341 //                                                     0.2578330049, 0.257829976764542
342                                                     0.2578199266368004, 0.25782997676455244
343                                                 }), false);
344     }
345 
346     @Test
347     void testMinpackBrownDennis() {
348         minpackTest(new BrownDennisFunction(20,
349                                             new double[] { 25.0, 5.0, -5.0, -1.0 },
350                                             2815.43839161816, 292.954288244866,
351                                             new double[] {
352                                                 -11.59125141003, 13.2024883984741,
353                                                 -0.403574643314272, 0.236736269844604
354                                             }), false);
355         minpackTest(new BrownDennisFunction(20,
356                                             new double[] { 250.0, 50.0, -50.0, -10.0 },
357                                             555073.354173069, 292.954270581415,
358                                             new double[] {
359                                                 -11.5959274272203, 13.2041866926242,
360                                                 -0.403417362841545, 0.236771143410386
361                                             }), false);
362         minpackTest(new BrownDennisFunction(20,
363                                             new double[] { 2500.0, 500.0, -500.0, -100.0 },
364                                             61211252.2338581, 292.954306151134,
365                                             new double[] {
366                                                 -11.5902596937374, 13.2020628854665,
367                                                 -0.403688070279258, 0.236665033746463
368                                             }), false);
369     }
370 
371     @Test
372     void testMinpackChebyquad() {
373         minpackTest(new ChebyquadFunction(1, 8, 1.0,
374                                           1.88623796907732, 1.88623796907732,
375                                           new double[] { 0.5 }), false);
376         minpackTest(new ChebyquadFunction(1, 8, 10.0,
377                                           5383344372.34005, 1.88424820499951,
378                                           new double[] { 0.9817314924684 }), false);
379         minpackTest(new ChebyquadFunction(1, 8, 100.0,
380                                           0.118088726698392e19, 1.88424820499347,
381                                           new double[] { 0.9817314852934 }), false);
382         minpackTest(new ChebyquadFunction(8, 8, 1.0,
383                                           0.196513862833975, 0.0593032355046727,
384                                           new double[] {
385                                               0.0431536648587336, 0.193091637843267,
386                                               0.266328593812698,  0.499999334628884,
387                                               0.500000665371116,  0.733671406187302,
388                                               0.806908362156733,  0.956846335141266
389                                           }), false);
390         minpackTest(new ChebyquadFunction(9, 9, 1.0,
391                                           0.16994993465202, 0.0,
392                                           new double[] {
393                                               0.0442053461357828, 0.199490672309881,
394                                               0.23561910847106,   0.416046907892598,
395                                               0.5,                0.583953092107402,
396                                               0.764380891528940,  0.800509327690119,
397                                               0.955794653864217
398                                           }), false);
399         minpackTest(new ChebyquadFunction(10, 10, 1.0,
400                                           0.183747831178711, 0.0806471004038253,
401                                           new double[] {
402                                               0.0596202671753563, 0.166708783805937,
403                                               0.239171018813509,  0.398885290346268,
404                                               0.398883667870681,  0.601116332129320,
405                                               0.60111470965373,   0.760828981186491,
406                                               0.833291216194063,  0.940379732824644
407                                           }), false);
408     }
409 
410     @Test
411     void testMinpackBrownAlmostLinear() {
412         minpackTest(new BrownAlmostLinearFunction(10, 0.5,
413                                                   16.5302162063499, 0.0,
414                                                   new double[] {
415                                                       0.979430303349862, 0.979430303349862,
416                                                       0.979430303349862, 0.979430303349862,
417                                                       0.979430303349862, 0.979430303349862,
418                                                       0.979430303349862, 0.979430303349862,
419                                                       0.979430303349862, 1.20569696650138
420                                                   }), false);
421         minpackTest(new BrownAlmostLinearFunction(10, 5.0,
422                                                   9765624.00089211, 0.0,
423                                                   new double[] {
424                                                       0.979430303349865, 0.979430303349865,
425                                                       0.979430303349865, 0.979430303349865,
426                                                       0.979430303349865, 0.979430303349865,
427                                                       0.979430303349865, 0.979430303349865,
428                                                       0.979430303349865, 1.20569696650135
429                                                   }), false);
430         minpackTest(new BrownAlmostLinearFunction(10, 50.0,
431                                                   0.9765625e17, 0.0,
432                                                   new double[] {
433                                                       1.0, 1.0, 1.0, 1.0, 1.0,
434                                                       1.0, 1.0, 1.0, 1.0, 1.0
435                                                   }), false);
436         minpackTest(new BrownAlmostLinearFunction(30, 0.5,
437                                                   83.476044467848, 0.0,
438                                                   new double[] {
439                                                       0.997754216442807, 0.997754216442807,
440                                                       0.997754216442807, 0.997754216442807,
441                                                       0.997754216442807, 0.997754216442807,
442                                                       0.997754216442807, 0.997754216442807,
443                                                       0.997754216442807, 0.997754216442807,
444                                                       0.997754216442807, 0.997754216442807,
445                                                       0.997754216442807, 0.997754216442807,
446                                                       0.997754216442807, 0.997754216442807,
447                                                       0.997754216442807, 0.997754216442807,
448                                                       0.997754216442807, 0.997754216442807,
449                                                       0.997754216442807, 0.997754216442807,
450                                                       0.997754216442807, 0.997754216442807,
451                                                       0.997754216442807, 0.997754216442807,
452                                                       0.997754216442807, 0.997754216442807,
453                                                       0.997754216442807, 1.06737350671578
454                                                   }), false);
455         minpackTest(new BrownAlmostLinearFunction(40, 0.5,
456                                                   128.026364472323, 0.0,
457                                                   new double[] {
458                                                       1.00000000000002, 1.00000000000002,
459                                                       1.00000000000002, 1.00000000000002,
460                                                       1.00000000000002, 1.00000000000002,
461                                                       1.00000000000002, 1.00000000000002,
462                                                       1.00000000000002, 1.00000000000002,
463                                                       1.00000000000002, 1.00000000000002,
464                                                       1.00000000000002, 1.00000000000002,
465                                                       1.00000000000002, 1.00000000000002,
466                                                       1.00000000000002, 1.00000000000002,
467                                                       1.00000000000002, 1.00000000000002,
468                                                       1.00000000000002, 1.00000000000002,
469                                                       1.00000000000002, 1.00000000000002,
470                                                       1.00000000000002, 1.00000000000002,
471                                                       1.00000000000002, 1.00000000000002,
472                                                       1.00000000000002, 1.00000000000002,
473                                                       1.00000000000002, 1.00000000000002,
474                                                       1.00000000000002, 1.00000000000002,
475                                                       0.999999999999121
476                                                   }), false);
477     }
478 
479     @Test
480     void testMinpackOsborne1() {
481         minpackTest(new Osborne1Function(new double[] { 0.5, 1.5, -1.0, 0.01, 0.02, },
482                                          0.937564021037838, 0.00739249260904843,
483                                          new double[] {
484                                              0.375410049244025, 1.93584654543108,
485                                              -1.46468676748716, 0.0128675339110439,
486                                              0.0221227011813076
487                                          }), false);
488     }
489 
490     @Test
491     void testMinpackOsborne2() {
492         minpackTest(new Osborne2Function(new double[] {
493                     1.3, 0.65, 0.65, 0.7, 0.6,
494                     3.0, 5.0, 7.0, 2.0, 4.5, 5.5
495                 },
496                 1.44686540984712, 0.20034404483314,
497                 new double[] {
498                     1.30997663810096,  0.43155248076,
499                     0.633661261602859, 0.599428560991695,
500                     0.754179768272449, 0.904300082378518,
501                     1.36579949521007, 4.82373199748107,
502                     2.39868475104871, 4.56887554791452,
503                     5.67534206273052
504                 }), false);
505     }
506 
507     private void minpackTest(MinpackFunction function, boolean exceptionExpected) {
508         final double tol = 2.22044604926e-16;
509         final double sqrtTol = FastMath.sqrt(tol);
510 
511         LevenbergMarquardtOptimizer optimizer = new LevenbergMarquardtOptimizer()
512             .withCostRelativeTolerance(sqrtTol)
513             .withParameterRelativeTolerance(sqrtTol)
514             .withOrthoTolerance(tol);
515 
516         LeastSquaresProblem problem = new LeastSquaresBuilder()
517                 .maxEvaluations(400 * (function.getN() + 1))
518                 .maxIterations(2000)
519                 .model(function.getModelFunction(), function.getModelFunctionJacobian())
520                 .target(function.getTarget())
521                 .weight(new DiagonalMatrix(function.getWeight()))
522                 .start(function.getStartPoint())
523                 .build();
524 
525         try {
526             final Optimum optimum = optimizer.optimize(problem);
527             assertFalse(exceptionExpected);
528             function.checkTheoreticalMinCost(optimum.getRMS());
529             function.checkTheoreticalMinParams(optimum.getPoint().toArray());
530         } catch (MathIllegalStateException e) {
531             assertTrue(exceptionExpected);
532         }
533     }
534 
535     private static abstract class MinpackFunction {
536         protected int      n;
537         protected int      m;
538         protected double[] startParams;
539         protected double   theoreticalMinCost;
540         protected double[] theoreticalMinParams;
541         protected double   costAccuracy;
542         protected double   paramsAccuracy;
543 
544         protected MinpackFunction(int m, double[] startParams,
545                                   double theoreticalMinCost,
546                                   double[] theoreticalMinParams) {
547             this.m = m;
548             this.n = startParams.length;
549             this.startParams          = startParams.clone();
550             this.theoreticalMinCost   = theoreticalMinCost;
551             this.theoreticalMinParams = theoreticalMinParams;
552             this.costAccuracy         = 1.0e-8;
553             this.paramsAccuracy       = 1.0e-5;
554         }
555 
556         protected static double[] buildArray(int n, double x) {
557             double[] array = new double[n];
558             Arrays.fill(array, x);
559             return array;
560         }
561 
562         public double[] getTarget() {
563             return buildArray(m, 0.0);
564         }
565 
566         public double[] getWeight() {
567             return buildArray(m, 1.0);
568         }
569 
570         public double[] getStartPoint() {
571             return startParams.clone();
572         }
573 
574         protected void setCostAccuracy(double costAccuracy) {
575             this.costAccuracy = costAccuracy;
576         }
577 
578         protected void setParamsAccuracy(double paramsAccuracy) {
579             this.paramsAccuracy = paramsAccuracy;
580         }
581 
582         public int getN() {
583             return startParams.length;
584         }
585 
586         public void checkTheoreticalMinCost(double rms) {
587             double threshold = costAccuracy * (1.0 + theoreticalMinCost);
588             assertEquals(theoreticalMinCost, FastMath.sqrt(m) * rms, threshold);
589         }
590 
591         public void checkTheoreticalMinParams(double[] params) {
592             if (theoreticalMinParams != null) {
593                 for (int i = 0; i < theoreticalMinParams.length; ++i) {
594                     double mi = theoreticalMinParams[i];
595                     double vi = params[i];
596                     assertEquals(mi, vi, paramsAccuracy * (1.0 + FastMath.abs(mi)));
597                 }
598             }
599         }
600 
601         public MultivariateVectorFunction getModelFunction() {
602             return new MultivariateVectorFunction() {
603                 public double[] value(double[] point) {
604                     return computeValue(point);
605                 }
606             };
607         }
608 
609         public MultivariateMatrixFunction getModelFunctionJacobian() {
610             return new MultivariateMatrixFunction() {
611                 public double[][] value(double[] point) {
612                     return computeJacobian(point);
613                 }
614             };
615         }
616 
617         public abstract double[][] computeJacobian(double[] variables);
618         public abstract double[] computeValue(double[] variables);
619     }
620 
621     private static class LinearFullRankFunction extends MinpackFunction {
622 
623         public LinearFullRankFunction(int m, int n, double x0,
624                                       double theoreticalStartCost,
625                                       double theoreticalMinCost) {
626             super(m, buildArray(n, x0), theoreticalMinCost,
627                   buildArray(n, -1.0));
628         }
629 
630         @Override
631         public double[][] computeJacobian(double[] variables) {
632             double t = 2.0 / m;
633             double[][] jacobian = new double[m][];
634             for (int i = 0; i < m; ++i) {
635                 jacobian[i] = new double[n];
636                 for (int j = 0; j < n; ++j) {
637                     jacobian[i][j] = (i == j) ? (1 - t) : -t;
638                 }
639             }
640             return jacobian;
641         }
642 
643         @Override
644         public double[] computeValue(double[] variables) {
645             double sum = 0;
646             for (int i = 0; i < n; ++i) {
647                 sum += variables[i];
648             }
649             double t  = 1 + 2 * sum / m;
650             double[] f = new double[m];
651             for (int i = 0; i < n; ++i) {
652                 f[i] = variables[i] - t;
653             }
654             Arrays.fill(f, n, m, -t);
655             return f;
656         }
657     }
658 
659     private static class LinearRank1Function extends MinpackFunction {
660 
661         public LinearRank1Function(int m, int n, double x0,
662                                    double theoreticalStartCost,
663                                    double theoreticalMinCost) {
664             super(m, buildArray(n, x0), theoreticalMinCost, null);
665         }
666 
667         @Override
668         public double[][] computeJacobian(double[] variables) {
669             double[][] jacobian = new double[m][];
670             for (int i = 0; i < m; ++i) {
671                 jacobian[i] = new double[n];
672                 for (int j = 0; j < n; ++j) {
673                     jacobian[i][j] = (i + 1) * (j + 1);
674                 }
675             }
676             return jacobian;
677         }
678 
679         @Override
680         public double[] computeValue(double[] variables) {
681             double[] f = new double[m];
682             double sum = 0;
683             for (int i = 0; i < n; ++i) {
684                 sum += (i + 1) * variables[i];
685             }
686             for (int i = 0; i < m; ++i) {
687                 f[i] = (i + 1) * sum - 1;
688             }
689             return f;
690         }
691     }
692 
693     private static class LinearRank1ZeroColsAndRowsFunction extends MinpackFunction {
694 
695         public LinearRank1ZeroColsAndRowsFunction(int m, int n, double x0) {
696             super(m, buildArray(n, x0),
697                   FastMath.sqrt((m * (m + 3) - 6) / (2.0 * (2 * m - 3))),
698                   null);
699         }
700 
701         @Override
702         public double[][] computeJacobian(double[] variables) {
703             double[][] jacobian = new double[m][];
704             for (int i = 0; i < m; ++i) {
705                 jacobian[i] = new double[n];
706                 jacobian[i][0] = 0;
707                 for (int j = 1; j < (n - 1); ++j) {
708                     if (i == 0) {
709                         jacobian[i][j] = 0;
710                     } else if (i != (m - 1)) {
711                         jacobian[i][j] = i * (j + 1);
712                     } else {
713                         jacobian[i][j] = 0;
714                     }
715                 }
716                 jacobian[i][n - 1] = 0;
717             }
718             return jacobian;
719         }
720 
721         @Override
722         public double[] computeValue(double[] variables) {
723             double[] f = new double[m];
724             double sum = 0;
725             for (int i = 1; i < (n - 1); ++i) {
726                 sum += (i + 1) * variables[i];
727             }
728             for (int i = 0; i < (m - 1); ++i) {
729                 f[i] = i * sum - 1;
730             }
731             f[m - 1] = -1;
732             return f;
733         }
734     }
735 
736     private static class RosenbrockFunction extends MinpackFunction {
737         public RosenbrockFunction(double[] startParams, double theoreticalStartCost) {
738             super(2, startParams, 0.0, buildArray(2, 1.0));
739         }
740 
741         @Override
742         public double[][] computeJacobian(double[] variables) {
743             double x1 = variables[0];
744             return new double[][] { { -20 * x1, 10 }, { -1, 0 } };
745         }
746 
747         @Override
748         public double[] computeValue(double[] variables) {
749             double x1 = variables[0];
750             double x2 = variables[1];
751             return new double[] { 10 * (x2 - x1 * x1), 1 - x1 };
752         }
753     }
754 
755     private static class HelicalValleyFunction extends MinpackFunction {
756         public HelicalValleyFunction(double[] startParams,
757                                      double theoreticalStartCost) {
758             super(3, startParams, 0.0, new double[] { 1.0, 0.0, 0.0 });
759         }
760 
761         @Override
762         public double[][] computeJacobian(double[] variables) {
763             double x1 = variables[0];
764             double x2 = variables[1];
765             double tmpSquare = x1 * x1 + x2 * x2;
766             double tmp1 = twoPi * tmpSquare;
767             double tmp2 = FastMath.sqrt(tmpSquare);
768             return new double[][] {
769                 {  100 * x2 / tmp1, -100 * x1 / tmp1, 10 },
770                 { 10 * x1 / tmp2, 10 * x2 / tmp2, 0 },
771                 { 0, 0, 1 }
772             };
773         }
774 
775         @Override
776         public double[] computeValue(double[] variables) {
777             double x1 = variables[0];
778             double x2 = variables[1];
779             double x3 = variables[2];
780             double tmp1;
781             if (x1 == 0) {
782                 tmp1 = (x2 >= 0) ? 0.25 : -0.25;
783             } else {
784                 tmp1 = FastMath.atan(x2 / x1) / twoPi;
785                 if (x1 < 0) {
786                     tmp1 += 0.5;
787                 }
788             }
789             double tmp2 = FastMath.sqrt(x1 * x1 + x2 * x2);
790             return new double[] {
791                 10.0 * (x3 - 10 * tmp1),
792                 10.0 * (tmp2 - 1),
793                 x3
794             };
795         }
796 
797         private static final double twoPi = 2.0 * FastMath.PI;
798     }
799 
800     private static class PowellSingularFunction extends MinpackFunction {
801 
802         public PowellSingularFunction(double[] startParams,
803                                       double theoreticalStartCost) {
804             super(4, startParams, 0.0, buildArray(4, 0.0));
805         }
806 
807         @Override
808         public double[][] computeJacobian(double[] variables) {
809             double x1 = variables[0];
810             double x2 = variables[1];
811             double x3 = variables[2];
812             double x4 = variables[3];
813             return new double[][] {
814                 { 1, 10, 0, 0 },
815                 { 0, 0, sqrt5, -sqrt5 },
816                 { 0, 2 * (x2 - 2 * x3), -4 * (x2 - 2 * x3), 0 },
817                 { 2 * sqrt10 * (x1 - x4), 0, 0, -2 * sqrt10 * (x1 - x4) }
818             };
819         }
820 
821         @Override
822         public double[] computeValue(double[] variables) {
823             double x1 = variables[0];
824             double x2 = variables[1];
825             double x3 = variables[2];
826             double x4 = variables[3];
827             return new double[] {
828                 x1 + 10 * x2,
829                 sqrt5 * (x3 - x4),
830                 (x2 - 2 * x3) * (x2 - 2 * x3),
831                 sqrt10 * (x1 - x4) * (x1 - x4)
832             };
833         }
834 
835         private static final double sqrt5  = FastMath.sqrt( 5.0);
836         private static final double sqrt10 = FastMath.sqrt(10.0);
837   }
838 
839     private static class FreudensteinRothFunction extends MinpackFunction {
840 
841         public FreudensteinRothFunction(double[] startParams,
842                                         double theoreticalStartCost,
843                                         double theoreticalMinCost,
844                                         double[] theoreticalMinParams) {
845             super(2, startParams, theoreticalMinCost,
846                   theoreticalMinParams);
847         }
848 
849         @Override
850         public double[][] computeJacobian(double[] variables) {
851             double x2 = variables[1];
852             return new double[][] {
853                 { 1, x2 * (10 - 3 * x2) -  2 },
854                 { 1, x2 * ( 2 + 3 * x2) - 14, }
855             };
856         }
857 
858         @Override
859         public double[] computeValue(double[] variables) {
860             double x1 = variables[0];
861             double x2 = variables[1];
862             return new double[] {
863                 -13.0 + x1 + ((5.0 - x2) * x2 -  2.0) * x2,
864                 -29.0 + x1 + ((1.0 + x2) * x2 - 14.0) * x2
865             };
866         }
867     }
868 
869     private static class BardFunction extends MinpackFunction {
870 
871         public BardFunction(double x0,
872                             double theoreticalStartCost,
873                             double theoreticalMinCost,
874                             double[] theoreticalMinParams) {
875             super(15, buildArray(3, x0), theoreticalMinCost,
876                   theoreticalMinParams);
877         }
878 
879         @Override
880         public double[][] computeJacobian(double[] variables) {
881             double   x2 = variables[1];
882             double   x3 = variables[2];
883             double[][] jacobian = new double[m][];
884             for (int i = 0; i < m; ++i) {
885                 double tmp1 = i  + 1;
886                 double tmp2 = 15 - i;
887                 double tmp3 = (i <= 7) ? tmp1 : tmp2;
888                 double tmp4 = x2 * tmp2 + x3 * tmp3;
889                 tmp4 *= tmp4;
890                 jacobian[i] = new double[] { -1, tmp1 * tmp2 / tmp4, tmp1 * tmp3 / tmp4 };
891             }
892             return jacobian;
893         }
894 
895         @Override
896         public double[] computeValue(double[] variables) {
897             double   x1 = variables[0];
898             double   x2 = variables[1];
899             double   x3 = variables[2];
900             double[] f = new double[m];
901             for (int i = 0; i < m; ++i) {
902                 double tmp1 = i + 1;
903                 double tmp2 = 15 - i;
904                 double tmp3 = (i <= 7) ? tmp1 : tmp2;
905                 f[i] = y[i] - (x1 + tmp1 / (x2 * tmp2 + x3 * tmp3));
906             }
907             return f;
908         }
909 
910         private static final double[] y = {
911             0.14, 0.18, 0.22, 0.25, 0.29,
912             0.32, 0.35, 0.39, 0.37, 0.58,
913             0.73, 0.96, 1.34, 2.10, 4.39
914         };
915     }
916 
917     private static class KowalikOsborneFunction extends MinpackFunction {
918 
919         public KowalikOsborneFunction(double[] startParams,
920                                       double theoreticalStartCost,
921                                       double theoreticalMinCost,
922                                       double[] theoreticalMinParams) {
923             super(11, startParams, theoreticalMinCost,
924                   theoreticalMinParams);
925             if (theoreticalStartCost > 20.0) {
926                 setCostAccuracy(2.0e-4);
927                 setParamsAccuracy(5.0e-3);
928             }
929         }
930 
931         @Override
932         public double[][] computeJacobian(double[] variables) {
933             double   x1 = variables[0];
934             double   x2 = variables[1];
935             double   x3 = variables[2];
936             double   x4 = variables[3];
937             double[][] jacobian = new double[m][];
938             for (int i = 0; i < m; ++i) {
939                 double tmp = v[i] * (v[i] + x3) + x4;
940                 double j1  = -v[i] * (v[i] + x2) / tmp;
941                 double j2  = -v[i] * x1 / tmp;
942                 double j3  = j1 * j2;
943                 double j4  = j3 / v[i];
944                 jacobian[i] = new double[] { j1, j2, j3, j4 };
945             }
946             return jacobian;
947         }
948 
949         @Override
950         public double[] computeValue(double[] variables) {
951             double x1 = variables[0];
952             double x2 = variables[1];
953             double x3 = variables[2];
954             double x4 = variables[3];
955             double[] f = new double[m];
956             for (int i = 0; i < m; ++i) {
957                 f[i] = y[i] - x1 * (v[i] * (v[i] + x2)) / (v[i] * (v[i] + x3) + x4);
958             }
959             return f;
960         }
961 
962         private static final double[] v = {
963             4.0, 2.0, 1.0, 0.5, 0.25, 0.167, 0.125, 0.1, 0.0833, 0.0714, 0.0625
964         };
965 
966         private static final double[] y = {
967             0.1957, 0.1947, 0.1735, 0.1600, 0.0844, 0.0627,
968             0.0456, 0.0342, 0.0323, 0.0235, 0.0246
969         };
970     }
971 
972     private static class MeyerFunction extends MinpackFunction {
973 
974         public MeyerFunction(double[] startParams,
975                              double theoreticalStartCost,
976                              double theoreticalMinCost,
977                              double[] theoreticalMinParams) {
978             super(16, startParams, theoreticalMinCost,
979                   theoreticalMinParams);
980             if (theoreticalStartCost > 1.0e6) {
981                 setCostAccuracy(7.0e-3);
982                 setParamsAccuracy(2.0e-2);
983             }
984         }
985 
986         @Override
987         public double[][] computeJacobian(double[] variables) {
988             double   x1 = variables[0];
989             double   x2 = variables[1];
990             double   x3 = variables[2];
991             double[][] jacobian = new double[m][];
992             for (int i = 0; i < m; ++i) {
993                 double temp = 5.0 * (i + 1) + 45.0 + x3;
994                 double tmp1 = x2 / temp;
995                 double tmp2 = FastMath.exp(tmp1);
996                 double tmp3 = x1 * tmp2 / temp;
997                 jacobian[i] = new double[] { tmp2, tmp3, -tmp1 * tmp3 };
998             }
999             return jacobian;
1000         }
1001 
1002         @Override
1003         public double[] computeValue(double[] variables) {
1004             double x1 = variables[0];
1005             double x2 = variables[1];
1006             double x3 = variables[2];
1007             double[] f = new double[m];
1008             for (int i = 0; i < m; ++i) {
1009                 f[i] = x1 * FastMath.exp(x2 / (5.0 * (i + 1) + 45.0 + x3)) - y[i];
1010             }
1011             return f;
1012         }
1013 
1014         private static final double[] y = {
1015             34780.0, 28610.0, 23650.0, 19630.0,
1016             16370.0, 13720.0, 11540.0,  9744.0,
1017             8261.0,  7030.0,  6005.0,  5147.0,
1018             4427.0,  3820.0,  3307.0,  2872.0
1019         };
1020     }
1021 
1022     private static class WatsonFunction extends MinpackFunction {
1023 
1024         public WatsonFunction(int n, double x0,
1025                               double theoreticalStartCost,
1026                               double theoreticalMinCost,
1027                               double[] theoreticalMinParams) {
1028             super(31, buildArray(n, x0), theoreticalMinCost,
1029                   theoreticalMinParams);
1030         }
1031 
1032         @Override
1033         public double[][] computeJacobian(double[] variables) {
1034             double[][] jacobian = new double[m][];
1035 
1036             for (int i = 0; i < (m - 2); ++i) {
1037                 double div = (i + 1) / 29.0;
1038                 double s2  = 0.0;
1039                 double dx  = 1.0;
1040                 for (int j = 0; j < n; ++j) {
1041                     s2 += dx * variables[j];
1042                     dx *= div;
1043                 }
1044                 double temp= 2 * div * s2;
1045                 dx = 1.0 / div;
1046                 jacobian[i] = new double[n];
1047                 for (int j = 0; j < n; ++j) {
1048                     jacobian[i][j] = dx * (j - temp);
1049                     dx *= div;
1050                 }
1051             }
1052 
1053             jacobian[m - 2]    = new double[n];
1054             jacobian[m - 2][0] = 1;
1055 
1056             jacobian[m - 1]   = new double[n];
1057             jacobian[m - 1][0]= -2 * variables[0];
1058             jacobian[m - 1][1]= 1;
1059 
1060             return jacobian;
1061         }
1062 
1063         @Override
1064         public double[] computeValue(double[] variables) {
1065             double[] f = new double[m];
1066             for (int i = 0; i < (m - 2); ++i) {
1067                 double div = (i + 1) / 29.0;
1068                 double s1 = 0;
1069                 double dx = 1;
1070                 for (int j = 1; j < n; ++j) {
1071                     s1 += j * dx * variables[j];
1072                     dx *= div;
1073                 }
1074                 double s2 = 0;
1075                 dx = 1;
1076                 for (int j = 0; j < n; ++j) {
1077                     s2 += dx * variables[j];
1078                     dx *= div;
1079                 }
1080                 f[i] = s1 - s2 * s2 - 1;
1081             }
1082 
1083             double x1 = variables[0];
1084             double x2 = variables[1];
1085             f[m - 2] = x1;
1086             f[m - 1] = x2 - x1 * x1 - 1;
1087 
1088             return f;
1089         }
1090     }
1091 
1092     private static class Box3DimensionalFunction extends MinpackFunction {
1093 
1094         public Box3DimensionalFunction(int m, double[] startParams,
1095                                        double theoreticalStartCost) {
1096             super(m, startParams, 0.0,
1097                   new double[] { 1.0, 10.0, 1.0 });
1098         }
1099 
1100         @Override
1101         public double[][] computeJacobian(double[] variables) {
1102             double   x1 = variables[0];
1103             double   x2 = variables[1];
1104             double[][] jacobian = new double[m][];
1105             for (int i = 0; i < m; ++i) {
1106                 double tmp = (i + 1) / 10.0;
1107                 jacobian[i] = new double[] {
1108                     -tmp * FastMath.exp(-tmp * x1),
1109                     tmp * FastMath.exp(-tmp * x2),
1110                     FastMath.exp(-i - 1) - FastMath.exp(-tmp)
1111                 };
1112             }
1113             return jacobian;
1114         }
1115 
1116         @Override
1117         public double[] computeValue(double[] variables) {
1118             double x1 = variables[0];
1119             double x2 = variables[1];
1120             double x3 = variables[2];
1121             double[] f = new double[m];
1122             for (int i = 0; i < m; ++i) {
1123                 double tmp = (i + 1) / 10.0;
1124                 f[i] = FastMath.exp(-tmp * x1) - FastMath.exp(-tmp * x2)
1125                     + (FastMath.exp(-i - 1) - FastMath.exp(-tmp)) * x3;
1126             }
1127             return f;
1128         }
1129     }
1130 
1131     private static class JennrichSampsonFunction extends MinpackFunction {
1132 
1133         public JennrichSampsonFunction(int m, double[] startParams,
1134                                        double theoreticalStartCost,
1135                                        double theoreticalMinCost,
1136                                        double[] theoreticalMinParams) {
1137             super(m, startParams, theoreticalMinCost,
1138                   theoreticalMinParams);
1139         }
1140 
1141         @Override
1142         public double[][] computeJacobian(double[] variables) {
1143             double   x1 = variables[0];
1144             double   x2 = variables[1];
1145             double[][] jacobian = new double[m][];
1146             for (int i = 0; i < m; ++i) {
1147                 double t = i + 1;
1148                 jacobian[i] = new double[] { -t * FastMath.exp(t * x1), -t * FastMath.exp(t * x2) };
1149             }
1150             return jacobian;
1151         }
1152 
1153         @Override
1154         public double[] computeValue(double[] variables) {
1155             double x1 = variables[0];
1156             double x2 = variables[1];
1157             double[] f = new double[m];
1158             for (int i = 0; i < m; ++i) {
1159                 double temp = i + 1;
1160                 f[i] = 2 + 2 * temp - FastMath.exp(temp * x1) - FastMath.exp(temp * x2);
1161             }
1162             return f;
1163         }
1164     }
1165 
1166     private static class BrownDennisFunction extends MinpackFunction {
1167 
1168         public BrownDennisFunction(int m, double[] startParams,
1169                                    double theoreticalStartCost,
1170                                    double theoreticalMinCost,
1171                                    double[] theoreticalMinParams) {
1172             super(m, startParams, theoreticalMinCost,
1173                 theoreticalMinParams);
1174             setCostAccuracy(2.5e-8);
1175         }
1176 
1177         @Override
1178         public double[][] computeJacobian(double[] variables) {
1179             double   x1 = variables[0];
1180             double   x2 = variables[1];
1181             double   x3 = variables[2];
1182             double   x4 = variables[3];
1183             double[][] jacobian = new double[m][];
1184             for (int i = 0; i < m; ++i) {
1185                 double temp = (i + 1) / 5.0;
1186                 double ti   = FastMath.sin(temp);
1187                 double tmp1 = x1 + temp * x2 - FastMath.exp(temp);
1188                 double tmp2 = x3 + ti   * x4 - FastMath.cos(temp);
1189                 jacobian[i] = new double[] {
1190                     2 * tmp1, 2 * temp * tmp1, 2 * tmp2, 2 * ti * tmp2
1191                 };
1192             }
1193             return jacobian;
1194         }
1195 
1196         @Override
1197         public double[] computeValue(double[] variables) {
1198             double x1 = variables[0];
1199             double x2 = variables[1];
1200             double x3 = variables[2];
1201             double x4 = variables[3];
1202             double[] f = new double[m];
1203             for (int i = 0; i < m; ++i) {
1204                 double temp = (i + 1) / 5.0;
1205                 double tmp1 = x1 + temp * x2 - FastMath.exp(temp);
1206                 double tmp2 = x3 + FastMath.sin(temp) * x4 - FastMath.cos(temp);
1207                 f[i] = tmp1 * tmp1 + tmp2 * tmp2;
1208             }
1209             return f;
1210         }
1211     }
1212 
1213     private static class ChebyquadFunction extends MinpackFunction {
1214 
1215         private static double[] buildChebyquadArray(int n, double factor) {
1216             double[] array = new double[n];
1217             double inv = factor / (n + 1);
1218             for (int i = 0; i < n; ++i) {
1219                 array[i] = (i + 1) * inv;
1220             }
1221             return array;
1222         }
1223 
1224         public ChebyquadFunction(int n, int m, double factor,
1225                                  double theoreticalStartCost,
1226                                  double theoreticalMinCost,
1227                                  double[] theoreticalMinParams) {
1228             super(m, buildChebyquadArray(n, factor), theoreticalMinCost,
1229                   theoreticalMinParams);
1230         }
1231 
1232         @Override
1233         public double[][] computeJacobian(double[] variables) {
1234             double[][] jacobian = new double[m][];
1235             for (int i = 0; i < m; ++i) {
1236                 jacobian[i] = new double[n];
1237             }
1238 
1239             double dx = 1.0 / n;
1240             for (int j = 0; j < n; ++j) {
1241                 double tmp1 = 1;
1242                 double tmp2 = 2 * variables[j] - 1;
1243                 double temp = 2 * tmp2;
1244                 double tmp3 = 0;
1245                 double tmp4 = 2;
1246                 for (int i = 0; i < m; ++i) {
1247                     jacobian[i][j] = dx * tmp4;
1248                     double ti = 4 * tmp2 + temp * tmp4 - tmp3;
1249                     tmp3 = tmp4;
1250                     tmp4 = ti;
1251                     ti   = temp * tmp2 - tmp1;
1252                     tmp1 = tmp2;
1253                     tmp2 = ti;
1254                 }
1255             }
1256 
1257             return jacobian;
1258         }
1259 
1260         @Override
1261         public double[] computeValue(double[] variables) {
1262             double[] f = new double[m];
1263 
1264             for (int j = 0; j < n; ++j) {
1265                 double tmp1 = 1;
1266                 double tmp2 = 2 * variables[j] - 1;
1267                 double temp = 2 * tmp2;
1268                 for (int i = 0; i < m; ++i) {
1269                     f[i] += tmp2;
1270                     double ti = temp * tmp2 - tmp1;
1271                     tmp1 = tmp2;
1272                     tmp2 = ti;
1273                 }
1274             }
1275 
1276             double dx = 1.0 / n;
1277             boolean iev = false;
1278             for (int i = 0; i < m; ++i) {
1279                 f[i] *= dx;
1280                 if (iev) {
1281                     f[i] += 1.0 / (i * (i + 2));
1282                 }
1283                 iev = ! iev;
1284             }
1285 
1286             return f;
1287         }
1288     }
1289 
1290     private static class BrownAlmostLinearFunction extends MinpackFunction {
1291 
1292         public BrownAlmostLinearFunction(int m, double factor,
1293                                          double theoreticalStartCost,
1294                                          double theoreticalMinCost,
1295                                          double[] theoreticalMinParams) {
1296             super(m, buildArray(m, factor), theoreticalMinCost,
1297                   theoreticalMinParams);
1298         }
1299 
1300         @Override
1301         public double[][] computeJacobian(double[] variables) {
1302             double[][] jacobian = new double[m][];
1303             for (int i = 0; i < m; ++i) {
1304                 jacobian[i] = new double[n];
1305             }
1306 
1307             double prod = 1;
1308             for (int j = 0; j < n; ++j) {
1309                 prod *= variables[j];
1310                 for (int i = 0; i < n; ++i) {
1311                     jacobian[i][j] = 1;
1312                 }
1313                 jacobian[j][j] = 2;
1314             }
1315 
1316             for (int j = 0; j < n; ++j) {
1317                 double temp = variables[j];
1318                 if (temp == 0) {
1319                     temp = 1;
1320                     prod = 1;
1321                     for (int k = 0; k < n; ++k) {
1322                         if (k != j) {
1323                             prod *= variables[k];
1324                         }
1325                     }
1326                 }
1327                 jacobian[n - 1][j] = prod / temp;
1328             }
1329 
1330             return jacobian;
1331         }
1332 
1333         @Override
1334         public double[] computeValue(double[] variables) {
1335             double[] f = new double[m];
1336             double sum  = -(n + 1);
1337             double prod = 1;
1338             for (int j = 0; j < n; ++j) {
1339                 sum  += variables[j];
1340                 prod *= variables[j];
1341             }
1342             for (int i = 0; i < n; ++i) {
1343                 f[i] = variables[i] + sum;
1344             }
1345             f[n - 1] = prod - 1;
1346             return f;
1347         }
1348     }
1349 
1350     private static class Osborne1Function extends MinpackFunction {
1351 
1352         public Osborne1Function(double[] startParams,
1353                                 double theoreticalStartCost,
1354                                 double theoreticalMinCost,
1355                                 double[] theoreticalMinParams) {
1356             super(33, startParams, theoreticalMinCost,
1357                   theoreticalMinParams);
1358         }
1359 
1360         @Override
1361         public double[][] computeJacobian(double[] variables) {
1362             double   x2 = variables[1];
1363             double   x3 = variables[2];
1364             double   x4 = variables[3];
1365             double   x5 = variables[4];
1366             double[][] jacobian = new double[m][];
1367             for (int i = 0; i < m; ++i) {
1368                 double temp = 10.0 * i;
1369                 double tmp1 = FastMath.exp(-temp * x4);
1370                 double tmp2 = FastMath.exp(-temp * x5);
1371                 jacobian[i] = new double[] {
1372                     -1, -tmp1, -tmp2, temp * x2 * tmp1, temp * x3 * tmp2
1373                 };
1374             }
1375             return jacobian;
1376         }
1377 
1378         @Override
1379         public double[] computeValue(double[] variables) {
1380             double x1 = variables[0];
1381             double x2 = variables[1];
1382             double x3 = variables[2];
1383             double x4 = variables[3];
1384             double x5 = variables[4];
1385             double[] f = new double[m];
1386             for (int i = 0; i < m; ++i) {
1387                 double temp = 10.0 * i;
1388                 double tmp1 = FastMath.exp(-temp * x4);
1389                 double tmp2 = FastMath.exp(-temp * x5);
1390                 f[i] = y[i] - (x1 + x2 * tmp1 + x3 * tmp2);
1391             }
1392             return f;
1393         }
1394 
1395         private static final double[] y = {
1396             0.844, 0.908, 0.932, 0.936, 0.925, 0.908, 0.881, 0.850, 0.818, 0.784, 0.751,
1397             0.718, 0.685, 0.658, 0.628, 0.603, 0.580, 0.558, 0.538, 0.522, 0.506, 0.490,
1398             0.478, 0.467, 0.457, 0.448, 0.438, 0.431, 0.424, 0.420, 0.414, 0.411, 0.406
1399         };
1400     }
1401 
1402     private static class Osborne2Function extends MinpackFunction {
1403 
1404         public Osborne2Function(double[] startParams,
1405                                 double theoreticalStartCost,
1406                                 double theoreticalMinCost,
1407                                 double[] theoreticalMinParams) {
1408             super(65, startParams, theoreticalMinCost,
1409                   theoreticalMinParams);
1410         }
1411 
1412         @Override
1413         public double[][] computeJacobian(double[] variables) {
1414             double   x01 = variables[0];
1415             double   x02 = variables[1];
1416             double   x03 = variables[2];
1417             double   x04 = variables[3];
1418             double   x05 = variables[4];
1419             double   x06 = variables[5];
1420             double   x07 = variables[6];
1421             double   x08 = variables[7];
1422             double   x09 = variables[8];
1423             double   x10 = variables[9];
1424             double   x11 = variables[10];
1425             double[][] jacobian = new double[m][];
1426             for (int i = 0; i < m; ++i) {
1427                 double temp = i / 10.0;
1428                 double tmp1 = FastMath.exp(-x05 * temp);
1429                 double tmp2 = FastMath.exp(-x06 * (temp - x09) * (temp - x09));
1430                 double tmp3 = FastMath.exp(-x07 * (temp - x10) * (temp - x10));
1431                 double tmp4 = FastMath.exp(-x08 * (temp - x11) * (temp - x11));
1432                 jacobian[i] = new double[] {
1433                     -tmp1,
1434                     -tmp2,
1435                     -tmp3,
1436                     -tmp4,
1437                     temp * x01 * tmp1,
1438                     x02 * (temp - x09) * (temp - x09) * tmp2,
1439                     x03 * (temp - x10) * (temp - x10) * tmp3,
1440                     x04 * (temp - x11) * (temp - x11) * tmp4,
1441                     -2 * x02 * x06 * (temp - x09) * tmp2,
1442                     -2 * x03 * x07 * (temp - x10) * tmp3,
1443                     -2 * x04 * x08 * (temp - x11) * tmp4
1444                 };
1445             }
1446             return jacobian;
1447         }
1448 
1449         @Override
1450         public double[] computeValue(double[] variables) {
1451             double x01 = variables[0];
1452             double x02 = variables[1];
1453             double x03 = variables[2];
1454             double x04 = variables[3];
1455             double x05 = variables[4];
1456             double x06 = variables[5];
1457             double x07 = variables[6];
1458             double x08 = variables[7];
1459             double x09 = variables[8];
1460             double x10 = variables[9];
1461             double x11 = variables[10];
1462             double[] f = new double[m];
1463             for (int i = 0; i < m; ++i) {
1464                 double temp = i / 10.0;
1465                 double tmp1 = FastMath.exp(-x05 * temp);
1466                 double tmp2 = FastMath.exp(-x06 * (temp - x09) * (temp - x09));
1467                 double tmp3 = FastMath.exp(-x07 * (temp - x10) * (temp - x10));
1468                 double tmp4 = FastMath.exp(-x08 * (temp - x11) * (temp - x11));
1469                 f[i] = y[i] - (x01 * tmp1 + x02 * tmp2 + x03 * tmp3 + x04 * tmp4);
1470             }
1471             return f;
1472         }
1473 
1474         private static final double[] y = {
1475             1.366, 1.191, 1.112, 1.013, 0.991,
1476             0.885, 0.831, 0.847, 0.786, 0.725,
1477             0.746, 0.679, 0.608, 0.655, 0.616,
1478             0.606, 0.602, 0.626, 0.651, 0.724,
1479             0.649, 0.649, 0.694, 0.644, 0.624,
1480             0.661, 0.612, 0.558, 0.533, 0.495,
1481             0.500, 0.423, 0.395, 0.375, 0.372,
1482             0.391, 0.396, 0.405, 0.428, 0.429,
1483             0.523, 0.562, 0.607, 0.653, 0.672,
1484             0.708, 0.633, 0.668, 0.645, 0.632,
1485             0.591, 0.559, 0.597, 0.625, 0.739,
1486             0.710, 0.729, 0.720, 0.636, 0.581,
1487             0.428, 0.292, 0.162, 0.098, 0.054
1488         };
1489     }
1490 }