From 7ecbbc9aa4d8f0fa3a99504898f21de77b9a8df4 Mon Sep 17 00:00:00 2001 From: Thomas Capricelli Date: Mon, 10 Aug 2009 12:34:51 +0200 Subject: another nist test with difficulty 'leverage', this one passes --- unsupported/test/NonLinear.cpp | 73 ++++++++++++++++++++++++++++++++++++++++++ 1 file changed, 73 insertions(+) diff --git a/unsupported/test/NonLinear.cpp b/unsupported/test/NonLinear.cpp index 72885cac7..d541b248c 100644 --- a/unsupported/test/NonLinear.cpp +++ b/unsupported/test/NonLinear.cpp @@ -1017,10 +1017,83 @@ void testNistHahn1(void) } +struct misra1d_functor { + static int f(void * /*p*/, int m, int n, const double *b, double *fvec, double *fjac, int ldfjac, int iflag) + { + static const double x[14] = { 77.6E0, 114.9E0, 141.1E0, 190.8E0, 239.9E0, 289.0E0, 332.8E0, 378.4E0, 434.8E0, 477.3E0, 536.8E0, 593.1E0, 689.1E0, 760.0E0}; + static const double y[14] = { 10.07E0, 14.73E0, 17.94E0, 23.93E0, 29.61E0, 35.18E0, 40.02E0, 44.82E0, 50.76E0, 55.05E0, 61.01E0, 66.40E0, 75.47E0, 81.78E0}; + int i; + + assert(m==14); + assert(n==2); + assert(ldfjac==14); + if (iflag != 2) {// compute fvec at b + for(i=0; i<14; i++) { + fvec[i] = b[0]*b[1]*x[i]/(1.+b[1]*x[i]) - y[i]; + } + } + else { // compute fjac at b + for(i=0; i<14; i++) { + double den = 1.+b[1]*x[i]; + fjac[i+ldfjac*0] = b[1]*x[i] / den; + fjac[i+ldfjac*1] = b[0]*x[i]*(den-b[1]*x[i])/den/den; + } + } + return 0; + } +}; + +// http://www.itl.nist.gov/div898/strd/nls/data/misra1d.shtml +void testNistMisra1d(void) +{ + const int m=14, n=2; + int info, nfev, njev; + + Eigen::VectorXd x(n), fvec(m), wa1, diag; + Eigen::MatrixXd fjac; + VectorXi ipvt; + + /* + * First try + */ + x<< 500., 0.0001; + // do the computation + info = ei_lmder(x, fvec, nfev, njev, fjac, ipvt, wa1, diag); + + // check return value + VERIFY( 3 == info); + VERIFY( 9 == nfev); + VERIFY( 7 == njev); + // check norm^2 + VERIFY_IS_APPROX(fvec.squaredNorm(), 5.6419295283E-02); + // check x + VERIFY_IS_APPROX(x[0], 4.3736970754E+02); + VERIFY_IS_APPROX(x[1], 3.0227324449E-04); + + /* + * Second try + */ + x<< 450., 0.0003; + // do the computation + info = ei_lmder(x, fvec, nfev, njev, fjac, ipvt, wa1, diag); + + // check return value + VERIFY( 1 == info); + VERIFY( 4 == nfev); + VERIFY( 3 == njev); + // check norm^2 + VERIFY_IS_APPROX(fvec.squaredNorm(), 5.6419295283E-02); + // check x + VERIFY_IS_APPROX(x[0], 4.3736970754E+02); + VERIFY_IS_APPROX(x[1], 3.0227324449E-04); +} + + void test_NonLinear() { CALL_SUBTEST(testNistMisra1a()); // CALL_SUBTEST(testNistHahn1()); + CALL_SUBTEST(testNistMisra1d()); CALL_SUBTEST(testChkder()); CALL_SUBTEST(testLmder1()); -- cgit v1.2.3