add analyzePattern/factorize API to iterative solvers and basic preconditioners

Eigen/src/IterativeLinearSolvers/BasicPreconditioners.h

@@ -64,7 +64,13 @@ class DiagonalPreconditioner
     Index cols() const { return m_invdiag.size(); }
     
     template<typename MatrixType>
-    DiagonalPreconditioner& compute(const MatrixType& mat)
+    DiagonalPreconditioner& analyzePattern(const MatrixType& )
+    {
+      return *this;
+    }
+    
+    template<typename MatrixType>
+    DiagonalPreconditioner& factorize(const MatrixType& mat)
     {
       m_invdiag.resize(mat.cols());
       for(int j=0; j<mat.outerSize(); ++j)
@@ -80,6 +86,12 @@ class DiagonalPreconditioner
       return *this;
     }
     
+    template<typename MatrixType>
+    DiagonalPreconditioner& compute(const MatrixType& mat)
+    {
+      return factorize(mat);
+    }
+
     template<typename Rhs, typename Dest>
     void _solve(const Rhs& b, Dest& x) const
     {
@@ -131,6 +143,12 @@ class IdentityPreconditioner
     template<typename MatrixType>
     IdentityPreconditioner(const MatrixType& ) {}
     
+    template<typename MatrixType>
+    IdentityPreconditioner& analyzePattern(const MatrixType& ) { return *this; }
+    
+    template<typename MatrixType>
+    IdentityPreconditioner& factorize(const MatrixType& ) { return *this; }
+
     template<typename MatrixType>
     IdentityPreconditioner& compute(const MatrixType& ) { return *this; }
     

Eigen/src/IterativeLinearSolvers/IncompleteLUT.h

@@ -97,10 +97,10 @@ class IncompleteLUT
     void factorize(const MatrixType& amat);
     
     /**
-    * Compute an incomplete LU factorization with dual threshold on the matrix mat
+      * Compute an incomplete LU factorization with dual threshold on the matrix mat
-    * No pivoting is done in this version
+      * No pivoting is done in this version
-    * 
+      * 
-    **/
+      **/
     template<typename MatrixType>
     IncompleteLUT<Scalar>& compute(const MatrixType& amat)
     {

Eigen/src/IterativeLinearSolvers/IterativeSolverBase.h

@@ -71,6 +71,39 @@ public:
 
   ~IterativeSolverBase() {}
   
+  /** Initializes the iterative solver for the sparcity pattern of the matrix \a A for further solving \c Ax=b problems.
+    *
+    * Currently, this function mostly call analyzePattern on the preconditioner. In the future
+    * we might, for instance, implement column reodering for faster matrix vector products.
+    */
+  Derived& analyzePattern(const MatrixType& A)
+  {
+    m_preconditioner.analyzePattern(A);
+    m_isInitialized = true;
+    m_analysisIsOk = true;
+    m_info = Success;
+    return derived();
+  }
+  
+  /** Initializes the iterative solver with the numerical values of the matrix \a A for further solving \c Ax=b problems.
+    *
+    * Currently, this function mostly call factorize on the preconditioner.
+    *
+    * \warning this class stores a reference to the matrix A as well as some
+    * precomputed values that depend on it. Therefore, if \a A is changed
+    * this class becomes invalid. Call compute() to update it with the new
+    * matrix A, or modify a copy of A.
+    */
+  Derived& factorize(const MatrixType& A)
+  {
+    eigen_assert(m_analysisIsOk && "You must first call analyzePattern()"); 
+    mp_matrix = &A;
+    m_preconditioner.factorize(A);
+    m_factorizationIsOk = true;
+    m_info = Success;
+    return derived();
+  }
+
   /** Initializes the iterative solver with the matrix \a A for further solving \c Ax=b problems.
     *
     * Currently, this function mostly initialized/compute the preconditioner. In the future
@@ -86,6 +119,8 @@ public:
     mp_matrix = &A;
     m_preconditioner.compute(A);
     m_isInitialized = true;
+    m_analysisIsOk = true;
+    m_factorizationIsOk = true;
     m_info = Success;
     return derived();
   }
@@ -194,6 +229,8 @@ protected:
   void init()
   {
     m_isInitialized = false;
+    m_analysisIsOk = false;
+    m_factorizationIsOk = false;
     m_maxIterations = -1;
     m_tolerance = NumTraits<Scalar>::epsilon();
   }
@@ -206,7 +243,7 @@ protected:
   mutable RealScalar m_error;
   mutable int m_iterations;
   mutable ComputationInfo m_info;
-  mutable bool m_isInitialized;
+  mutable bool m_isInitialized, m_analysisIsOk, m_factorizationIsOk;
 };
 
 namespace internal {