get rid of ei_qform + lot of other cleaning, now that we do not depend on

minpack qr factorization.

2 files changed, 10 insertions, 101 deletions

diff --git a/unsupported/Eigen/src/NonLinearOptimization/HybridNonLinearSolver.h b/unsupported/Eigen/src/NonLinearOptimization/HybridNonLinearSolver.h
index ee5b63c44..4b985db74 100644
--- a/unsupported/Eigen/src/NonLinearOptimization/HybridNonLinearSolver.h
+++ b/unsupported/Eigen/src/NonLinearOptimization/HybridNonLinearSolver.h
@@ -217,7 +217,7 @@ HybridNonLinearSolver<FunctorType,Scalar>::solveOneStep(
         const int mode
         )
 {
-    int i, j;
+    int j;
     jeval = true;
 
     /* calculate the jacobian matrix. */
@@ -246,30 +246,15 @@ HybridNonLinearSolver<FunctorType,Scalar>::solveOneStep(
     /* compute the qr factorization of the jacobian. */
     wa2 = fjac.colwise().blueNorm();
     HouseholderQR<JacobianType> qrfac(fjac); // no pivoting:
-    fjac = qrfac.matrixQR();
-    wa1 = fjac.diagonal();
-    fjac.diagonal() = qrfac.hCoeffs();
-    // TODO : avoid this:
-    for(int ii=0; ii< fjac.cols(); ii++) fjac.col(ii).segment(ii+1, fjac.rows()-ii-1) *= fjac(ii,ii); // rescale vectors
-
-    /* form (q transpose)*fvec and store in qtf. */
-    qtf = fvec;
-    for (j = 0; j < n; ++j)
-        if (fjac(j,j) != 0.) {
-            sum = 0.;
-            for (i = j; i < n; ++i)
-                sum += fjac(i,j) * qtf[i];
-            temp = -sum / fjac(j,j);
-            for (i = j; i < n; ++i)
-                qtf[i] += fjac(i,j) * temp;
-        }
 
     /* copy the triangular factor of the qr factorization into r. */
     R = qrfac.matrixQR();
-    sing = wa1.cwiseAbs().minCoeff()==0.;
 
     /* accumulate the orthogonal factor in fjac. */
-    ei_qform<Scalar>(n, n, fjac.data(), fjac.rows(), wa1.data());
+    fjac = qrfac.householderQ();
+
+    /* form (q transpose)*fvec and store in qtf. */
+    qtf = fjac.transpose() * fvec;
 
     /* rescale if necessary. */
     if (mode != 2)
@@ -480,13 +465,12 @@ HybridNonLinearSolver<FunctorType,Scalar>::solveNumericalDiffOneStep(
         const int mode
         )
 {
-    int i, j;
+    int j;
     jeval = true;
     if (parameters.nb_of_subdiagonals<0) parameters.nb_of_subdiagonals= n-1;
     if (parameters.nb_of_superdiagonals<0) parameters.nb_of_superdiagonals= n-1;
 
     /* calculate the jacobian matrix. */
-
     if (ei_fdjac1(functor, x, fvec, fjac, parameters.nb_of_subdiagonals, parameters.nb_of_superdiagonals, parameters.epsfcn) <0)
         return UserAksed;
     nfev += std::min(parameters.nb_of_subdiagonals+parameters.nb_of_superdiagonals+ 1, n);
@@ -512,30 +496,15 @@ HybridNonLinearSolver<FunctorType,Scalar>::solveNumericalDiffOneStep(
     /* compute the qr factorization of the jacobian. */
     wa2 = fjac.colwise().blueNorm();
     HouseholderQR<JacobianType> qrfac(fjac); // no pivoting:
-    fjac = qrfac.matrixQR();
-    wa1 = fjac.diagonal();
-    fjac.diagonal() = qrfac.hCoeffs();
-    // TODO : avoid this:
-    for(int ii=0; ii< fjac.cols(); ii++) fjac.col(ii).segment(ii+1, fjac.rows()-ii-1) *= fjac(ii,ii); // rescale vectors
-
-    /* form (q transpose)*fvec and store in qtf. */
-    qtf = fvec;
-    for (j = 0; j < n; ++j)
-        if (fjac(j,j) != 0.) {
-            sum = 0.;
-            for (i = j; i < n; ++i)
-                sum += fjac(i,j) * qtf[i];
-            temp = -sum / fjac(j,j);
-            for (i = j; i < n; ++i)
-                qtf[i] += fjac(i,j) * temp;
-        }
 
     /* copy the triangular factor of the qr factorization into r. */
     R = qrfac.matrixQR();
-    sing = wa1.cwiseAbs().minCoeff()==0.;
 
     /* accumulate the orthogonal factor in fjac. */
-    ei_qform<Scalar>(n, n, fjac.data(), fjac.rows(), wa1.data());
+    fjac = qrfac.householderQ();
+
+    /* form (q transpose)*fvec and store in qtf. */
+    qtf = fjac.transpose() * fvec;
 
     /* rescale if necessary. */
     if (mode != 2)
diff --git a/unsupported/Eigen/src/NonLinearOptimization/qform.h b/unsupported/Eigen/src/NonLinearOptimization/qform.h
deleted file mode 100644
index 1c3e3e267..000000000
--- a/unsupported/Eigen/src/NonLinearOptimization/qform.h
+++ /dev/null
@@ -1,60 +0,0 @@
-
-template <typename Scalar>
-void ei_qform(int m, int n, Scalar *q, int
-        ldq, Scalar *wa)
-{
-    /* System generated locals */
-    int q_dim1, q_offset;
-
-    /* Local variables */
-    int i, j, k, l;
-    Scalar sum, temp;
-    int minmn;
-
-    /* Parameter adjustments */
-    --wa;
-    q_dim1 = ldq;
-    q_offset = 1 + q_dim1 * 1;
-    q -= q_offset;
-
-    /* Function Body */
-
-    /*     zero out upper triangle of q in the first min(m,n) columns. */
-
-    minmn = std::min(m,n);
-    for (j = 2; j <= minmn; ++j) {
-        for (i = 1; i <= j-1; ++i)
-            q[i + j * q_dim1] = 0.;
-    }
-
-    /*     initialize remaining columns to those of the identity matrix. */
-
-    for (j = n+1; j <= m; ++j) {
-        for (i = 1; i <= m; ++i)
-            q[i + j * q_dim1] = 0.;
-        q[j + j * q_dim1] = 1.;
-    }
-
-    /*     accumulate q from its factored form. */
-
-    for (l = 1; l <= minmn; ++l) {
-        k = minmn - l + 1;
-        for (i = k; i <= m; ++i) {
-            wa[i] = q[i + k * q_dim1];
-            q[i + k * q_dim1] = 0.;
-        }
-        q[k + k * q_dim1] = 1.;
-        if (wa[k] == 0.)
-            continue;
-        for (j = k; j <= m; ++j) {
-            sum = 0.;
-            for (i = k; i <= m; ++i)
-                sum += q[i + j * q_dim1] * wa[i];
-            temp = sum / wa[k];
-            for (i = k; i <= m; ++i)
-                q[i + j * q_dim1] -= temp * wa[i];
-        }
-    }
-
-}
-