diff --git a/unsupported/Eigen/src/NonLinear/LevenbergMarquardt.h b/unsupported/Eigen/src/NonLinear/LevenbergMarquardt.h
index 4e0565237..f541c385b 100644
--- a/unsupported/Eigen/src/NonLinear/LevenbergMarquardt.h
+++ b/unsupported/Eigen/src/NonLinear/LevenbergMarquardt.h
@@ -55,7 +55,7 @@ public:
             );
 
     static Status lmdif1(
-            FunctorType &_functor,
+            FunctorType &functor,
             Matrix< Scalar, Dynamic, 1 >  &x,
             int *nfev,
             const Scalar tol = ei_sqrt(epsilon<Scalar>())
@@ -200,9 +200,13 @@ LevenbergMarquardt<FunctorType,Scalar>::minimizeOneStep(
 
     /* calculate the jacobian matrix. */
 
-    if (functor.df(x, fjac) < 0)
+    int df_ret = functor.df(x, fjac);
+    if (df_ret<0)
         return UserAsked;
-    ++njev;
+    if (df_ret>0)
+        // numerical diff, we evaluated the function df_ret times
+        nfev += df_ret;
+    else njev++;
 
     /* compute the qr factorization of the jacobian. */
 
@@ -702,10 +706,8 @@ LevenbergMarquardt<FunctorType,Scalar>::lmdif1(
     lm.parameters.maxfev = 200*(n+1);
 
     Status info = Status(lm.minimize(x));
-
     if (nfev)
         * nfev = lm.nfev;
-
     return info;
 }
 
diff --git a/unsupported/Eigen/src/NumericalDiff/NumericalDiff.h b/unsupported/Eigen/src/NumericalDiff/NumericalDiff.h
index 276b315f8..dcd435864 100644
--- a/unsupported/Eigen/src/NumericalDiff/NumericalDiff.h
+++ b/unsupported/Eigen/src/NumericalDiff/NumericalDiff.h
@@ -74,6 +74,7 @@ public:
         val1.resize(Functor::values());
         val2.resize(Functor::values());
 
+        // initialization
         switch(mode) {
             case Forward:
                 // compute f(x)
@@ -86,8 +87,7 @@ public:
                 assert(false);
         };
 
-        /* Function Body */
-
+        // Function Body
         for (int j = 0; j < n; ++j) {
             h = eps * ei_abs(x[j]);
             if (h == 0.) {
diff --git a/unsupported/test/NonLinear.cpp b/unsupported/test/NonLinear.cpp
index ee081cfe7..6403bd457 100644
--- a/unsupported/test/NonLinear.cpp
+++ b/unsupported/test/NonLinear.cpp
@@ -561,7 +561,7 @@ void testLmdif1()
 
   // check return value
   VERIFY( 1 == info);
-  VERIFY(nfev==21);
+  VERIFY(nfev==26);
 
   // check norm
   functor(x, fvec);
@@ -592,7 +592,7 @@ void testLmdif()
 
   // check return values
   VERIFY( 1 == info);
-  VERIFY(lm.nfev==21);
+  VERIFY(lm.nfev==26);
 
   // check norm
   fnorm = lm.fvec.blueNorm();