6#ifndef DUNE_ISTL_SOLVERS_HH
7#define DUNE_ISTL_SOLVERS_HH
17#include <dune/common/exceptions.hh>
18#include <dune/common/math.hh>
19#include <dune/common/simd/io.hh>
20#include <dune/common/simd/simd.hh>
21#include <dune/common/std/type_traits.hh>
22#include <dune/common/timer.hh>
79 _op->applyscaleadd(-1,x,b);
83 if(iteration.step(0, def)){
99 _op->applyscaleadd(-1,v,b);
101 if(iteration.step(i, def))
147 _op->applyscaleadd(-1,x,b);
150 if(iteration.step(0, def)){
165 auto alpha =
_sp->dot(q,p);
168 _sp->dot(p,b)/alpha);
173 if(iteration.step(i, def))
224 condition_estimate_(condition_estimate)
227 condition_estimate_ =
false;
228 std::cerr <<
"WARNING: Condition estimate was disabled. It is only available for double and float field types!" << std::endl;
240 scalar_real_type reduction,
int maxit,
int verbose,
bool condition_estimate) :
IterativeSolver<X,X>(op, sp, prec, reduction, maxit, verbose),
241 condition_estimate_(condition_estimate)
243 if (condition_estimate && !(std::is_same<field_type,float>::value || std::is_same<field_type,double>::value)) {
244 condition_estimate_ =
false;
245 std::cerr <<
"WARNING: Condition estimate was disabled. It is only available for double and float field types!" << std::endl;
258 scalar_real_type reduction,
int maxit,
int verbose,
bool condition_estimate)
260 condition_estimate_(condition_estimate)
262 if (condition_estimate && !(std::is_same<field_type,float>::value || std::is_same<field_type,double>::value)) {
263 condition_estimate_ =
false;
264 std::cerr <<
"WARNING: Condition estimate was disabled. It is only available for double and float field types!" << std::endl;
284 _op->applyscaleadd(-1,x,b);
287 if(iteration.step(0, def)){
296 std::vector<real_type> lambdas(0);
297 std::vector<real_type> betas(0);
305 rholast =
_sp->dot(p,b);
313 alpha =
_sp->dot(p,q);
316 if (condition_estimate_)
317 lambdas.push_back(std::real(lambda));
323 if(iteration.step(i, def))
332 if (condition_estimate_)
333 betas.push_back(std::real(beta));
341 if (condition_estimate_) {
352 for (
auto row = T.createbegin(); row != T.createend(); ++row) {
354 row.insert(row.index()-1);
355 row.insert(row.index());
356 if (row.index() < T.N() - 1)
357 row.insert(row.index()+1);
359 for (
int row = 0; row < i; ++row) {
361 T[row][row-1] = sqrt(betas[row-1]) / lambdas[row-1];
364 T[row][row] = 1.0 / lambdas[row];
366 T[row][row] += betas[row-1] / lambdas[row-1];
370 T[row][row+1] = sqrt(betas[row]) / lambdas[row];
386 std::cout <<
"Min eigv estimate: " << Simd::io(min_eigv) <<
'\n';
387 std::cout <<
"Max eigv estimate: " << Simd::io(max_eigv) <<
'\n';
388 std::cout <<
"Condition estimate: "
389 << Simd::io(max_eigv / min_eigv) << std::endl;
393 std::cerr <<
"WARNING: Condition estimate was requested. This requires ARPACK, but ARPACK was not found!" << std::endl;
399 bool condition_estimate_ =
false;
442 const Simd::Scalar<real_type> EPSILON=1e-80;
445 field_type rho, rho_new, alpha, beta, h, omega;
466 _op->applyscaleadd(-1,x,r);
471 if(iteration.step(0, norm)){
486 for (it = 0.5; it <
_maxit; it+=.5)
493 rho_new =
_sp->dot(rt,r);
496 if (Simd::allTrue(abs(rho) <= EPSILON))
497 DUNE_THROW(
SolverAbort,
"breakdown in BiCGSTAB - rho "
498 << Simd::io(rho) <<
" <= EPSILON " << EPSILON
499 <<
" after " << it <<
" iterations");
500 if (Simd::allTrue(abs(omega) <= EPSILON))
501 DUNE_THROW(
SolverAbort,
"breakdown in BiCGSTAB - omega "
502 << Simd::io(omega) <<
" <= EPSILON " << EPSILON
503 <<
" after " << it <<
" iterations");
512 ( rho_new / rho ) * ( alpha / omega ));
528 if ( Simd::allTrue(abs(h) < EPSILON) )
529 DUNE_THROW(
SolverAbort,
"abs(h) < EPSILON in BiCGSTAB - abs(h) "
530 << Simd::io(abs(h)) <<
" < EPSILON " << EPSILON
531 <<
" after " << it <<
" iterations");
549 if(iteration.step(it, norm)){
582 if(iteration.step(it, norm)){
597 template<
class CountType>
636 _op->applyscaleadd(-1.0,x,b);
645 if (iteration.step(0, def)){
653 std::array<real_type,2> c{{0.0,0.0}};
655 std::array<field_type,2> s{{0.0,0.0}};
658 std::array<field_type,3> T{{0.0,0.0,0.0}};
661 std::array<field_type,2> xi{{1.0,0.0}};
665 beta = sqrt(
_sp->dot(b,z));
669 std::array<X,3> p{{b,b,b}};
675 std::array<X,3> q{{b,b,b}};
697 q[i2].axpy(-beta,q[i0]);
701 alpha =
_sp->dot(z,q[i2]);
702 q[i2].axpy(-alpha,q[i1]);
705 _prec->apply(z,q[i2]);
709 beta = sqrt(
_sp->dot(q[i2],z));
726 T[2] = c[(i+1)%2]*alpha - s[(i+1)%2]*T[1];
727 T[1] = c[(i+1)%2]*T[1] + s[(i+1)%2]*alpha;
733 generateGivensRotation(T[2],beta,c[i%2],s[i%2]);
735 T[2] = c[i%2]*T[2] + s[i%2]*beta;
737 xi[i%2] = -s[i%2]*xi[(i+1)%2];
738 xi[(i+1)%2] *= c[i%2];
742 p[i2].axpy(-T[1],p[i1]);
743 p[i2].axpy(-T[0],p[i0]);
747 x.axpy(beta0*xi[(i+1)%2],p[i2]);
754 def = abs(beta0*xi[i%2]);
755 if(iteration.step(i, def)){
775 real_type norm_max = max(norm_dx, norm_dy);
776 real_type norm_min = min(norm_dx, norm_dy);
779 cs = Simd::cond(norm_dy < eps,
781 Simd::cond(norm_dx < eps,
783 Simd::cond(norm_dy > norm_dx,
787 sn = Simd::cond(norm_dy < eps,
789 Simd::cond(norm_dx < eps,
791 Simd::cond(norm_dy > norm_dx,
825 template<
class X,
class Y=X,
class F = Y>
926 const Simd::Scalar<real_type> EPSILON = 1e-80;
930 std::vector<field_type,fAlloc> s(m+1), sn(m);
931 std::vector<real_type,rAlloc> cs(m);
936 std::vector< std::vector<field_type,fAlloc> > H(m+1,s);
937 std::vector<F> v(m+1,b);
945 _op->applyscaleadd(-1.0,x,b);
947 v[0] = 0.0;
_prec->apply(v[0],b);
948 norm =
_sp->norm(v[0]);
949 if(iteration.step(0, norm)){
969 _op->apply(v[i],v[i+1]);
970 _prec->apply(w,v[i+1]);
971 for(
int k=0; k<i+1; k++) {
976 H[k][i] =
_sp->dot(v[k],w);
978 w.axpy(-H[k][i],v[k]);
980 H[i+1][i] =
_sp->norm(w);
981 if(Simd::allTrue(abs(H[i+1][i]) < EPSILON))
983 "breakdown in GMRes - |w| == 0.0 after " << j <<
" iterations");
987 v[i+1] *= Simd::cond(norm==
real_type(0.),
992 for(
int k=0; k<i; k++)
1004 iteration.step(j, norm);
1020 std::cout <<
"=== GMRes::restart" << std::endl;
1024 _op->applyscaleadd(-1.0,x,b);
1027 _prec->apply(v[0],b);
1028 norm =
_sp->norm(v[0]);
1040 const std::vector<std::vector<field_type,fAlloc> >& H,
1041 const std::vector<field_type,fAlloc>& s,
1042 const std::vector<X>& v) {
1044 std::vector<field_type,fAlloc> y(s);
1047 for(
int a=i-1; a>=0; a--) {
1049 for(
int b=a+1; b<i; b++)
1050 rhs -= H[a][b]*y[b];
1061 template<
typename T>
1062 typename std::enable_if<std::is_same<field_type,real_type>::value,T>::type
conjugate(
const T& t) {
1066 template<
typename T>
1067 typename std::enable_if<!std::is_same<field_type,real_type>::value,T>::type
conjugate(
const T& t) {
1082 real_type norm_max = max(norm_dx, norm_dy);
1083 real_type norm_min = min(norm_dx, norm_dy);
1086 cs = Simd::cond(norm_dy < eps,
1088 Simd::cond(norm_dx < eps,
1090 Simd::cond(norm_dy > norm_dx,
1094 sn = Simd::cond(norm_dy < eps,
1096 Simd::cond(norm_dx < eps,
1098 Simd::cond(norm_dy > norm_dx,
1137 template<
class X,
class Y=X,
class F = Y>
1172 const Simd::Scalar<real_type> EPSILON = 1e-80;
1176 std::vector<field_type,fAlloc> s(m+1), sn(m);
1177 std::vector<real_type,rAlloc> cs(m);
1180 std::vector< std::vector<field_type,fAlloc> > H(m+1,s);
1181 std::vector<F> v(m+1,b);
1182 std::vector<X> w(m+1,b);
1184 Iteration iteration(*
this,res);
1190 _op->applyscaleadd(-1.0, x, v[0]);
1192 norm =
_sp->norm(v[0]);
1193 if(iteration.step(0, norm)){
1202 v[0] *= (1.0 / norm);
1204 for(i=1; i<m+1; ++i)
1212 _prec->apply(w[i], v[i]);
1215 _op->apply(w[i], v[i+1]);
1217 for(
int kk=0; kk<i+1; kk++)
1223 H[kk][i] =
_sp->dot(v[kk],v[i+1]);
1225 v[i+1].axpy(-H[kk][i], v[kk]);
1227 H[i+1][i] =
_sp->norm(v[i+1]);
1228 if(Simd::allTrue(abs(H[i+1][i]) < EPSILON))
1229 DUNE_THROW(
SolverAbort,
"breakdown in fGMRes - |w| (-> "
1230 << w[i] <<
") == 0.0 after "
1231 << j <<
" iterations");
1250 iteration.step(j, norm);
1255 this->
update(tmp, i, H, s, w);
1265 std::cout <<
"=== fGMRes::restart" << std::endl;
1269 _op->applyscaleadd(-1.0, x,v[0]);
1271 norm =
_sp->norm(v[0]);
1363 _restart(configuration.
get<int>(
"restart"))
1368 _restart(configuration.
get<int>(
"restart"))
1393 Iteration iteration(*
this, res);
1395 _op->applyscaleadd(-1,x,b);
1397 std::vector<std::shared_ptr<X> > p(_restart);
1398 std::vector<field_type,fAlloc> pp(_restart);
1402 p[0].reset(
new X(x));
1405 if(iteration.step(0, def)){
1416 _prec->apply(*(p[0]),b);
1417 rho =
_sp->dot(*(p[0]),b);
1418 _op->apply(*(p[0]),q);
1419 pp[0] =
_sp->dot(*(p[0]),q);
1421 x.axpy(lambda,*(p[0]));
1427 if(iteration.step(i, def)){
1434 int end=std::min(_restart,
_maxit-i+1);
1435 for (ii=1; ii<end; ++ii )
1440 _prec->apply(prec_res,b);
1442 p[ii].reset(
new X(prec_res));
1443 _op->apply(prec_res, q);
1445 for(
int j=0; j<ii; ++j) {
1446 rho =
_sp->dot(q,*(p[j]))/pp[j];
1447 p[ii]->axpy(-rho, *(p[j]));
1451 _op->apply(*(p[ii]),q);
1452 pp[ii] =
_sp->dot(*(p[ii]),q);
1453 rho =
_sp->dot(*(p[ii]),b);
1454 lambda = rho/pp[ii];
1455 x.axpy(lambda,*(p[ii]));
1462 iteration.step(i, def);
1467 *(p[0])=*(p[_restart-1]);
1468 pp[0]=pp[_restart-1];
1530 scalar_real_type reduction,
int maxit,
int verbose,
int mmax = 10) :
IterativeSolver<X,X>(op, sp, prec, reduction, maxit, verbose),
_mmax(mmax)
1561 const ParameterTree& config)
1568 const ParameterTree& config)
1590 _op->applyscaleadd(-1,x,b);
1593 std::vector<X> d(
_mmax+1, x);
1594 std::vector<X> Ad(
_mmax+1, x);
1595 std::vector<field_type,rAlloc> ddotAd(
_mmax+1,0);
1599 if(iteration.step(0, def)){
1611 for (; i_bounded <=
_mmax && i<=
_maxit; i_bounded++) {
1613 _prec->apply(d[i_bounded], b);
1616 orthogonalizations(i_bounded,Ad,w,ddotAd,d);
1619 _op->apply(d[i_bounded], Ad[i_bounded]);
1620 ddotAd[i_bounded]=
_sp->dot(d[i_bounded],Ad[i_bounded]);
1621 alpha =
_sp->dot(d[i_bounded], b)/ddotAd[i_bounded];
1624 x.axpy(alpha, d[i_bounded]);
1625 b.axpy(-alpha, Ad[i_bounded]);
1630 iteration.step(i, def);
1634 cycle(Ad,d,ddotAd,i_bounded);
1647 for (
int k = 0; k < i_bounded; k++) {
1648 d[i_bounded].axpy(-
_sp->dot(Ad[k], w) / ddotAd[k], d[k]);
1653 virtual void cycle(std::vector<X>& Ad,std::vector<X>& d,std::vector<
field_type,ReboundAllocatorType<X,field_type> >& ddotAd,
int& i_bounded) {
1656 std::swap(Ad[0], Ad[
_mmax]);
1657 std::swap(d[0], d[
_mmax]);
1658 std::swap(ddotAd[0], ddotAd[
_mmax]);
1704 for (
int k = 0; k < _k_limit; k++) {
1706 d[i_bounded].axpy(-
_sp->dot(Ad[k], w) / ddotAd[k], d[k]);
1709 if(_k_limit<=i_bounded)
1715 virtual void cycle(std::vector<X>& Ad, [[maybe_unused]] std::vector<X>& d, [[maybe_unused]] std::vector<
field_type,ReboundAllocatorType<X,field_type> >& ddotAd,
int& i_bounded)
override {
1719 _k_limit = Ad.size();
Define general, extensible interface for inverse operators.
#define DUNE_REGISTER_ITERATIVE_SOLVER(name,...)
Definition solverregistry.hh:19
Implementation of the BCRSMatrix class.
Define base class for scalar product and norm.
Define general, extensible interface for operators. The available implementation wraps a matrix.
Definition allocator.hh:11
PropertyMapTypeSelector< Amg::VertexVisitedTag, Amg::PropertiesGraph< G, Amg::VertexProperties, EP, VM, EM > >::Type get(const Amg::VertexVisitedTag &tag, Amg::PropertiesGraph< G, Amg::VertexProperties, EP, VM, EM > &graph)
Definition dependency.hh:293
typename std::allocator_traits< typename AllocatorTraits< T >::type >::template rebind_alloc< X > ReboundAllocatorType
Definition allocator.hh:37
A sparse block matrix with compressed row storage.
Definition bcrsmatrix.hh:466
@ row_wise
Build in a row-wise manner.
Definition bcrsmatrix.hh:521
A vector of blocks with memory management.
Definition bvector.hh:395
Wrapper to use a range of ARPACK++ eigenvalue solvers.
Definition arpackpp.hh:245
void computeSymMaxMagnitude(const Real &epsilon, BlockVector &x, Real &lambda) const
Assume the matrix to be square, symmetric and perform IRLM to compute an approximation lambda of its ...
Definition arpackpp.hh:289
void computeSymMinMagnitude(const Real &epsilon, BlockVector &x, Real &lambda) const
Assume the matrix to be square, symmetric and perform IRLM to compute an approximation lambda of its ...
Definition arpackpp.hh:391
Thrown when a solver aborts due to some problem.
Definition istlexception.hh:46
A linear operator.
Definition operators.hh:67
Base class for matrix free definition of preconditioners.
Definition preconditioner.hh:32
Base class for scalar product and norm computation.
Definition scalarproducts.hh:52
Statistics about the application of an inverse operator.
Definition solver.hh:48
double condition_estimate
Estimate of condition number.
Definition solver.hh:79
void clear()
Resets all data.
Definition solver.hh:56
bool converged
True if convergence criterion has been met.
Definition solver.hh:73
Simd::Scalar< real_type > scalar_real_type
scalar type underlying the field_type
Definition solver.hh:114
Y range_type
Type of the range of the operator to be inverted.
Definition solver.hh:105
X domain_type
Type of the domain of the operator to be inverted.
Definition solver.hh:102
X::field_type field_type
The field type of the operator.
Definition solver.hh:108
FieldTraits< field_type >::real_type real_type
The real type of the field type (is the same if using real numbers, but differs for std::complex)
Definition solver.hh:111
Base class for all implementations of iterative solvers.
Definition solver.hh:203
std::shared_ptr< const ScalarProduct< X > > _sp
Definition solver.hh:506
IterativeSolver(const LinearOperator< X, X > &op, Preconditioner< X, X > &prec, scalar_real_type reduction, int maxit, int verbose)
General constructor to initialize an iterative solver.
Definition solver.hh:230
std::shared_ptr< const LinearOperator< X, X > > _op
Definition solver.hh:504
int _maxit
Definition solver.hh:508
int _verbose
Definition solver.hh:509
scalar_real_type _reduction
Definition solver.hh:507
std::shared_ptr< Preconditioner< X, X > > _prec
Definition solver.hh:505
Preconditioned loop solver.
Definition solvers.hh:59
virtual void apply(X &x, X &b, InverseOperatorResult &res)
Apply inverse operator,.
Definition solvers.hh:73
typename IterativeSolver< X, X >::template Iteration< unsigned int > Iteration
Definition solvers.hh:116
gradient method
Definition solvers.hh:124
virtual void apply(X &x, X &b, InverseOperatorResult &res)
Apply inverse operator.
Definition solvers.hh:142
typename IterativeSolver< X, X >::template Iteration< unsigned int > Iteration
Definition solvers.hh:187
conjugate gradient method
Definition solvers.hh:193
CGSolver(std::shared_ptr< const LinearOperator< X, X > > op, std::shared_ptr< ScalarProduct< X > > sp, std::shared_ptr< Preconditioner< X, X > > prec, scalar_real_type reduction, int maxit, int verbose, bool condition_estimate)
Constructor to initialize a CG solver.
Definition solvers.hh:256
static constexpr bool enableConditionEstimate
Definition solvers.hh:208
CGSolver(const LinearOperator< X, X > &op, const ScalarProduct< X > &sp, Preconditioner< X, X > &prec, scalar_real_type reduction, int maxit, int verbose, bool condition_estimate)
Constructor to initialize a CG solver.
Definition solvers.hh:239
virtual void apply(X &x, X &b, InverseOperatorResult &res)
Apply inverse operator.
Definition solvers.hh:279
CGSolver(const LinearOperator< X, X > &op, Preconditioner< X, X > &prec, scalar_real_type reduction, int maxit, int verbose, bool condition_estimate)
Constructor to initialize a CG solver.
Definition solvers.hh:222
typename IterativeSolver< X, X >::template Iteration< unsigned int > Iteration
Definition solvers.hh:412
Bi-conjugate Gradient Stabilized (BiCG-STAB)
Definition solvers.hh:419
typename IterativeSolver< X, X >::template Iteration< CountType > Iteration
Definition solvers.hh:598
virtual void apply(X &x, X &b, InverseOperatorResult &res)
Apply inverse operator.
Definition solvers.hh:439
Minimal Residual Method (MINRES)
Definition solvers.hh:609
virtual void apply(X &x, X &b, InverseOperatorResult &res)
Apply inverse operator.
Definition solvers.hh:627
typename IterativeSolver< X, X >::template Iteration< unsigned int > Iteration
Definition solvers.hh:808
implements the Generalized Minimal Residual (GMRes) method
Definition solvers.hh:827
RestartedGMResSolver(std::shared_ptr< const LinearOperator< X, Y > > op, std::shared_ptr< Preconditioner< X, X > > prec, const ParameterTree &configuration)
Constructor.
Definition solvers.hh:878
std::enable_if<!std::is_same< field_type, real_type >::value, T >::type conjugate(const T &t)
Definition solvers.hh:1067
RestartedGMResSolver(std::shared_ptr< const LinearOperator< X, Y > > op, std::shared_ptr< const ScalarProduct< X > > sp, std::shared_ptr< Preconditioner< X, X > > prec, const ParameterTree &configuration)
Definition solvers.hh:883
void update(X &w, int i, const std::vector< std::vector< field_type, fAlloc > > &H, const std::vector< field_type, fAlloc > &s, const std::vector< X > &v)
Definition solvers.hh:1039
std::enable_if< std::is_same< field_type, real_type >::value, T >::type conjugate(const T &t)
Definition solvers.hh:1062
ReboundAllocatorType< X, field_type > fAlloc
field_type Allocator retrieved from domain type
Definition solvers.hh:838
int _restart
Definition solvers.hh:1120
ReboundAllocatorType< X, real_type > rAlloc
real_type Allocator retrieved from domain type
Definition solvers.hh:840
virtual void apply(X &x, Y &b, double reduction, InverseOperatorResult &res)
Apply inverse operator.
Definition solvers.hh:923
RestartedGMResSolver(const LinearOperator< X, Y > &op, Preconditioner< X, Y > &prec, scalar_real_type reduction, int restart, int maxit, int verbose)
Set up RestartedGMResSolver solver.
Definition solvers.hh:850
RestartedGMResSolver(const LinearOperator< X, Y > &op, const ScalarProduct< X > &sp, Preconditioner< X, Y > &prec, scalar_real_type reduction, int restart, int maxit, int verbose)
Set up RestartedGMResSolver solver.
Definition solvers.hh:861
void generatePlaneRotation(field_type &dx, field_type &dy, real_type &cs, field_type &sn)
Definition solvers.hh:1073
RestartedGMResSolver(std::shared_ptr< const LinearOperator< X, Y > > op, std::shared_ptr< const ScalarProduct< X > > sp, std::shared_ptr< Preconditioner< X, Y > > prec, scalar_real_type reduction, int restart, int maxit, int verbose)
Set up RestartedGMResSolver solver.
Definition solvers.hh:894
virtual void apply(X &x, Y &b, InverseOperatorResult &res)
Apply inverse operator.
Definition solvers.hh:910
typename IterativeSolver< X, X >::template Iteration< unsigned int > Iteration
Definition solvers.hh:1119
void applyPlaneRotation(field_type &dx, field_type &dy, real_type &cs, field_type &sn)
Definition solvers.hh:1106
implements the Flexible Generalized Minimal Residual (FGMRes) method (right preconditioned)
Definition solvers.hh:1139
void apply(X &x, Y &b, double reduction, InverseOperatorResult &res) override
Apply inverse operator.
Definition solvers.hh:1169
Generalized preconditioned conjugate gradient solver.
Definition solvers.hh:1307
GeneralizedPCGSolver(const LinearOperator< X, X > &op, const ScalarProduct< X > &sp, Preconditioner< X, X > &prec, scalar_real_type reduction, int maxit, int verbose, int restart=10)
Set up nonlinear preconditioned conjugate gradient solver.
Definition solvers.hh:1343
GeneralizedPCGSolver(const LinearOperator< X, X > &op, Preconditioner< X, X > &prec, scalar_real_type reduction, int maxit, int verbose, int restart=10)
Set up nonlinear preconditioned conjugate gradient solver.
Definition solvers.hh:1331
virtual void apply(X &x, X &b, InverseOperatorResult &res)
Apply inverse operator.
Definition solvers.hh:1391
GeneralizedPCGSolver(std::shared_ptr< const LinearOperator< X, X > > op, std::shared_ptr< Preconditioner< X, X > > prec, const ParameterTree &configuration)
Constructor.
Definition solvers.hh:1361
GeneralizedPCGSolver(std::shared_ptr< const LinearOperator< X, X > > op, std::shared_ptr< const ScalarProduct< X > > sp, std::shared_ptr< Preconditioner< X, X > > prec, scalar_real_type reduction, int maxit, int verbose, int restart=10)
Set up nonlinear preconditioned conjugate gradient solver.
Definition solvers.hh:1377
GeneralizedPCGSolver(std::shared_ptr< const LinearOperator< X, X > > op, std::shared_ptr< const ScalarProduct< X > > sp, std::shared_ptr< Preconditioner< X, X > > prec, const ParameterTree &configuration)
Definition solvers.hh:1366
Accelerated flexible conjugate gradient method.
Definition solvers.hh:1501
RestartedFCGSolver(const LinearOperator< X, X > &op, Preconditioner< X, X > &prec, scalar_real_type reduction, int maxit, int verbose, int mmax=10)
Constructor to initialize a RestartedFCG solver.
Definition solvers.hh:1519
RestartedFCGSolver(std::shared_ptr< const LinearOperator< X, X > > op, std::shared_ptr< Preconditioner< X, X > > prec, const ParameterTree &config)
Constructor.
Definition solvers.hh:1559
RestartedFCGSolver(std::shared_ptr< const LinearOperator< X, X > > op, std::shared_ptr< const ScalarProduct< X > > sp, std::shared_ptr< Preconditioner< X, X > > prec, scalar_real_type reduction, int maxit, int verbose, int mmax=10)
Constructor to initialize a RestartedFCG solver.
Definition solvers.hh:1539
int _mmax
Definition solvers.hh:1662
RestartedFCGSolver(std::shared_ptr< const LinearOperator< X, X > > op, std::shared_ptr< const ScalarProduct< X > > sp, std::shared_ptr< Preconditioner< X, X > > prec, const ParameterTree &config)
Definition solvers.hh:1565
typename IterativeSolver< X, X >::template Iteration< unsigned int > Iteration
Definition solvers.hh:1669
virtual void apply(X &x, X &b, InverseOperatorResult &res)
Apply inverse operator.
Definition solvers.hh:1584
RestartedFCGSolver(const LinearOperator< X, X > &op, const ScalarProduct< X > &sp, Preconditioner< X, X > &prec, scalar_real_type reduction, int maxit, int verbose, int mmax=10)
Constructor to initialize a RestartedFCG solver.
Definition solvers.hh:1529
Complete flexible conjugate gradient method.
Definition solvers.hh:1680
virtual void apply(X &x, X &b, InverseOperatorResult &res) override
Apply inverse operator.
Definition solvers.hh:1694