4 #ifndef DUNE_ISTL_SOLVERS_HH
5 #define DUNE_ISTL_SOLVERS_HH
12 #include <type_traits>
15 #include <dune/common/exceptions.hh>
16 #include <dune/common/hybridutilities.hh>
17 #include <dune/common/math.hh>
18 #include <dune/common/simd/io.hh>
19 #include <dune/common/simd/simd.hh>
20 #include <dune/common/std/type_traits.hh>
21 #include <dune/common/timer.hh>
78 _op->applyscaleadd(-1,x,b);
82 if(iteration.step(0, def)){
98 _op->applyscaleadd(-1,v,b);
100 if(iteration.step(i, def))
146 _op->applyscaleadd(-1,x,b);
149 if(iteration.step(0, def)){
164 lambda =
_sp->dot(p,b)/
_sp->dot(q,p);
169 if(iteration.step(i, def))
204 using enableConditionEstimate_t = Dune::Std::bool_constant<(std::is_same<field_type,float>::value || std::is_same<field_type,double>::value)>;
220 condition_estimate_(condition_estimate)
223 condition_estimate_ =
false;
224 std::cerr <<
"WARNING: Condition estimate was disabled. It is only available for double and float field types!" << std::endl;
236 scalar_real_type reduction,
int maxit,
int verbose,
bool condition_estimate) :
IterativeSolver<X,X>(op, sp, prec, reduction, maxit, verbose),
237 condition_estimate_(condition_estimate)
239 if (condition_estimate && !(std::is_same<field_type,float>::value || std::is_same<field_type,double>::value)) {
240 condition_estimate_ =
false;
241 std::cerr <<
"WARNING: Condition estimate was disabled. It is only available for double and float field types!" << std::endl;
254 scalar_real_type reduction,
int maxit,
int verbose,
bool condition_estimate)
256 condition_estimate_(condition_estimate)
258 if (condition_estimate && !(std::is_same<field_type,float>::value || std::is_same<field_type,double>::value)) {
259 condition_estimate_ =
false;
260 std::cerr <<
"WARNING: Condition estimate was disabled. It is only available for double and float field types!" << std::endl;
280 _op->applyscaleadd(-1,x,b);
283 if(iteration.step(0, def)){
292 std::vector<real_type> lambdas(0);
293 std::vector<real_type> betas(0);
301 rholast =
_sp->dot(p,b);
309 alpha =
_sp->dot(p,q);
310 lambda = rholast/alpha;
312 if (condition_estimate_)
313 lambdas.push_back(std::real(
id(lambda)));
320 if(iteration.step(i, def))
329 if (condition_estimate_)
330 betas.push_back(std::real(
id(beta)));
339 if (condition_estimate_) {
352 row.insert(row.index()-1);
353 row.insert(row.index());
354 if (row.index() < T.
N() - 1)
355 row.insert(row.index()+1);
357 for (
int row = 0; row < i; ++row) {
359 T[row][row-1] = sqrt(betas[row-1]) / lambdas[row-1];
362 T[row][row] = 1.0 / id(lambdas[row]);
364 T[row][row] += betas[row-1] / lambdas[row-1];
368 T[row][row+1] = sqrt(betas[row]) / lambdas[row];
384 std::cout <<
"Min eigv estimate: " << Simd::io(min_eigv) <<
'\n';
385 std::cout <<
"Max eigv estimate: " << Simd::io(max_eigv) <<
'\n';
386 std::cout <<
"Condition estimate: "
387 << Simd::io(max_eigv / min_eigv) << std::endl;
391 std::cerr <<
"WARNING: Condition estimate was requested. This requires ARPACK, but ARPACK was not found!" << std::endl;
397 bool condition_estimate_ =
false;
440 const Simd::Scalar<real_type> EPSILON=1e-80;
443 field_type rho, rho_new, alpha, beta, h, omega;
464 _op->applyscaleadd(-1,x,r);
469 if(iteration.step(0, norm)){
484 for (it = 0.5; it <
_maxit; it+=.5)
491 rho_new =
_sp->dot(rt,r);
494 if (Simd::allTrue(abs(rho) <= EPSILON))
495 DUNE_THROW(
SolverAbort,
"breakdown in BiCGSTAB - rho "
496 << Simd::io(rho) <<
" <= EPSILON " << EPSILON
497 <<
" after " << it <<
" iterations");
498 if (Simd::allTrue(abs(omega) <= EPSILON))
499 DUNE_THROW(
SolverAbort,
"breakdown in BiCGSTAB - omega "
500 << Simd::io(omega) <<
" <= EPSILON " << EPSILON
501 <<
" after " << it <<
" iterations");
508 beta = ( rho_new / rho ) * ( alpha / omega );
524 if ( Simd::allTrue(abs(h) < EPSILON) )
525 DUNE_THROW(
SolverAbort,
"abs(h) < EPSILON in BiCGSTAB - abs(h) "
526 << Simd::io(abs(h)) <<
" < EPSILON " << EPSILON
527 <<
" after " << it <<
" iterations");
543 if(iteration.step(it, norm)){
557 omega =
_sp->dot(t,r)/
_sp->dot(t,t);
573 if(iteration.step(it, norm)){
588 template<
class CountType>
627 _op->applyscaleadd(-1,x,b);
631 if(iteration.step(0, def)){
639 std::array<real_type,2> c{{0.0,0.0}};
641 std::array<field_type,2> s{{0.0,0.0}};
644 std::array<field_type,3> T{{0.0,0.0,0.0}};
647 std::array<field_type,2> xi{{1.0,0.0}};
658 beta = sqrt(
_sp->dot(b,z));
662 std::array<X,3> p{{b,b,b}};
668 std::array<X,3> q{{b,b,b}};
686 q[i2].axpy(-beta,q[i0]);
690 alpha =
_sp->dot(z,q[i2]);
691 q[i2].axpy(-alpha,q[i1]);
694 _prec->apply(z,q[i2]);
698 beta = sqrt(
_sp->dot(q[i2],z));
711 T[2] = c[(i+1)%2]*alpha - s[(i+1)%2]*T[1];
712 T[1] = c[(i+1)%2]*T[1] + s[(i+1)%2]*alpha;
718 generateGivensRotation(T[2],beta,c[i%2],s[i%2]);
720 T[2] = c[i%2]*T[2] + s[i%2]*beta;
722 xi[i%2] = -s[i%2]*xi[(i+1)%2];
723 xi[(i+1)%2] *= c[i%2];
727 p[i2].axpy(-T[1],p[i1]);
728 p[i2].axpy(-T[0],p[i0]);
732 x.axpy(beta0*xi[(i+1)%2],p[i2]);
739 def = abs(beta0*xi[i%2]);
740 if(iteration.step(i, def)){
760 real_type norm_max = max(norm_dx, norm_dy);
761 real_type norm_min = min(norm_dx, norm_dy);
764 cs = Simd::cond(norm_dy < eps,
766 Simd::cond(norm_dx < eps,
768 Simd::cond(norm_dy > norm_dx,
772 sn = Simd::cond(norm_dy < eps,
774 Simd::cond(norm_dx < eps,
776 Simd::cond(norm_dy > norm_dx,
810 template<
class X,
class Y=X,
class F = Y>
911 const Simd::Scalar<real_type> EPSILON = 1e-80;
915 std::vector<field_type,fAlloc> s(m+1), sn(m);
916 std::vector<real_type,rAlloc> cs(m);
921 std::vector< std::vector<field_type,fAlloc> > H(m+1,s);
922 std::vector<F> v(m+1,b);
930 _op->applyscaleadd(-1.0,x,b);
932 v[0] = 0.0;
_prec->apply(v[0],b);
933 norm =
_sp->norm(v[0]);
934 if(iteration.step(0, norm)){
952 _op->apply(v[i],v[i+1]);
953 _prec->apply(w,v[i+1]);
954 for(
int k=0; k<i+1; k++) {
959 H[k][i] =
_sp->dot(v[k],w);
961 w.axpy(-H[k][i],v[k]);
963 H[i+1][i] =
_sp->norm(w);
964 if(Simd::allTrue(abs(H[i+1][i]) < EPSILON))
966 "breakdown in GMRes - |w| == 0.0 after " << j <<
" iterations");
969 v[i+1] = w; v[i+1] *=
real_type(1.0)/H[i+1][i];
972 for(
int k=0; k<i; k++)
984 iteration.step(j, norm);
1000 std::cout <<
"=== GMRes::restart" << std::endl;
1004 _op->applyscaleadd(-1.0,x,b);
1007 _prec->apply(v[0],b);
1008 norm =
_sp->norm(v[0]);
1020 const std::vector<std::vector<field_type,fAlloc> >& H,
1021 const std::vector<field_type,fAlloc>& s,
1022 const std::vector<X>& v) {
1024 std::vector<field_type,fAlloc> y(s);
1027 for(
int a=i-1; a>=0; a--) {
1029 for(
int b=a+1; b<i; b++)
1030 rhs -= H[a][b]*y[b];
1039 template<
typename T>
1040 typename std::enable_if<std::is_same<field_type,real_type>::value,T>::type
conjugate(
const T& t) {
1044 template<
typename T>
1045 typename std::enable_if<!std::is_same<field_type,real_type>::value,T>::type
conjugate(
const T& t) {
1060 real_type norm_max = max(norm_dx, norm_dy);
1061 real_type norm_min = min(norm_dx, norm_dy);
1064 cs = Simd::cond(norm_dy < eps,
1066 Simd::cond(norm_dx < eps,
1068 Simd::cond(norm_dy > norm_dx,
1072 sn = Simd::cond(norm_dy < eps,
1074 Simd::cond(norm_dx < eps,
1076 Simd::cond(norm_dy > norm_dx,
1115 template<
class X,
class Y=X,
class F = Y>
1150 const Simd::Scalar<real_type> EPSILON = 1e-80;
1154 std::vector<field_type,fAlloc> s(m+1), sn(m);
1155 std::vector<real_type,rAlloc> cs(m);
1158 std::vector< std::vector<field_type,fAlloc> > H(m+1,s);
1159 std::vector<F> v(m+1,b);
1160 std::vector<X> w(m+1,b);
1162 Iteration iteration(*
this,res);
1168 _op->applyscaleadd(-1.0, x, v[0]);
1170 norm =
_sp->norm(v[0]);
1171 if(iteration.step(0, norm)){
1180 v[0] *= (1.0 / norm);
1182 for(i=1; i<m+1; ++i)
1190 _prec->apply(w[i], v[i]);
1193 _op->apply(w[i], v[i+1]);
1195 for(
int k=0; k<i+1; k++)
1201 H[k][i] =
_sp->dot(v[k],v[i+1]);
1203 v[i+1].axpy(-H[k][i], v[k]);
1205 H[i+1][i] =
_sp->norm(v[i+1]);
1206 if(Simd::allTrue(abs(H[i+1][i]) < EPSILON))
1207 DUNE_THROW(
SolverAbort,
"breakdown in fGMRes - |w| (-> "
1208 << w[i] <<
") == 0.0 after "
1209 << j <<
" iterations");
1228 iteration.step(j, norm);
1233 this->
update(tmp, i, H, s, w);
1243 std::cout <<
"=== fGMRes::restart" << std::endl;
1247 _op->applyscaleadd(-1.0, x,v[0]);
1249 norm =
_sp->norm(v[0]);
1341 _restart(configuration.
get<int>(
"restart"))
1346 _restart(configuration.
get<int>(
"restart"))
1371 Iteration iteration(*
this, res);
1373 _op->applyscaleadd(-1,x,b);
1375 std::vector<std::shared_ptr<X> > p(_restart);
1376 std::vector<field_type,fAlloc> pp(_restart);
1380 p[0].reset(
new X(x));
1383 if(iteration.step(0, def)){
1394 _prec->apply(*(p[0]),b);
1395 rho =
_sp->dot(*(p[0]),b);
1396 _op->apply(*(p[0]),q);
1397 pp[0] =
_sp->dot(*(p[0]),q);
1399 x.axpy(lambda,*(p[0]));
1405 if(iteration.step(i, def)){
1412 int end=std::min(_restart,
_maxit-i+1);
1413 for (ii=1; ii<end; ++ii )
1418 _prec->apply(prec_res,b);
1420 p[ii].reset(
new X(prec_res));
1421 _op->apply(prec_res, q);
1423 for(
int j=0; j<ii; ++j) {
1424 rho =
_sp->dot(q,*(p[j]))/pp[j];
1425 p[ii]->axpy(-rho, *(p[j]));
1429 _op->apply(*(p[ii]),q);
1430 pp[ii] =
_sp->dot(*(p[ii]),q);
1431 rho =
_sp->dot(*(p[ii]),b);
1432 lambda = rho/pp[ii];
1433 x.axpy(lambda,*(p[ii]));
1440 iteration.step(i, def);
1445 *(p[0])=*(p[_restart-1]);
1446 pp[0]=pp[_restart-1];
1508 scalar_real_type reduction,
int maxit,
int verbose,
int mmax = 10) :
IterativeSolver<X,X>(op, sp, prec, reduction, maxit, verbose),
_mmax(mmax)
1540 const ParameterTree& config)
1562 _op->applyscaleadd(-1,x,b);
1565 std::vector<X> d(
_mmax+1, x);
1566 std::vector<X> Ad(
_mmax+1, x);
1567 std::vector<field_type,rAlloc> ddotAd(
_mmax+1,0);
1571 if(iteration.step(0, def)){
1583 for (; i_bounded <=
_mmax && i<=
_maxit; i_bounded++) {
1585 _prec->apply(d[i_bounded], b);
1588 orthogonalizations(i_bounded,Ad,w,ddotAd,d);
1591 _op->apply(d[i_bounded], Ad[i_bounded]);
1592 ddotAd[i_bounded]=
_sp->dot(d[i_bounded],Ad[i_bounded]);
1593 alpha =
_sp->dot(d[i_bounded], b)/ddotAd[i_bounded];
1596 x.axpy(alpha, d[i_bounded]);
1597 b.axpy(-alpha, Ad[i_bounded]);
1602 iteration.step(i, def);
1606 cycle(Ad,d,ddotAd,i_bounded);
1619 for (
int k = 0; k < i_bounded; k++) {
1620 d[i_bounded].axpy(-
_sp->dot(Ad[k], w) / ddotAd[k], d[k]);
1625 virtual void cycle(std::vector<X>& Ad,std::vector<X>& d,std::vector<
field_type,ReboundAllocatorType<X,field_type> >& ddotAd,
int& i_bounded) {
1628 std::swap(Ad[0], Ad[
_mmax]);
1629 std::swap(d[0], d[
_mmax]);
1630 std::swap(ddotAd[0], ddotAd[
_mmax]);
1676 for (
int k = 0; k < _k_limit; k++) {
1678 d[i_bounded].axpy(-
_sp->dot(Ad[k], w) / ddotAd[k], d[k]);
1681 if(_k_limit<=i_bounded)
1687 virtual void cycle(std::vector<X>& Ad,std::vector<X>& d,std::vector<
field_type,ReboundAllocatorType<X,field_type> >& ddotAd,
int& i_bounded)
override {
1691 _k_limit = Ad.size();
Implementation of the BCRSMatrix class.
Define general, extensible interface for operators. The available implementation wraps a matrix.
Define base class for scalar product and norm.
Define general, extensible interface for inverse operators.
DUNE_REGISTER_ITERATIVE_SOLVER("loopsolver", defaultIterativeSolverCreator< Dune::LoopSolver >())
Definition: allocator.hh:7
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:292
typename std::allocator_traits< typename AllocatorTraits< T >::type >::template rebind_alloc< X > ReboundAllocatorType
Definition: allocator.hh:33
std::enable_if_t<!is_complex< T >::value, T > conj(const T &r)
Definition: matrixmarket.hh:732
A sparse block matrix with compressed row storage.
Definition: bcrsmatrix.hh:426
@ row_wise
Build in a row-wise manner.
Definition: bcrsmatrix.hh:480
CreateIterator createend()
get create iterator pointing to one after the last block
Definition: bcrsmatrix.hh:1062
CreateIterator createbegin()
get initial create iterator
Definition: bcrsmatrix.hh:1056
size_type N() const
number of rows (counted in blocks)
Definition: bcrsmatrix.hh:1931
A vector of blocks with memory management.
Definition: bvector.hh:403
Wrapper to use a range of ARPACK++ eigenvalue solvers.
Definition: arpackpp.hh:242
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:286
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:388
Thrown when a solver aborts due to some problem.
Definition: istlexception.hh:43
Base class for scalar product and norm computation.
Definition: scalarproducts.hh:48
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:58
virtual void apply(X &x, X &b, InverseOperatorResult &res)
Apply inverse operator,.
Definition: solvers.hh:72
typename IterativeSolver< X, X >::template Iteration< unsigned int > Iteration
Definition: solvers.hh:115
gradient method
Definition: solvers.hh:123
virtual void apply(X &x, X &b, InverseOperatorResult &res)
Apply inverse operator.
Definition: solvers.hh:141
typename IterativeSolver< X, X >::template Iteration< unsigned int > Iteration
Definition: solvers.hh:183
conjugate gradient method
Definition: solvers.hh:189
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:218
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:235
virtual void apply(X &x, X &b, InverseOperatorResult &res)
Apply inverse operator.
Definition: solvers.hh:275
Dune::Std::bool_constant<(std::is_same< field_type, float >::value||std::is_same< field_type, double >::value)> enableConditionEstimate_t
Definition: solvers.hh:204
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:252
typename IterativeSolver< X, X >::template Iteration< unsigned int > Iteration
Definition: solvers.hh:410
Bi-conjugate Gradient Stabilized (BiCG-STAB)
Definition: solvers.hh:417
typename IterativeSolver< X, X >::template Iteration< CountType > Iteration
Definition: solvers.hh:589
virtual void apply(X &x, X &b, InverseOperatorResult &res)
Apply inverse operator.
Definition: solvers.hh:437
Minimal Residual Method (MINRES)
Definition: solvers.hh:600
virtual void apply(X &x, X &b, InverseOperatorResult &res)
Apply inverse operator.
Definition: solvers.hh:618
typename IterativeSolver< X, X >::template Iteration< unsigned int > Iteration
Definition: solvers.hh:793
implements the Generalized Minimal Residual (GMRes) method
Definition: solvers.hh:812
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:835
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:846
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:1019
RestartedGMResSolver(std::shared_ptr< LinearOperator< X, Y > > op, std::shared_ptr< Preconditioner< X, X > > prec, const ParameterTree &configuration)
Constructor.
Definition: solvers.hh:863
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:868
ReboundAllocatorType< X, field_type > fAlloc
field_type Allocator retrieved from domain type
Definition: solvers.hh:823
std::enable_if<!std::is_same< field_type, real_type >::value, T >::type conjugate(const T &t)
Definition: solvers.hh:1045
int _restart
Definition: solvers.hh:1098
ReboundAllocatorType< X, real_type > rAlloc
real_type Allocator retrieved from domain type
Definition: solvers.hh:825
virtual void apply(X &x, Y &b, double reduction, InverseOperatorResult &res)
Apply inverse operator.
Definition: solvers.hh:908
void generatePlaneRotation(field_type &dx, field_type &dy, real_type &cs, field_type &sn)
Definition: solvers.hh:1051
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:879
virtual void apply(X &x, Y &b, InverseOperatorResult &res)
Apply inverse operator.
Definition: solvers.hh:895
std::enable_if< std::is_same< field_type, real_type >::value, T >::type conjugate(const T &t)
Definition: solvers.hh:1040
typename IterativeSolver< X, X >::template Iteration< unsigned int > Iteration
Definition: solvers.hh:1097
void applyPlaneRotation(field_type &dx, field_type &dy, real_type &cs, field_type &sn)
Definition: solvers.hh:1084
implements the Flexible Generalized Minimal Residual (FGMRes) method (right preconditioned)
Definition: solvers.hh:1117
void apply(X &x, Y &b, double reduction, InverseOperatorResult &res) override
Apply inverse operator.
Definition: solvers.hh:1147
Generalized preconditioned conjugate gradient solver.
Definition: solvers.hh:1285
GeneralizedPCGSolver(std::shared_ptr< LinearOperator< X, X > > op, std::shared_ptr< Preconditioner< X, X > > prec, const ParameterTree &configuration)
Constructor.
Definition: solvers.hh:1339
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:1309
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:1344
virtual void apply(X &x, X &b, InverseOperatorResult &res)
Apply inverse operator.
Definition: solvers.hh:1369
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:1321
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:1355
Accelerated flexible conjugate gradient method.
Definition: solvers.hh:1479
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:1507
int _mmax
Definition: solvers.hh:1634
typename IterativeSolver< X, X >::template Iteration< unsigned int > Iteration
Definition: solvers.hh:1641
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:1497
virtual void apply(X &x, X &b, InverseOperatorResult &res)
Apply inverse operator.
Definition: solvers.hh:1556
RestartedFCGSolver(std::shared_ptr< LinearOperator< X, X >> op, std::shared_ptr< ScalarProduct< X >> sp, std::shared_ptr< Preconditioner< X, X >> prec, const ParameterTree &config)
Constructor.
Definition: solvers.hh:1537
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:1517
Complete flexible conjugate gradient method.
Definition: solvers.hh:1652
virtual void apply(X &x, X &b, InverseOperatorResult &res) override
Apply inverse operator.
Definition: solvers.hh:1666