3 #ifndef DUNE_PDELAB_BACKEND_EIGEN_SOLVERS_HH 4 #define DUNE_PDELAB_BACKEND_EIGEN_SOLVERS_HH 6 #include <dune/common/deprecated.hh> 7 #include <dune/common/timer.hh> 8 #include <dune/common/parallel/mpihelper.hh> 12 #include "../solver.hh" 16 #include <Eigen/Eigen> 17 #include <Eigen/Sparse> 27 template<
class PreconditionerImp>
28 class EigenBackend_BiCGSTAB_Base
29 :
public LinearResultStorage
37 explicit EigenBackend_BiCGSTAB_Base(
unsigned maxiter_=5000)
48 template<
class M,
class V,
class W>
49 void apply(M& A, V& z, W& r,
typename W::field_type reduction)
54 using Mat = Native<M>;
55 ::Eigen::BiCGSTAB<Mat, PreconditionerImp> solver;
56 solver.setMaxIterations(maxiter);
57 solver.setTolerance(reduction);
62 double elapsed = watch.elapsed();
64 res.converged = solver.info() == ::Eigen::ComputationInfo::Success;
65 res.iterations = solver.iterations();
66 res.elapsed = elapsed;
67 res.reduction = solver.error();
73 typename Dune::template FieldTraits<typename V::ElementType >::real_type norm(
const V& v)
const 83 class EigenBackend_BiCGSTAB_IILU
84 :
public EigenBackend_BiCGSTAB_Base<::Eigen::IncompleteLUT<double> >
87 explicit EigenBackend_BiCGSTAB_IILU(
unsigned maxiter_=5000)
88 : EigenBackend_BiCGSTAB_Base(maxiter_)
93 class EigenBackend_BiCGSTAB_Diagonal
94 :
public EigenBackend_BiCGSTAB_Base<::Eigen::DiagonalPreconditioner<double> >
97 explicit EigenBackend_BiCGSTAB_Diagonal(
unsigned maxiter_=5000)
98 : EigenBackend_BiCGSTAB_Base(maxiter_)
102 template<
class Preconditioner,
int UpLo >
103 class EigenBackend_CG_Base
104 :
public SequentialNorm,
public LinearResultStorage
112 explicit EigenBackend_CG_Base(
unsigned maxiter_=5000)
123 template<
class M,
class V,
class W>
124 void apply(M& A, V& z, W& r,
typename W::field_type reduction)
127 using Mat = Backend::Native<M>;
128 ::Eigen::ConjugateGradient<Mat, UpLo, Preconditioner> solver;
129 solver.setMaxIterations(maxiter);
130 solver.setTolerance(reduction);
133 solver.compute(
native(A));
135 double elapsed = watch.elapsed();
138 res.converged = solver.info() == ::Eigen::ComputationInfo::Success;
139 res.iterations = solver.iterations();
140 res.elapsed = elapsed;
141 res.reduction = solver.error();
151 class EigenBackend_CG_IILU_Up
152 :
public EigenBackend_CG_Base<::Eigen::IncompleteLUT<double>, ::Eigen::Upper >
155 explicit EigenBackend_CG_IILU_Up(
unsigned maxiter_=5000)
156 : EigenBackend_CG_Base(maxiter_)
160 class EigenBackend_CG_Diagonal_Up
161 :
public EigenBackend_CG_Base<::Eigen::DiagonalPreconditioner<double>, ::Eigen::Upper >
164 explicit EigenBackend_CG_Diagonal_Up(
unsigned maxiter_=5000)
165 : EigenBackend_CG_Base(maxiter_)
169 class EigenBackend_CG_IILU_Lo
170 :
public EigenBackend_CG_Base<::Eigen::IncompleteLUT<double>, ::Eigen::Lower >
173 explicit EigenBackend_CG_IILU_Lo(
unsigned maxiter_=5000)
174 : EigenBackend_CG_Base(maxiter_)
178 class EigenBackend_CG_Diagonal_Lo
179 :
public EigenBackend_CG_Base<::Eigen::DiagonalPreconditioner<double>, ::Eigen::Lower >
182 explicit EigenBackend_CG_Diagonal_Lo(
unsigned maxiter_=5000)
183 : EigenBackend_CG_Base(maxiter_)
187 #if EIGEN_VERSION_AT_LEAST(3,2,2) 188 template<
template<
class,
int,
class>
class Solver,
int UpLo>
190 template<
template<
class,
int>
class Solver,
int UpLo>
192 class EigenBackend_SPD_Base
193 :
public SequentialNorm,
public LinearResultStorage
201 explicit EigenBackend_SPD_Base()
211 template<
class M,
class V,
class W>
212 void apply(M& A, V& z, W& r,
typename W::field_type reduction)
215 using Mat = Backend::Native<M>;
216 #if EIGEN_VERSION_AT_LEAST(3,2,2) 220 Solver<Mat, UpLo, ::Eigen::AMDOrdering<typename Mat::Index> > solver;
222 Solver<Mat, UpLo> solver;
226 solver.compute(
native(A));
228 double elapsed = watch.elapsed();
230 res.converged = solver.info() == ::Eigen::ComputationInfo::Success;
231 res.iterations = solver.iterations();
232 res.elapsed = elapsed;
233 res.reduction = solver.error();
239 class EigenBackend_SimplicialCholesky_Up
240 :
public EigenBackend_SPD_Base<::Eigen::SimplicialCholesky, ::Eigen::Upper >
243 explicit EigenBackend_SimplicialCholesky_Up()
247 class EigenBackend_SimplicialCholesky_Lo
248 :
public EigenBackend_SPD_Base<::Eigen::SimplicialCholesky, ::Eigen::Lower >
251 explicit EigenBackend_SimplicialCholesky_Lo()
303 class EigenBackend_LeastSquares
304 :
public SequentialNorm,
public LinearResultStorage
307 const static unsigned int defaultFlag = ::Eigen::ComputeThinU | ::Eigen::ComputeThinV;
314 explicit EigenBackend_LeastSquares(
unsigned int flags = defaultFlag)
325 template<
class M,
class V,
class W>
326 void apply(M& A, V& z, W& r,
typename W::field_type reduction)
331 using Mat = Backend::Native<M>;
332 ::Eigen::JacobiSVD<Mat, ::Eigen::ColPivHouseholderQRPreconditioner> solver(A, flags_);
334 double elapsed = watch.elapsed();
336 res.converged = solver.info() == ::Eigen::ComputationInfo::Success;
337 res.iterations = solver.iterations();
338 res.elapsed = elapsed;
339 res.reduction = solver.error();
351 #endif // DUNE_PDELAB_BACKEND_EIGEN_SOLVERS_HH std::enable_if< std::is_base_of< impl::WrapperBase, T >::value, Native< T > &>::type native(T &t)
Definition: backend/interface.hh:192
typename native_type< T >::type Native
Alias of the native container type associated with T or T itself if it is not a backend wrapper...
Definition: backend/interface.hh:176
For backward compatibility – Do not use this!
Definition: adaptivity.hh:27