Generalize IterativeSolverBase::solve to hanlde any sparse type as the results (instead of SparseMatrix only)
diff --git a/Eigen/src/IterativeLinearSolvers/IterativeSolverBase.h b/Eigen/src/IterativeLinearSolvers/IterativeSolverBase.h index 3d62fef..7c2326e 100644 --- a/Eigen/src/IterativeLinearSolvers/IterativeSolverBase.h +++ b/Eigen/src/IterativeLinearSolvers/IterativeSolverBase.h
@@ -330,25 +330,27 @@ } /** \internal */ - template<typename Rhs, typename DestScalar, int DestOptions, typename DestIndex> - void _solve_impl(const Rhs& b, SparseMatrix<DestScalar,DestOptions,DestIndex> &dest) const + template<typename Rhs, typename DestDerived> + void _solve_impl(const Rhs& b, SparseMatrixBase<DestDerived> &aDest) const { eigen_assert(rows()==b.rows()); Index rhsCols = b.cols(); Index size = b.rows(); + DestDerived& dest(aDest.derived()); + typedef typename DestDerived::Scalar DestScalar; Eigen::Matrix<DestScalar,Dynamic,1> tb(size); Eigen::Matrix<DestScalar,Dynamic,1> tx(cols()); // We do not directly fill dest because sparse expressions have to be free of aliasing issue. // For non square least-square problems, b and dest might not have the same size whereas they might alias each-other. - SparseMatrix<DestScalar,DestOptions,DestIndex> tmp(cols(),rhsCols); + typename DestDerived::PlainObject tmp(cols(),rhsCols); for(Index k=0; k<rhsCols; ++k) { tb = b.col(k); tx = derived().solve(tb); tmp.col(k) = tx.sparseView(0); } - tmp.swap(dest); + dest.swap(tmp); } protected: