Lines Matching +full:- +full:lm
22 --g_test_level; \
35 static const double y[15]={1.4e-1, 1.8e-1, 2.2e-1, 2.5e-1, 2.9e-1, 3.2e-1, 3.5e-1, in fcn_chkder()
36 3.9e-1, 3.7e-1, 5.8e-1, 7.3e-1, 9.6e-1, 1.34, 2.1, 4.39}; in fcn_chkder()
45 tmp2 = 16-i-1; in fcn_chkder()
48 fvec[i] = y[i] - (x[0] + tmp1/(x[1]*tmp2 + x[2]*tmp3)); in fcn_chkder()
53 tmp2 = 16-i-1; in fcn_chkder()
61 fjac(i,0) = -1.; in fcn_chkder()
79 x << 9.2e-1, 1.3e-1, 5.4e-1; in testChkder()
87 fvecp -= fvec; in testChkder()
92 -1.181606, -1.429655, -1.606344, in testChkder()
93 -1.745269, -1.840654, -1.921586, in testChkder()
94 -1.984141, -2.022537, -2.468977, in testChkder()
95 -2.827562, -3.473582, -4.437612, in testChkder()
96 -6.047662, -9.267761, -18.91806; in testChkder()
98 -7.724666e-09, -3.432406e-09, -2.034843e-10, in testChkder()
99 2.313685e-09, 4.331078e-09, 5.984096e-09, in testChkder()
100 7.363281e-09, 8.53147e-09, 1.488591e-08, in testChkder()
101 2.33585e-08, 3.522012e-08, 5.301255e-08, in testChkder()
102 8.26666e-08, 1.419747e-07, 3.19899e-07; in testChkder()
146 static const double y[15] = {1.4e-1, 1.8e-1, 2.2e-1, 2.5e-1, 2.9e-1, 3.2e-1, 3.5e-1, in operator ()()
147 3.9e-1, 3.7e-1, 5.8e-1, 7.3e-1, 9.6e-1, 1.34, 2.1, 4.39}; in operator ()()
152 tmp2 = 16 - i - 1; in operator ()()
154 fvec[i] = y[i] - (x[0] + tmp1/(x[1]*tmp2 + x[2]*tmp3)); in operator ()()
165 tmp2 = 16 - i - 1; in df()
168 fjac(i,0) = -1; in df()
187 LevenbergMarquardt<lmder_functor> lm(functor); in testLmder1() local
188 info = lm.lmder1(x); in testLmder1()
192 LM_CHECK_N_ITERS(lm, 6, 5); in testLmder1()
195 VERIFY_IS_APPROX(lm.fvec.blueNorm(), 0.09063596); in testLmder1()
215 LevenbergMarquardt<lmder_functor> lm(functor); in testLmder() local
216 info = lm.minimize(x); in testLmder()
220 LM_CHECK_N_ITERS(lm, 6, 5); in testLmder()
223 fnorm = lm.fvec.blueNorm(); in testLmder()
232 covfac = fnorm*fnorm/(m-n); in testLmder()
233 internal::covar(lm.fjac, lm.permutation.indices()); // TODO : move this as a function of lm in testLmder()
237 0.0001531202, 0.002869941, -0.002656662, in testLmder()
238 0.002869941, 0.09480935, -0.09098995, in testLmder()
239 -0.002656662, -0.09098995, 0.08778727; in testLmder()
244 cov = covfac*lm.fjac.topLeftCorner<n,n>(); in testLmder()
261 temp = (3. - 2.*x[k])*x[k]; in operator ()()
263 if (k) temp1 = x[k-1]; in operator ()()
265 if (k != n-1) temp2 = x[k+1]; in operator ()()
266 fvec[k] = temp - temp1 - 2.*temp2 + 1.; in operator ()()
279 fjac(k,k) = 3.- 4.*x[k]; in df()
280 if (k) fjac(k,k-1) = -1.; in df()
281 if (k != n-1) fjac(k,k+1) = -2.; in df()
295 x.setConstant(n, -1.); in testHybrj1()
307 VERIFY_IS_APPROX(solver.fvec.blueNorm(), 1.192636e-08); in testHybrj1()
313 -0.5706545, -0.6816283, -0.7017325, in testHybrj1()
314 -0.7042129, -0.701369, -0.6918656, in testHybrj1()
315 -0.665792, -0.5960342, -0.4164121; in testHybrj1()
326 x.setConstant(n, -1.); in testHybrj()
341 VERIFY_IS_APPROX(solver.fvec.blueNorm(), 1.192636e-08); in testHybrj()
347 -0.5706545, -0.6816283, -0.7017325, in testHybrj()
348 -0.7042129, -0.701369, -0.6918656, in testHybrj()
349 -0.665792, -0.5960342, -0.4164121; in testHybrj()
365 temp = (3. - 2.*x[k])*x[k]; in operator ()()
367 if (k) temp1 = x[k-1]; in operator ()()
369 if (k != n-1) temp2 = x[k+1]; in operator ()()
370 fvec[k] = temp - temp1 - 2.*temp2 + 1.; in operator ()()
382 x.setConstant(n, -1.); in testHybrd1()
394 VERIFY_IS_APPROX(solver.fvec.blueNorm(), 1.192636e-08); in testHybrd1()
398 …x_ref << -0.5706545, -0.6816283, -0.7017325, -0.7042129, -0.701369, -0.6918656, -0.665792, -0.5960… in testHybrd1()
409 x.setConstant(n, -1.); in testHybrd()
425 VERIFY_IS_APPROX(solver.fvec.blueNorm(), 1.192636e-08); in testHybrd()
430 -0.5706545, -0.6816283, -0.7017325, in testHybrd()
431 -0.7042129, -0.701369, -0.6918656, in testHybrd()
432 -0.665792, -0.5960342, -0.4164121; in testHybrd()
443 static const double y[15]={1.4e-1, 1.8e-1, 2.2e-1, 2.5e-1, 2.9e-1, 3.2e-1, 3.5e-1, in operator ()()
444 3.9e-1, 3.7e-1, 5.8e-1, 7.3e-1, 9.6e-1, 1.34, 2.1, 4.39}; in operator ()()
452 tmp2 = 16 - i - 1; in operator ()()
454 fvec[i] = y[i] - (x[0] + tmp1/(x[1]*tmp2 + x[2]*tmp3)); in operator ()()
464 VectorXd::Index i = rownb-2; in df()
466 tmp2 = 16 - i - 1; in df()
469 jac_row[0] = -1; in df()
488 LevenbergMarquardt<lmstr_functor> lm(functor); in testLmstr1() local
489 info = lm.lmstr1(x); in testLmstr1()
493 LM_CHECK_N_ITERS(lm, 6, 5); in testLmstr1()
496 VERIFY_IS_APPROX(lm.fvec.blueNorm(), 0.09063596); in testLmstr1()
516 LevenbergMarquardt<lmstr_functor> lm(functor); in testLmstr() local
517 info = lm.minimizeOptimumStorage(x); in testLmstr()
521 LM_CHECK_N_ITERS(lm, 6, 5); in testLmstr()
524 fnorm = lm.fvec.blueNorm(); in testLmstr()
541 static const double y[15]={1.4e-1,1.8e-1,2.2e-1,2.5e-1,2.9e-1,3.2e-1,3.5e-1,3.9e-1, in operator ()()
542 3.7e-1,5.8e-1,7.3e-1,9.6e-1,1.34e0,2.1e0,4.39e0}; in operator ()()
549 tmp2 = 15 - i; in operator ()()
553 fvec[i] = y[i] - (x[0] + tmp1/(x[1]*tmp2 + x[2]*tmp3)); in operator ()()
571 DenseIndex nfev = -1; // initialize to avoid maybe-uninitialized warning in testLmdif1()
602 LevenbergMarquardt<NumericalDiff<lmdif_functor> > lm(numDiff); in testLmdif() local
603 info = lm.minimize(x); in testLmdif()
607 VERIFY_IS_EQUAL(lm.nfev, 26); in testLmdif()
610 fnorm = lm.fvec.blueNorm(); in testLmdif()
619 covfac = fnorm*fnorm/(m-n); in testLmdif()
620 internal::covar(lm.fjac, lm.permutation.indices()); // TODO : move this as a function of lm in testLmdif()
624 0.0001531202, 0.002869942, -0.002656662, in testLmdif()
625 0.002869942, 0.09480937, -0.09098997, in testLmdif()
626 -0.002656662, -0.09098997, 0.08778729; in testLmdif()
631 cov = covfac*lm.fjac.topLeftCorner<n,n>(); in testLmdif()
650 fvec[i] = exp(-b[0]*x)/(b[1]+b[2]*x) - m_y[i]; in operator ()()
662 double e = exp(-b[0]*x); in df()
663 fjac(i,0) = -x*e*factor; in df()
664 fjac(i,1) = -e*factor*factor; in df()
665 fjac(i,2) = -x*e*factor*factor; in df()
687 LevenbergMarquardt<chwirut2_functor> lm(functor); in testNistChwirut2() local
688 info = lm.minimize(x); in testNistChwirut2()
692 LM_CHECK_N_ITERS(lm, 10, 8); in testNistChwirut2()
694 VERIFY_IS_APPROX(lm.fvec.squaredNorm(), 5.1304802941E+02); in testNistChwirut2()
696 VERIFY_IS_APPROX(x[0], 1.6657666537E-01); in testNistChwirut2()
697 VERIFY_IS_APPROX(x[1], 5.1653291286E-03); in testNistChwirut2()
698 VERIFY_IS_APPROX(x[2], 1.2150007096E-02); in testNistChwirut2()
705 lm.resetParameters(); in testNistChwirut2()
706 lm.parameters.ftol = 1.E6*NumTraits<double>::epsilon(); in testNistChwirut2()
707 lm.parameters.xtol = 1.E6*NumTraits<double>::epsilon(); in testNistChwirut2()
708 info = lm.minimize(x); in testNistChwirut2()
712 LM_CHECK_N_ITERS(lm, 7, 6); in testNistChwirut2()
714 VERIFY_IS_APPROX(lm.fvec.squaredNorm(), 5.1304802941E+02); in testNistChwirut2()
716 VERIFY_IS_APPROX(x[0], 1.6657666537E-01); in testNistChwirut2()
717 VERIFY_IS_APPROX(x[1], 5.1653291286E-03); in testNistChwirut2()
718 VERIFY_IS_APPROX(x[2], 1.2150007096E-02); in testNistChwirut2()
732 fvec[i] = b[0]*(1.-exp(-b[1]*m_x[i])) - m_y[i] ; in operator ()()
742 fjac(i,0) = (1.-exp(-b[1]*m_x[i])); in df()
743 fjac(i,1) = (b[0]*m_x[i]*exp(-b[1]*m_x[i])); in df()
765 LevenbergMarquardt<misra1a_functor> lm(functor); in testNistMisra1a() local
766 info = lm.minimize(x); in testNistMisra1a()
770 LM_CHECK_N_ITERS(lm, 19, 15); in testNistMisra1a()
772 VERIFY_IS_APPROX(lm.fvec.squaredNorm(), 1.2455138894E-01); in testNistMisra1a()
775 VERIFY_IS_APPROX(x[1], 5.5015643181E-04); in testNistMisra1a()
782 info = lm.minimize(x); in testNistMisra1a()
786 LM_CHECK_N_ITERS(lm, 5, 4); in testNistMisra1a()
788 VERIFY_IS_APPROX(lm.fvec.squaredNorm(), 1.2455138894E-01); in testNistMisra1a()
791 VERIFY_IS_APPROX(x[1], 5.5015643181E-04); in testNistMisra1a()
810 fvec[i] = (b[0]+b[1]*x+b[2]*xx+b[3]*xxx) / (1.+b[4]*x+b[5]*xx+b[6]*xxx) - m_y[i]; in operator ()()
827 fact = - (b[0]+b[1]*x+b[2]*xx+b[3]*xxx) * fact * fact; in df()
850 x<< 10., -1., .05, -.00001, -.05, .001, -.000001; in testNistHahn1()
853 LevenbergMarquardt<hahn1_functor> lm(functor); in testNistHahn1() local
854 info = lm.minimize(x); in testNistHahn1()
858 LM_CHECK_N_ITERS(lm, 11, 10); in testNistHahn1()
860 VERIFY_IS_APPROX(lm.fvec.squaredNorm(), 1.5324382854E+00); in testNistHahn1()
863 VERIFY_IS_APPROX(x[1],-1.2269296921E-01); in testNistHahn1()
864 VERIFY_IS_APPROX(x[2], 4.0863750610E-03); in testNistHahn1()
865 VERIFY_IS_APPROX(x[3],-1.426264e-06); // shoulde be : -1.4262662514E-06 in testNistHahn1()
866 VERIFY_IS_APPROX(x[4],-5.7609940901E-03); in testNistHahn1()
867 VERIFY_IS_APPROX(x[5], 2.4053735503E-04); in testNistHahn1()
868 VERIFY_IS_APPROX(x[6],-1.2314450199E-07); in testNistHahn1()
873 x<< .1, -.1, .005, -.000001, -.005, .0001, -.0000001; in testNistHahn1()
875 info = lm.minimize(x); in testNistHahn1()
879 LM_CHECK_N_ITERS(lm, 11, 10); in testNistHahn1()
881 VERIFY_IS_APPROX(lm.fvec.squaredNorm(), 1.5324382854E+00); in testNistHahn1()
884 VERIFY_IS_APPROX(x[1], -0.1226933); // should be : -1.2269296921E-01 in testNistHahn1()
885 VERIFY_IS_APPROX(x[2], 0.004086383); // should be : 4.0863750610E-03 in testNistHahn1()
886 VERIFY_IS_APPROX(x[3], -1.426277e-06); // shoulde be : -1.4262662514E-06 in testNistHahn1()
887 VERIFY_IS_APPROX(x[4],-5.7609940901E-03); in testNistHahn1()
888 VERIFY_IS_APPROX(x[5], 0.00024053772); // should be : 2.4053735503E-04 in testNistHahn1()
889 VERIFY_IS_APPROX(x[6], -1.231450e-07); // should be : -1.2314450199E-07 in testNistHahn1()
903 fvec[i] = b[0]*b[1]*x[i]/(1.+b[1]*x[i]) - y[i]; in operator ()()
915 fjac(i,1) = b[0]*x[i]*(den-b[1]*x[i])/den/den; in df()
937 LevenbergMarquardt<misra1d_functor> lm(functor); in testNistMisra1d() local
938 info = lm.minimize(x); in testNistMisra1d()
942 LM_CHECK_N_ITERS(lm, 9, 7); in testNistMisra1d()
944 VERIFY_IS_APPROX(lm.fvec.squaredNorm(), 5.6419295283E-02); in testNistMisra1d()
947 VERIFY_IS_APPROX(x[1], 3.0227324449E-04); in testNistMisra1d()
954 info = lm.minimize(x); in testNistMisra1d()
958 LM_CHECK_N_ITERS(lm, 4, 3); in testNistMisra1d()
960 VERIFY_IS_APPROX(lm.fvec.squaredNorm(), 5.6419295283E-02); in testNistMisra1d()
963 VERIFY_IS_APPROX(x[1], 3.0227324449E-04); in testNistMisra1d()
977 fvec[i] = b[0]*exp(-b[1]*x[i]) + b[2]*exp(-b[3]*x[i]) + b[4]*exp(-b[5]*x[i]) - y[i]; in operator ()()
986 fjac(i,0) = exp(-b[1]*x[i]); in df()
987 fjac(i,1) = -b[0]*x[i]*exp(-b[1]*x[i]); in df()
988 fjac(i,2) = exp(-b[3]*x[i]); in df()
989 fjac(i,3) = -b[2]*x[i]*exp(-b[3]*x[i]); in df()
990 fjac(i,4) = exp(-b[5]*x[i]); in df()
991 fjac(i,5) = -b[4]*x[i]*exp(-b[5]*x[i]); in df()
996 …-02, 1.000000000000E-01, 1.500000000000E-01, 2.000000000000E-01, 2.500000000000E-01, 3.00000000000…
997 …-01 ,7.679338563728E-01 ,6.388775523106E-01 ,5.337835317402E-01 ,4.479363617347E-01 ,3.77584788435…
1013 LevenbergMarquardt<lanczos1_functor> lm(functor); in testNistLanczos1() local
1014 info = lm.minimize(x); in testNistLanczos1()
1018 LM_CHECK_N_ITERS(lm, 79, 72); in testNistLanczos1()
1021 std::cout << lm.fvec.squaredNorm() << "\n"; in testNistLanczos1()
1022 VERIFY(lm.fvec.squaredNorm() <= 1.4307867721E-25); in testNistLanczos1()
1024 VERIFY_IS_APPROX(x[0], 9.5100000027E-02); in testNistLanczos1()
1026 VERIFY_IS_APPROX(x[2], 8.6070000013E-01); in testNistLanczos1()
1036 info = lm.minimize(x); in testNistLanczos1()
1040 LM_CHECK_N_ITERS(lm, 9, 8); in testNistLanczos1()
1042 VERIFY(lm.fvec.squaredNorm() <= 1.4307867721E-25); in testNistLanczos1()
1044 VERIFY_IS_APPROX(x[0], 9.5100000027E-02); in testNistLanczos1()
1046 VERIFY_IS_APPROX(x[2], 8.6070000013E-01); in testNistLanczos1()
1063 fvec[i] = b[0] / (1.+exp(b[1]-b[2]*x[i])) - y[i]; in operator ()()
1074 double e = exp(b[1]-b[2]*x[i]); in df()
1076 fjac(i,1) = -b[0]*e/(1.+e)/(1.+e); in df()
1099 LevenbergMarquardt<rat42_functor> lm(functor); in testNistRat42() local
1100 info = lm.minimize(x); in testNistRat42()
1104 LM_CHECK_N_ITERS(lm, 10, 8); in testNistRat42()
1106 VERIFY_IS_APPROX(lm.fvec.squaredNorm(), 8.0565229338E+00); in testNistRat42()
1110 VERIFY_IS_APPROX(x[2], 6.7359200066E-02); in testNistRat42()
1117 info = lm.minimize(x); in testNistRat42()
1121 LM_CHECK_N_ITERS(lm, 6, 5); in testNistRat42()
1123 VERIFY_IS_APPROX(lm.fvec.squaredNorm(), 8.0565229338E+00); in testNistRat42()
1127 VERIFY_IS_APPROX(x[2], 6.7359200066E-02); in testNistRat42()
1140 fvec[i] = b[0] * exp(b[1]/(x[i]+b[2])) - y[i]; in operator ()()
1153 fjac(i,2) = -b[1]*b[0]*factor*factor*e; in df()
1175 LevenbergMarquardt<MGH10_functor> lm(functor); in testNistMGH10() local
1176 info = lm.minimize(x); in testNistMGH10()
1180 LM_CHECK_N_ITERS(lm, 284, 249); in testNistMGH10()
1182 VERIFY_IS_APPROX(lm.fvec.squaredNorm(), 8.7945855171E+01); in testNistMGH10()
1184 VERIFY_IS_APPROX(x[0], 5.6096364710E-03); in testNistMGH10()
1193 info = lm.minimize(x); in testNistMGH10()
1197 LM_CHECK_N_ITERS(lm, 126, 116); in testNistMGH10()
1199 VERIFY_IS_APPROX(lm.fvec.squaredNorm(), 8.7945855171E+01); in testNistMGH10()
1201 VERIFY_IS_APPROX(x[0], 5.6096364710E-03); in testNistMGH10()
1217 fvec[i] = b[0]*(1.-exp(-b[1]*x[i])) - y[i]; in operator ()()
1226 double e = exp(-b[1]*x[i]); in df()
1227 fjac(i,0) = 1.-e; in df()
1249 LevenbergMarquardt<BoxBOD_functor> lm(functor); in testNistBoxBOD() local
1250 lm.parameters.ftol = 1.E6*NumTraits<double>::epsilon(); in testNistBoxBOD()
1251 lm.parameters.xtol = 1.E6*NumTraits<double>::epsilon(); in testNistBoxBOD()
1252 lm.parameters.factor = 10.; in testNistBoxBOD()
1253 info = lm.minimize(x); in testNistBoxBOD()
1257 LM_CHECK_N_ITERS(lm, 31, 25); in testNistBoxBOD()
1259 VERIFY_IS_APPROX(lm.fvec.squaredNorm(), 1.1680088766E+03); in testNistBoxBOD()
1262 VERIFY_IS_APPROX(x[1], 5.4723748542E-01); in testNistBoxBOD()
1269 lm.resetParameters(); in testNistBoxBOD()
1270 lm.parameters.ftol = NumTraits<double>::epsilon(); in testNistBoxBOD()
1271 lm.parameters.xtol = NumTraits<double>::epsilon(); in testNistBoxBOD()
1272 info = lm.minimize(x); in testNistBoxBOD()
1276 LM_CHECK_N_ITERS(lm, 15, 14); in testNistBoxBOD()
1278 VERIFY_IS_APPROX(lm.fvec.squaredNorm(), 1.1680088766E+03); in testNistBoxBOD()
1281 VERIFY_IS_APPROX(x[1], 5.4723748542E-01); in testNistBoxBOD()
1294 fvec[i] = b[0] + b[1]*exp(-b[3]*x[i]) + b[2]*exp(-b[4]*x[i]) - y[i]; in operator ()()
1304 fjac(i,1) = exp(-b[3]*x[i]); in df()
1305 fjac(i,2) = exp(-b[4]*x[i]); in df()
1306 fjac(i,3) = -x[i]*b[1]*exp(-b[3]*x[i]); in df()
1307 fjac(i,4) = -x[i]*b[2]*exp(-b[4]*x[i]); in df()
1313 …-01, 9.080000E-01, 9.320000E-01, 9.360000E-01, 9.250000E-01, 9.080000E-01, 8.810000E-01, 8.500000E…
1326 x<< 50., 150., -100., 1., 2.; in testNistMGH17()
1329 LevenbergMarquardt<MGH17_functor> lm(functor); in testNistMGH17() local
1330 lm.parameters.ftol = NumTraits<double>::epsilon(); in testNistMGH17()
1331 lm.parameters.xtol = NumTraits<double>::epsilon(); in testNistMGH17()
1332 lm.parameters.maxfev = 1000; in testNistMGH17()
1333 info = lm.minimize(x); in testNistMGH17()
1336 VERIFY_IS_APPROX(lm.fvec.squaredNorm(), 5.4648946975E-05); in testNistMGH17()
1338 VERIFY_IS_APPROX(x[0], 3.7541005211E-01); in testNistMGH17()
1340 VERIFY_IS_APPROX(x[2], -1.4646871366E+00); in testNistMGH17()
1341 VERIFY_IS_APPROX(x[3], 1.2867534640E-02); in testNistMGH17()
1342 VERIFY_IS_APPROX(x[4], 2.2122699662E-02); in testNistMGH17()
1346 LM_CHECK_N_ITERS(lm, 602, 545); in testNistMGH17()
1351 x<< 0.5 ,1.5 ,-1 ,0.01 ,0.02; in testNistMGH17()
1353 lm.resetParameters(); in testNistMGH17()
1354 info = lm.minimize(x); in testNistMGH17()
1358 LM_CHECK_N_ITERS(lm, 18, 15); in testNistMGH17()
1360 VERIFY_IS_APPROX(lm.fvec.squaredNorm(), 5.4648946975E-05); in testNistMGH17()
1362 VERIFY_IS_APPROX(x[0], 3.7541005211E-01); in testNistMGH17()
1364 VERIFY_IS_APPROX(x[2], -1.4646871366E+00); in testNistMGH17()
1365 VERIFY_IS_APPROX(x[3], 1.2867534640E-02); in testNistMGH17()
1366 VERIFY_IS_APPROX(x[4], 2.2122699662E-02); in testNistMGH17()
1380 fvec[i] = b[0]*(xx+x*b[1])/(xx+x*b[2]+b[3]) - y[i]; in operator ()()
1394 fjac(i,2) = - b[0]*(xx+x*b[1]) * x * factor * factor; in df()
1395 fjac(i,3) = - b[0]*(xx+x*b[1]) * factor * factor; in df()
1400 …ctor::_x[11] = { 4., 2., 1., 5.E-1 , 2.5E-01, 1.670000E-01, 1.250000E-01, 1.E-01, 8.330000E-02, 7…
1401 …57000E-01, 1.947000E-01, 1.735000E-01, 1.600000E-01, 8.440000E-02, 6.270000E-02, 4.560000E-02, 3.4…
1417 LevenbergMarquardt<MGH09_functor> lm(functor); in testNistMGH09() local
1418 lm.parameters.maxfev = 1000; in testNistMGH09()
1419 info = lm.minimize(x); in testNistMGH09()
1423 LM_CHECK_N_ITERS(lm, 490, 376); in testNistMGH09()
1425 VERIFY_IS_APPROX(lm.fvec.squaredNorm(), 3.0750560385E-04); in testNistMGH09()
1427 VERIFY_IS_APPROX(x[0], 0.1928077089); // should be 1.9280693458E-01 in testNistMGH09()
1428 VERIFY_IS_APPROX(x[1], 0.19126423573); // should be 1.9128232873E-01 in testNistMGH09()
1429 VERIFY_IS_APPROX(x[2], 0.12305309914); // should be 1.2305650693E-01 in testNistMGH09()
1430 VERIFY_IS_APPROX(x[3], 0.13605395375); // should be 1.3606233068E-01 in testNistMGH09()
1437 lm.resetParameters(); in testNistMGH09()
1438 info = lm.minimize(x); in testNistMGH09()
1442 LM_CHECK_N_ITERS(lm, 18, 16); in testNistMGH09()
1444 VERIFY_IS_APPROX(lm.fvec.squaredNorm(), 3.0750560385E-04); in testNistMGH09()
1446 VERIFY_IS_APPROX(x[0], 0.19280781); // should be 1.9280693458E-01 in testNistMGH09()
1447 VERIFY_IS_APPROX(x[1], 0.19126265); // should be 1.9128232873E-01 in testNistMGH09()
1448 VERIFY_IS_APPROX(x[2], 0.12305280); // should be 1.2305650693E-01 in testNistMGH09()
1449 VERIFY_IS_APPROX(x[3], 0.13605322); // should be 1.3606233068E-01 in testNistMGH09()
1464 fvec[i] = b[0]* pow(b[1]+x[i],-1./b[2]) - y[i]; in operator ()()
1473 double e = pow(b[1]+x[i],-1./b[2]); in df()
1475 fjac(i,1) = - b[0]*e/b[2]/(b[1]+x[i]); in df()
1483 …-34.834702E0 ,-34.393200E0 ,-34.152901E0 ,-33.979099E0 ,-33.845901E0 ,-33.732899E0 ,-33.640301E0 ,…
1484 …-32.238602E0 ,-32.229900E0 ,-32.224098E0 ,-32.215401E0 ,-32.203800E0 ,-32.198002E0 ,-32.189400E0 ,…
1485 -31.827299E0 ,-31.821600E0 ,-31.818701E0 ,-31.812901E0 ,-31.809999E0 ,-31.807100E0 ,-31.801300E0 ,-…
1498 x<< -2000., 50., 0.8; in testNistBennett5()
1501 LevenbergMarquardt<Bennett5_functor> lm(functor); in testNistBennett5() local
1502 lm.parameters.maxfev = 1000; in testNistBennett5()
1503 info = lm.minimize(x); in testNistBennett5()
1507 LM_CHECK_N_ITERS(lm, 758, 744); in testNistBennett5()
1509 VERIFY_IS_APPROX(lm.fvec.squaredNorm(), 5.2404744073E-04); in testNistBennett5()
1511 VERIFY_IS_APPROX(x[0], -2.5235058043E+03); in testNistBennett5()
1513 VERIFY_IS_APPROX(x[2], 9.3218483193E-01); in testNistBennett5()
1517 x<< -1500., 45., 0.85; in testNistBennett5()
1519 lm.resetParameters(); in testNistBennett5()
1520 info = lm.minimize(x); in testNistBennett5()
1524 LM_CHECK_N_ITERS(lm, 203, 192); in testNistBennett5()
1526 VERIFY_IS_APPROX(lm.fvec.squaredNorm(), 5.2404744073E-04); in testNistBennett5()
1528 VERIFY_IS_APPROX(x[0], -2523.3007865); // should be -2.5235058043E+03 in testNistBennett5()
1530 VERIFY_IS_APPROX(x[2], 0.93219881891); // should be 9.3218483193E-01); in testNistBennett5()
1545 fvec[i] = (b[0]+b[1]*x+b[2]*xx+b[3]*xxx) / (1.+b[4]*x+b[5]*xx+b[6]*xxx) - _y[i]; in operator ()()
1561 fact = - (b[0]+b[1]*x+b[2]*xx+b[3]*xxx) * fact * fact; in df()
1569 …-3.067E0, -2.981E0, -2.921E0, -2.912E0, -2.840E0, -2.797E0, -2.702E0, -2.699E0, -2.633E0, -2.481E0…
1586 LevenbergMarquardt<thurber_functor> lm(functor); in testNistThurber() local
1587 lm.parameters.ftol = 1.E4*NumTraits<double>::epsilon(); in testNistThurber()
1588 lm.parameters.xtol = 1.E4*NumTraits<double>::epsilon(); in testNistThurber()
1589 info = lm.minimize(x); in testNistThurber()
1593 LM_CHECK_N_ITERS(lm, 39,36); in testNistThurber()
1595 VERIFY_IS_APPROX(lm.fvec.squaredNorm(), 5.6427082397E+03); in testNistThurber()
1601 VERIFY_IS_APPROX(x[4], 9.6629502864E-01); in testNistThurber()
1602 VERIFY_IS_APPROX(x[5], 3.9797285797E-01); in testNistThurber()
1603 VERIFY_IS_APPROX(x[6], 4.9727297349E-02); in testNistThurber()
1610 lm.resetParameters(); in testNistThurber()
1611 lm.parameters.ftol = 1.E4*NumTraits<double>::epsilon(); in testNistThurber()
1612 lm.parameters.xtol = 1.E4*NumTraits<double>::epsilon(); in testNistThurber()
1613 info = lm.minimize(x); in testNistThurber()
1617 LM_CHECK_N_ITERS(lm, 29, 28); in testNistThurber()
1619 VERIFY_IS_APPROX(lm.fvec.squaredNorm(), 5.6427082397E+03); in testNistThurber()
1625 VERIFY_IS_APPROX(x[4], 9.6629502864E-01); in testNistThurber()
1626 VERIFY_IS_APPROX(x[5], 3.9797285797E-01); in testNistThurber()
1627 VERIFY_IS_APPROX(x[6], 4.9727297349E-02); in testNistThurber()
1640 fvec[i] = b[0] * pow(1.+exp(b[1]-b[2]*x[i]),-1./b[3]) - y[i]; in operator ()()
1649 double e = exp(b[1]-b[2]*x[i]); in df()
1650 double power = -1./b[3]; in df()
1652 fjac(i,1) = power*b[0]*e*pow(1.+e, power-1.); in df()
1653 fjac(i,2) = -power*b[0]*e*x[i]*pow(1.+e, power-1.); in df()
1676 LevenbergMarquardt<rat43_functor> lm(functor); in testNistRat43() local
1677 lm.parameters.ftol = 1.E6*NumTraits<double>::epsilon(); in testNistRat43()
1678 lm.parameters.xtol = 1.E6*NumTraits<double>::epsilon(); in testNistRat43()
1679 info = lm.minimize(x); in testNistRat43()
1683 LM_CHECK_N_ITERS(lm, 27, 20); in testNistRat43()
1685 VERIFY_IS_APPROX(lm.fvec.squaredNorm(), 8.7864049080E+03); in testNistRat43()
1689 VERIFY_IS_APPROX(x[2], 7.5962938329E-01); in testNistRat43()
1697 lm.resetParameters(); in testNistRat43()
1698 lm.parameters.ftol = 1.E5*NumTraits<double>::epsilon(); in testNistRat43()
1699 lm.parameters.xtol = 1.E5*NumTraits<double>::epsilon(); in testNistRat43()
1700 info = lm.minimize(x); in testNistRat43()
1704 LM_CHECK_N_ITERS(lm, 9, 8); in testNistRat43()
1706 VERIFY_IS_APPROX(lm.fvec.squaredNorm(), 8.7864049080E+03); in testNistRat43()
1710 VERIFY_IS_APPROX(x[2], 7.5962938329E-01); in testNistRat43()
1726 fvec[i] = b[0]/b[1] * exp(-0.5*(x[i]-b[2])*(x[i]-b[2])/(b[1]*b[1])) - y[i]; in operator ()()
1736 double e = exp(-0.5*(x[i]-b[2])*(x[i]-b[2])/b12); in df()
1738 fjac(i,1) = ((x[i]-b[2])*(x[i]-b[2])/b12-1.) * b[0]*e/b12; in df()
1739 fjac(i,2) = (x[i]-b[2])*e*b[0]/b[1]/b12; in df()
1761 LevenbergMarquardt<eckerle4_functor> lm(functor); in testNistEckerle4() local
1762 info = lm.minimize(x); in testNistEckerle4()
1766 LM_CHECK_N_ITERS(lm, 18, 15); in testNistEckerle4()
1768 VERIFY_IS_APPROX(lm.fvec.squaredNorm(), 1.4635887487E-03); in testNistEckerle4()
1779 info = lm.minimize(x); in testNistEckerle4()
1783 LM_CHECK_N_ITERS(lm, 7, 6); in testNistEckerle4()
1785 VERIFY_IS_APPROX(lm.fvec.squaredNorm(), 1.4635887487E-03); in testNistEckerle4()
1830 printf("info, nfev : %d, %d\n", info, lm.nfev);
1838 printf("fvec.blueNorm() : %.32g\n", lm.fvec.blueNorm());
1840 printf("info, nfev, njev : %d, %d, %d\n", info, lm.nfev, lm.njev);
1841 printf("fvec.squaredNorm() : %.13g\n", lm.fvec.squaredNorm());