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