another nist test with difficulty 'leverage', this one passes

1 files changed, 73 insertions, 0 deletions

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<misra1d_functor, double>(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<misra1d_functor, double>(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());