4#ifndef DUNE_ISTL_SOLVERS_HH
5#define DUNE_ISTL_SOLVERS_HH
15#include <dune/common/exceptions.hh>
16#include <dune/common/math.hh>
17#include <dune/common/simd/io.hh>
18#include <dune/common/simd/simd.hh>
19#include <dune/common/std/type_traits.hh>
20#include <dune/common/timer.hh>
77 _op->applyscaleadd(-1,x,b);
81 if(iteration.step(0, def)){
97 _op->applyscaleadd(-1,v,b);
99 if(iteration.step(i, def))
145 _op->applyscaleadd(-1,x,b);
148 if(iteration.step(0, def)){
163 lambda =
_sp->dot(p,b)/
_sp->dot(q,p);
168 if(iteration.step(i, def))
205 = std::bool_constant<enableConditionEstimate>;
221 condition_estimate_(condition_estimate)
224 condition_estimate_ =
false;
225 std::cerr <<
"WARNING: Condition estimate was disabled. It is only available for double and float field types!" << std::endl;
237 scalar_real_type reduction,
int maxit,
int verbose,
bool condition_estimate) :
IterativeSolver<X,X>(op, sp, prec, reduction, maxit, verbose),
238 condition_estimate_(condition_estimate)
240 if (condition_estimate && !(std::is_same<field_type,float>::value || std::is_same<field_type,double>::value)) {
241 condition_estimate_ =
false;
242 std::cerr <<
"WARNING: Condition estimate was disabled. It is only available for double and float field types!" << std::endl;
255 scalar_real_type reduction,
int maxit,
int verbose,
bool condition_estimate)
257 condition_estimate_(condition_estimate)
259 if (condition_estimate && !(std::is_same<field_type,float>::value || std::is_same<field_type,double>::value)) {
260 condition_estimate_ =
false;
261 std::cerr <<
"WARNING: Condition estimate was disabled. It is only available for double and float field types!" << std::endl;
281 _op->applyscaleadd(-1,x,b);
284 if(iteration.step(0, def)){
293 std::vector<real_type> lambdas(0);
294 std::vector<real_type> betas(0);
302 rholast =
_sp->dot(p,b);
310 alpha =
_sp->dot(p,q);
311 lambda = rholast/alpha;
313 if (condition_estimate_)
314 lambdas.push_back(std::real(lambda));
320 if(iteration.step(i, def))
329 if (condition_estimate_)
330 betas.push_back(std::real(beta));
338 if (condition_estimate_) {
351 row.insert(row.index()-1);
352 row.insert(row.index());
353 if (row.index() < T.
N() - 1)
354 row.insert(row.index()+1);
356 for (
int row = 0; row < i; ++row) {
358 T[row][row-1] = sqrt(betas[row-1]) / lambdas[row-1];
361 T[row][row] = 1.0 / lambdas[row];
363 T[row][row] += betas[row-1] / lambdas[row-1];
367 T[row][row+1] = sqrt(betas[row]) / lambdas[row];
383 std::cout <<
"Min eigv estimate: " << Simd::io(min_eigv) <<
'\n';
384 std::cout <<
"Max eigv estimate: " << Simd::io(max_eigv) <<
'\n';
385 std::cout <<
"Condition estimate: "
386 << Simd::io(max_eigv / min_eigv) << std::endl;
390 std::cerr <<
"WARNING: Condition estimate was requested. This requires ARPACK, but ARPACK was not found!" << std::endl;
396 bool condition_estimate_ =
false;
439 const Simd::Scalar<real_type> EPSILON=1e-80;
442 field_type rho, rho_new, alpha, beta, h, omega;
463 _op->applyscaleadd(-1,x,r);
468 if(iteration.step(0, norm)){
483 for (it = 0.5; it <
_maxit; it+=.5)
490 rho_new =
_sp->dot(rt,r);
493 if (Simd::allTrue(abs(rho) <= EPSILON))
494 DUNE_THROW(
SolverAbort,
"breakdown in BiCGSTAB - rho "
495 << Simd::io(rho) <<
" <= EPSILON " << EPSILON
496 <<
" after " << it <<
" iterations");
497 if (Simd::allTrue(abs(omega) <= EPSILON))
498 DUNE_THROW(
SolverAbort,
"breakdown in BiCGSTAB - omega "
499 << Simd::io(omega) <<
" <= EPSILON " << EPSILON
500 <<
" after " << it <<
" iterations");
507 beta = ( rho_new / rho ) * ( alpha / omega );
523 if ( Simd::allTrue(abs(h) < EPSILON) )
524 DUNE_THROW(
SolverAbort,
"abs(h) < EPSILON in BiCGSTAB - abs(h) "
525 << Simd::io(abs(h)) <<
" < EPSILON " << EPSILON
526 <<
" after " << it <<
" iterations");
542 if(iteration.step(it, norm)){
556 omega =
_sp->dot(t,r)/
_sp->dot(t,t);
572 if(iteration.step(it, norm)){
587 template<
class CountType>
626 _op->applyscaleadd(-1,x,b);
630 if(iteration.step(0, def)){
638 std::array<real_type,2> c{{0.0,0.0}};
640 std::array<field_type,2> s{{0.0,0.0}};
643 std::array<field_type,3> T{{0.0,0.0,0.0}};
646 std::array<field_type,2> xi{{1.0,0.0}};
657 beta = sqrt(
_sp->dot(b,z));
661 std::array<X,3> p{{b,b,b}};
667 std::array<X,3> q{{b,b,b}};
685 q[i2].axpy(-beta,q[i0]);
689 alpha =
_sp->dot(z,q[i2]);
690 q[i2].axpy(-alpha,q[i1]);
693 _prec->apply(z,q[i2]);
697 beta = sqrt(
_sp->dot(q[i2],z));
710 T[2] = c[(i+1)%2]*alpha - s[(i+1)%2]*T[1];
711 T[1] = c[(i+1)%2]*T[1] + s[(i+1)%2]*alpha;
717 generateGivensRotation(T[2],beta,c[i%2],s[i%2]);
719 T[2] = c[i%2]*T[2] + s[i%2]*beta;
721 xi[i%2] = -s[i%2]*xi[(i+1)%2];
722 xi[(i+1)%2] *= c[i%2];
726 p[i2].axpy(-T[1],p[i1]);
727 p[i2].axpy(-T[0],p[i0]);
731 x.axpy(beta0*xi[(i+1)%2],p[i2]);
738 def = abs(beta0*xi[i%2]);
739 if(iteration.step(i, def)){
759 real_type norm_max = max(norm_dx, norm_dy);
760 real_type norm_min = min(norm_dx, norm_dy);
763 cs = Simd::cond(norm_dy < eps,
765 Simd::cond(norm_dx < eps,
767 Simd::cond(norm_dy > norm_dx,
771 sn = Simd::cond(norm_dy < eps,
773 Simd::cond(norm_dx < eps,
775 Simd::cond(norm_dy > norm_dx,
786 using IterativeSolver<X,X>
::_op;
787 using IterativeSolver<X,X>
::_prec;
788 using IterativeSolver<X,X>
::_sp;
809 template<
class X,
class Y=X,
class F = Y>
910 const Simd::Scalar<real_type> EPSILON = 1e-80;
914 std::vector<field_type,fAlloc> s(m+1), sn(m);
915 std::vector<real_type,rAlloc> cs(m);
920 std::vector< std::vector<field_type,fAlloc> > H(m+1,s);
921 std::vector<F> v(m+1,b);
929 _op->applyscaleadd(-1.0,x,b);
931 v[0] = 0.0;
_prec->apply(v[0],b);
932 norm =
_sp->norm(v[0]);
933 if(iteration.step(0, norm)){
951 _op->apply(v[i],v[i+1]);
952 _prec->apply(w,v[i+1]);
953 for(
int k=0; k<i+1; k++) {
958 H[k][i] =
_sp->dot(v[k],w);
960 w.axpy(-H[k][i],v[k]);
962 H[i+1][i] =
_sp->norm(w);
963 if(Simd::allTrue(abs(H[i+1][i]) < EPSILON))
965 "breakdown in GMRes - |w| == 0.0 after " << j <<
" iterations");
968 v[i+1] = w; v[i+1] *=
real_type(1.0)/H[i+1][i];
971 for(
int k=0; k<i; k++)
983 iteration.step(j, norm);
999 std::cout <<
"=== GMRes::restart" << std::endl;
1003 _op->applyscaleadd(-1.0,x,b);
1006 _prec->apply(v[0],b);
1007 norm =
_sp->norm(v[0]);
1019 const std::vector<std::vector<field_type,fAlloc> >& H,
1020 const std::vector<field_type,fAlloc>& s,
1021 const std::vector<X>& v) {
1023 std::vector<field_type,fAlloc> y(s);
1026 for(
int a=i-1; a>=0; a--) {
1028 for(
int b=a+1; b<i; b++)
1029 rhs -= H[a][b]*y[b];
1038 template<
typename T>
1039 typename std::enable_if<std::is_same<field_type,real_type>::value,T>::type
conjugate(
const T& t) {
1043 template<
typename T>
1044 typename std::enable_if<!std::is_same<field_type,real_type>::value,T>::type
conjugate(
const T& t) {
1059 real_type norm_max = max(norm_dx, norm_dy);
1060 real_type norm_min = min(norm_dx, norm_dy);
1063 cs = Simd::cond(norm_dy < eps,
1065 Simd::cond(norm_dx < eps,
1067 Simd::cond(norm_dy > norm_dx,
1071 sn = Simd::cond(norm_dy < eps,
1073 Simd::cond(norm_dx < eps,
1075 Simd::cond(norm_dy > norm_dx,
1114 template<
class X,
class Y=X,
class F = Y>
1149 const Simd::Scalar<real_type> EPSILON = 1e-80;
1153 std::vector<field_type,fAlloc> s(m+1), sn(m);
1154 std::vector<real_type,rAlloc> cs(m);
1157 std::vector< std::vector<field_type,fAlloc> > H(m+1,s);
1158 std::vector<F> v(m+1,b);
1159 std::vector<X> w(m+1,b);
1161 Iteration iteration(*
this,res);
1167 _op->applyscaleadd(-1.0, x, v[0]);
1169 norm =
_sp->norm(v[0]);
1170 if(iteration.step(0, norm)){
1179 v[0] *= (1.0 / norm);
1181 for(i=1; i<m+1; ++i)
1189 _prec->apply(w[i], v[i]);
1192 _op->apply(w[i], v[i+1]);
1194 for(
int kk=0; kk<i+1; kk++)
1200 H[kk][i] =
_sp->dot(v[kk],v[i+1]);
1202 v[i+1].axpy(-H[kk][i], v[kk]);
1204 H[i+1][i] =
_sp->norm(v[i+1]);
1205 if(Simd::allTrue(abs(H[i+1][i]) < EPSILON))
1206 DUNE_THROW(
SolverAbort,
"breakdown in fGMRes - |w| (-> "
1207 << w[i] <<
") == 0.0 after "
1208 << j <<
" iterations");
1227 iteration.step(j, norm);
1232 this->
update(tmp, i, H, s, w);
1242 std::cout <<
"=== fGMRes::restart" << std::endl;
1246 _op->applyscaleadd(-1.0, x,v[0]);
1248 norm =
_sp->norm(v[0]);
1340 _restart(configuration.
get<int>(
"restart"))
1345 _restart(configuration.
get<int>(
"restart"))
1370 Iteration iteration(*
this, res);
1372 _op->applyscaleadd(-1,x,b);
1374 std::vector<std::shared_ptr<X> > p(_restart);
1375 std::vector<field_type,fAlloc> pp(_restart);
1379 p[0].reset(
new X(x));
1382 if(iteration.step(0, def)){
1393 _prec->apply(*(p[0]),b);
1394 rho =
_sp->dot(*(p[0]),b);
1395 _op->apply(*(p[0]),q);
1396 pp[0] =
_sp->dot(*(p[0]),q);
1398 x.axpy(lambda,*(p[0]));
1404 if(iteration.step(i, def)){
1411 int end=std::min(_restart,
_maxit-i+1);
1412 for (ii=1; ii<end; ++ii )
1417 _prec->apply(prec_res,b);
1419 p[ii].reset(
new X(prec_res));
1420 _op->apply(prec_res, q);
1422 for(
int j=0; j<ii; ++j) {
1423 rho =
_sp->dot(q,*(p[j]))/pp[j];
1424 p[ii]->axpy(-rho, *(p[j]));
1428 _op->apply(*(p[ii]),q);
1429 pp[ii] =
_sp->dot(*(p[ii]),q);
1430 rho =
_sp->dot(*(p[ii]),b);
1431 lambda = rho/pp[ii];
1432 x.axpy(lambda,*(p[ii]));
1439 iteration.step(i, def);
1444 *(p[0])=*(p[_restart-1]);
1445 pp[0]=pp[_restart-1];
1507 scalar_real_type reduction,
int maxit,
int verbose,
int mmax = 10) :
IterativeSolver<X,X>(op, sp, prec, reduction, maxit, verbose),
_mmax(mmax)
1538 const ParameterTree& config)
1545 const ParameterTree& config)
1567 _op->applyscaleadd(-1,x,b);
1570 std::vector<X> d(
_mmax+1, x);
1571 std::vector<X> Ad(
_mmax+1, x);
1572 std::vector<field_type,rAlloc> ddotAd(
_mmax+1,0);
1576 if(iteration.step(0, def)){
1588 for (; i_bounded <=
_mmax && i<=
_maxit; i_bounded++) {
1590 _prec->apply(d[i_bounded], b);
1593 orthogonalizations(i_bounded,Ad,w,ddotAd,d);
1596 _op->apply(d[i_bounded], Ad[i_bounded]);
1597 ddotAd[i_bounded]=
_sp->dot(d[i_bounded],Ad[i_bounded]);
1598 alpha =
_sp->dot(d[i_bounded], b)/ddotAd[i_bounded];
1601 x.axpy(alpha, d[i_bounded]);
1602 b.axpy(-alpha, Ad[i_bounded]);
1607 iteration.step(i, def);
1611 cycle(Ad,d,ddotAd,i_bounded);
1624 for (
int k = 0; k < i_bounded; k++) {
1625 d[i_bounded].axpy(-
_sp->dot(Ad[k], w) / ddotAd[k], d[k]);
1630 virtual void cycle(std::vector<X>& Ad,std::vector<X>& d,std::vector<
field_type,ReboundAllocatorType<X,field_type> >& ddotAd,
int& i_bounded) {
1633 std::swap(Ad[0], Ad[
_mmax]);
1634 std::swap(d[0], d[
_mmax]);
1635 std::swap(ddotAd[0], ddotAd[
_mmax]);
1681 for (
int k = 0; k < _k_limit; k++) {
1683 d[i_bounded].axpy(-
_sp->dot(Ad[k], w) / ddotAd[k], d[k]);
1686 if(_k_limit<=i_bounded)
1692 virtual void cycle(std::vector<X>& Ad,std::vector<X>& d,std::vector<
field_type,ReboundAllocatorType<X,field_type> >& ddotAd,
int& i_bounded)
override {
1696 _k_limit = Ad.size();
Define general, extensible interface for inverse operators.
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.
#define DUNE_REGISTER_ITERATIVE_SOLVER(name,...)
Definition: solverregistry.hh:17
Definition: allocator.hh:9
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:291
typename std::allocator_traits< typename AllocatorTraits< T >::type >::template rebind_alloc< X > ReboundAllocatorType
Definition: allocator.hh:35
A sparse block matrix with compressed row storage.
Definition: bcrsmatrix.hh:464
@ row_wise
Build in a row-wise manner.
Definition: bcrsmatrix.hh:519
CreateIterator createend()
get create iterator pointing to one after the last block
Definition: bcrsmatrix.hh:1101
CreateIterator createbegin()
get initial create iterator
Definition: bcrsmatrix.hh:1095
size_type N() const
number of rows (counted in blocks)
Definition: bcrsmatrix.hh:1970
A vector of blocks with memory management.
Definition: bvector.hh:393
Wrapper to use a range of ARPACK++ eigenvalue solvers.
Definition: arpackpp.hh:243
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:287
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:389
Thrown when a solver aborts due to some problem.
Definition: istlexception.hh:53
A linear operator.
Definition: operators.hh:65
Base class for matrix free definition of preconditioners.
Definition: preconditioner.hh:30
Base class for scalar product and norm computation.
Definition: scalarproducts.hh:50
Statistics about the application of an inverse operator.
Definition: solver.hh:46
double condition_estimate
Estimate of condition number.
Definition: solver.hh:77
void clear()
Resets all data.
Definition: solver.hh:54
bool converged
True if convergence criterion has been met.
Definition: solver.hh:71
Simd::Scalar< real_type > scalar_real_type
scalar type underlying the field_type
Definition: solver.hh:112
Y range_type
Type of the range of the operator to be inverted.
Definition: solver.hh:103
X domain_type
Type of the domain of the operator to be inverted.
Definition: solver.hh:100
X::field_type field_type
The field type of the operator.
Definition: solver.hh:106
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:109
Base class for all implementations of iterative solvers.
Definition: solver.hh:201
std::shared_ptr< LinearOperator< X, X > > _op
Definition: solver.hh:502
std::shared_ptr< ScalarProduct< X > > _sp
Definition: solver.hh:504
int _maxit
Definition: solver.hh:506
int _verbose
Definition: solver.hh:507
scalar_real_type _reduction
Definition: solver.hh:505
std::shared_ptr< Preconditioner< X, X > > _prec
Definition: solver.hh:503
Preconditioned loop solver.
Definition: solvers.hh:57
virtual void apply(X &x, X &b, InverseOperatorResult &res)
Apply inverse operator,.
Definition: solvers.hh:71
typename IterativeSolver< X, X >::template Iteration< unsigned int > Iteration
Definition: solvers.hh:114
gradient method
Definition: solvers.hh:122
virtual void apply(X &x, X &b, InverseOperatorResult &res)
Apply inverse operator.
Definition: solvers.hh:140
typename IterativeSolver< X, X >::template Iteration< unsigned int > Iteration
Definition: solvers.hh:182
conjugate gradient method
Definition: solvers.hh:188
static constexpr bool enableConditionEstimate
Definition: solvers.hh:203
CGSolver(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:219
std::bool_constant< enableConditionEstimate > enableConditionEstimate_t
Definition: solvers.hh:205
CGSolver(LinearOperator< X, X > &op, 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:236
virtual void apply(X &x, X &b, InverseOperatorResult &res)
Apply inverse operator.
Definition: solvers.hh:276
CGSolver(std::shared_ptr< 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:253
typename IterativeSolver< X, X >::template Iteration< unsigned int > Iteration
Definition: solvers.hh:409
Bi-conjugate Gradient Stabilized (BiCG-STAB)
Definition: solvers.hh:416
typename IterativeSolver< X, X >::template Iteration< CountType > Iteration
Definition: solvers.hh:588
virtual void apply(X &x, X &b, InverseOperatorResult &res)
Apply inverse operator.
Definition: solvers.hh:436
Minimal Residual Method (MINRES)
Definition: solvers.hh:599
virtual void apply(X &x, X &b, InverseOperatorResult &res)
Apply inverse operator.
Definition: solvers.hh:617
typename IterativeSolver< X, X >::template Iteration< unsigned int > Iteration
Definition: solvers.hh:792
implements the Generalized Minimal Residual (GMRes) method
Definition: solvers.hh:811
RestartedGMResSolver(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:834
std::enable_if<!std::is_same< field_type, real_type >::value, T >::type conjugate(const T &t)
Definition: solvers.hh:1044
RestartedGMResSolver(LinearOperator< X, Y > &op, ScalarProduct< X > &sp, Preconditioner< X, Y > &prec, scalar_real_type reduction, int restart, int maxit, int verbose)
Set up RestartedGMResSolver solver.
Definition: solvers.hh:845
RestartedGMResSolver(std::shared_ptr< LinearOperator< X, Y > > op, std::shared_ptr< 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:878
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:1018
RestartedGMResSolver(std::shared_ptr< LinearOperator< X, Y > > op, std::shared_ptr< Preconditioner< X, X > > prec, const ParameterTree &configuration)
Constructor.
Definition: solvers.hh:862
std::enable_if< std::is_same< field_type, real_type >::value, T >::type conjugate(const T &t)
Definition: solvers.hh:1039
RestartedGMResSolver(std::shared_ptr< LinearOperator< X, Y > > op, std::shared_ptr< ScalarProduct< X > > sp, std::shared_ptr< Preconditioner< X, X > > prec, const ParameterTree &configuration)
Definition: solvers.hh:867
ReboundAllocatorType< X, field_type > fAlloc
field_type Allocator retrieved from domain type
Definition: solvers.hh:822
int _restart
Definition: solvers.hh:1097
ReboundAllocatorType< X, real_type > rAlloc
real_type Allocator retrieved from domain type
Definition: solvers.hh:824
virtual void apply(X &x, Y &b, double reduction, InverseOperatorResult &res)
Apply inverse operator.
Definition: solvers.hh:907
void generatePlaneRotation(field_type &dx, field_type &dy, real_type &cs, field_type &sn)
Definition: solvers.hh:1050
virtual void apply(X &x, Y &b, InverseOperatorResult &res)
Apply inverse operator.
Definition: solvers.hh:894
typename IterativeSolver< X, X >::template Iteration< unsigned int > Iteration
Definition: solvers.hh:1096
void applyPlaneRotation(field_type &dx, field_type &dy, real_type &cs, field_type &sn)
Definition: solvers.hh:1083
implements the Flexible Generalized Minimal Residual (FGMRes) method (right preconditioned)
Definition: solvers.hh:1116
void apply(X &x, Y &b, double reduction, InverseOperatorResult &res) override
Apply inverse operator.
Definition: solvers.hh:1146
Generalized preconditioned conjugate gradient solver.
Definition: solvers.hh:1284
GeneralizedPCGSolver(std::shared_ptr< LinearOperator< X, X > > op, std::shared_ptr< Preconditioner< X, X > > prec, const ParameterTree &configuration)
Constructor.
Definition: solvers.hh:1338
GeneralizedPCGSolver(std::shared_ptr< 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, int restart=10)
Set up nonlinear preconditioned conjugate gradient solver.
Definition: solvers.hh:1354
GeneralizedPCGSolver(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:1308
GeneralizedPCGSolver(std::shared_ptr< LinearOperator< X, X > > op, std::shared_ptr< ScalarProduct< X > > sp, std::shared_ptr< Preconditioner< X, X > > prec, const ParameterTree &configuration)
Definition: solvers.hh:1343
virtual void apply(X &x, X &b, InverseOperatorResult &res)
Apply inverse operator.
Definition: solvers.hh:1368
GeneralizedPCGSolver(LinearOperator< X, X > &op, 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:1320
Accelerated flexible conjugate gradient method.
Definition: solvers.hh:1478
RestartedFCGSolver(LinearOperator< X, X > &op, 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:1506
RestartedFCGSolver(std::shared_ptr< 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, int mmax=10)
Constructor to initialize a RestartedFCG solver.
Definition: solvers.hh:1516
int _mmax
Definition: solvers.hh:1639
typename IterativeSolver< X, X >::template Iteration< unsigned int > Iteration
Definition: solvers.hh:1646
RestartedFCGSolver(std::shared_ptr< LinearOperator< X, X > > op, std::shared_ptr< ScalarProduct< X > > sp, std::shared_ptr< Preconditioner< X, X > > prec, const ParameterTree &config)
Definition: solvers.hh:1542
RestartedFCGSolver(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:1496
virtual void apply(X &x, X &b, InverseOperatorResult &res)
Apply inverse operator.
Definition: solvers.hh:1561
RestartedFCGSolver(std::shared_ptr< LinearOperator< X, X > > op, std::shared_ptr< Preconditioner< X, X > > prec, const ParameterTree &config)
Constructor.
Definition: solvers.hh:1536
Complete flexible conjugate gradient method.
Definition: solvers.hh:1657
virtual void apply(X &x, X &b, InverseOperatorResult &res) override
Apply inverse operator.
Definition: solvers.hh:1671