47 #ifndef EWOMS_SOLVERS_HH
48 #define EWOMS_SOLVERS_HH
62 #include <dune/common/version.hh>
63 #include <dune/istl/istlexception.hh>
64 #include <dune/istl/operators.hh>
65 #include <dune/istl/scalarproducts.hh>
66 #if DUNE_VERSION_NEWER(DUNE_COMMON, 2,3)
67 #include <dune/istl/preconditioner.hh>
68 #include <dune/istl/solver.hh>
70 #include <dune/istl/preconditioners.hh>
71 #include <dune/istl/solvers.hh>
73 #include <dune/common/array.hh>
74 #include <dune/common/timer.hh>
75 #include <dune/common/ftraits.hh>
77 #include <type_traits>
101 template <
class X,
class Y>
119 {
return *convergenceCriterion_; }
126 {
return *convergenceCriterion_; }
133 { convergenceCriterion_ = convCrit; }
144 virtual void apply(X& x, Y& b, Dune::InverseOperatorResult& res) = 0;
151 std::shared_ptr<ConvergenceCriterion<X> > convergenceCriterion_;
178 #if DUNE_VERSION_NEWER(DUNE_COMMON, 2,4)
179 typedef typename Dune::FieldTraits<field_type>::real_type
real_type;
204 template <
class L,
class P>
206 real_type reduction,
int maxit,
int verbose) :
207 ssp(), _op(op), _prec(prec), _sp(ssp), _maxit(maxit), _verbose(verbose)
209 static_assert(static_cast<int>(L::category) == static_cast<int>(P::category),
210 "L and P have to have the same category!");
211 static_assert(static_cast<int>(L::category) == static_cast<int>(Dune::SolverCategory::sequential),
212 "L has to be sequential!");
214 auto crit = std::make_shared<ResidReductionCriterion<X>>(_sp, reduction);
238 template <
class L,
class S,
class P>
240 real_type reduction,
int maxit,
int verbose) :
241 _op(op), _prec(prec), _sp(sp), _maxit(maxit), _verbose(verbose)
243 static_assert(static_cast<int>(L::category) == static_cast<int>(P::category),
244 "L and P must have the same category!");
245 static_assert(static_cast<int>(L::category) == static_cast<int>(S::category),
246 "L and S must have the same category!");
248 auto crit = std::make_shared<ResidReductionCriterion<X>>(_sp, reduction);
254 virtual void apply(X& x, X& b, Dune::InverseOperatorResult& res)
266 _op.applyscaleadd(-1, x, b);
270 std::cout <<
"=== LoopSolver" << std::endl << std::flush;
277 res.converged =
true;
290 for (; i <= _maxit; i++) {
294 _op.applyscaleadd(-1, v, b);
300 res.converged =
true;
306 i = std::min(_maxit, i);
318 res.conv_rate = std::pow(res.reduction, 1.0/res.iterations);
319 res.elapsed = watch.elapsed();
324 std::cout <<
"=== rate=" << res.conv_rate
325 <<
", T=" << res.elapsed
326 <<
", TIT=" << res.elapsed/i
327 <<
", IT=" << i << std::endl;
332 Dune::SeqScalarProduct<X> ssp;
333 Dune::LinearOperator<X, X> &_op;
334 Dune::Preconditioner<X, X> &_prec;
335 Dune::ScalarProduct<X> &_sp;
336 std::shared_ptr<ConvergenceCriterion> convergenceCriterion_;
355 #if DUNE_VERSION_NEWER(DUNE_COMMON, 2,4)
356 typedef typename Dune::FieldTraits<field_type>::real_type
real_type;
368 template <
class L,
class P>
370 real_type reduction,
int maxit,
int verbose) :
371 ssp(), _op(op), _prec(prec), _sp(ssp), _maxit(maxit), _verbose(verbose)
373 static_assert(static_cast<int>(L::category) == static_cast<int>(P::category),
374 "L and P have to have the same category!");
375 static_assert(static_cast<int>(L::category) == static_cast<int>(Dune::SolverCategory::sequential),
376 "L has to be sequential!");
378 auto crit = std::make_shared<ResidReductionCriterion<X>>(_sp, reduction);
386 template <
class L,
class S,
class P>
388 real_type reduction,
int maxit,
int verbose) :
389 _op(op), _prec(prec), _sp(sp), _maxit(maxit), _verbose(verbose)
391 static_assert(static_cast<int>(L::category) == static_cast<int>(P::category),
392 "L and P have to have the same category!");
393 static_assert(static_cast<int>(L::category) == static_cast<int>(S::category),
394 "L and S have to have the same category!");
396 auto crit = std::make_shared<ResidReductionCriterion<X>>(_sp, reduction);
405 virtual void apply(X& x, X& b, Dune::InverseOperatorResult& res)
410 _op.applyscaleadd(-1, x, b);
417 std::cout <<
"=== GradientSolver" << std::endl << std::flush;
423 res.converged =
true;
433 for (; i <=_maxit; i++)
438 lambda = _sp.dot(p, b)/_sp.dot(q, p);
446 res.converged =
true;
452 i = std::min(_maxit, i);
460 res.conv_rate = std::pow(res.reduction, 1.0/res.iterations);
461 res.elapsed = watch.elapsed();
465 Dune::SeqScalarProduct<X> ssp;
466 Dune::LinearOperator<X, X> &_op;
467 Dune::Preconditioner<X, X> &_prec;
468 Dune::ScalarProduct<X> &_sp;
469 std::shared_ptr<ConvergenceCriterion> convergenceCriterion_;
490 typedef typename Dune::FieldTraits<field_type>::real_type
real_type;
497 template <
class L,
class P>
498 CGSolver(L& op, P& prec, real_type reduction,
int maxit,
int verbose) :
499 ssp(), _op(op), _prec(prec), _sp(ssp), _maxit(maxit), _verbose(verbose)
501 static_assert(static_cast<int>(L::category) == static_cast<int>(P::category),
502 "L and P must have the same category!");
503 static_assert(static_cast<int>(L::category) == static_cast<int>(Dune::SolverCategory::sequential),
504 "L must be sequential!");
506 auto crit = std::make_shared<ResidReductionCriterion<X>>(_sp, reduction);
514 template <
class L,
class S,
class P>
515 CGSolver(L& op, S& sp, P& prec, real_type reduction,
int maxit,
int verbose) :
516 _op(op), _prec(prec), _sp(sp), _maxit(maxit), _verbose(verbose)
518 static_assert(static_cast<int>(L::category) == static_cast<int>(P::category),
519 "L and P must have the same category!");
520 static_assert(static_cast<int>(L::category) == static_cast<int>(S::category),
521 "L and S must have the same category!");
523 auto crit = std::make_shared<ResidReductionCriterion<X>>(_sp, reduction);
532 virtual void apply(X& x, X& b, Dune::InverseOperatorResult& res)
537 _op.applyscaleadd(-1, x, b);
545 std::cout <<
"=== CGSolver" << std::endl << std::flush;
551 res.converged =
true;
560 field_type rho, rholast, lambda, alpha, beta;
565 rholast = _sp.dot(p, b);
569 for (; i <= _maxit; i++) {
572 alpha = _sp.dot(p, q);
573 lambda = rholast/alpha;
582 res.converged =
true;
597 i = std::min(_maxit, i);
605 res.conv_rate = std::pow(res.reduction, 1.0/res.iterations);
606 res.elapsed = watch.elapsed();
610 Dune::SeqScalarProduct<X> ssp;
611 Dune::LinearOperator<X, X> &_op;
612 Dune::Preconditioner<X, X> &_prec;
613 Dune::ScalarProduct<X> &_sp;
614 std::shared_ptr<ConvergenceCriterion> convergenceCriterion_;
634 #if DUNE_VERSION_NEWER(DUNE_COMMON, 2,4)
635 typedef typename Dune::FieldTraits<field_type>::real_type
real_type;
646 template <
class L,
class P>
648 real_type reduction,
int maxit,
int verbose) :
649 ssp(), _op(op), _prec(prec), _sp(ssp), _maxit(maxit), _verbose(verbose)
651 static_assert(static_cast<int>(L::category) == static_cast<int>(P::category),
652 "L and P must be of the same category!");
653 static_assert(static_cast<int>(L::category) == static_cast<int>(Dune::SolverCategory::sequential),
654 "L must be sequential!");
656 auto crit = std::make_shared<ResidReductionCriterion<X>>(_sp, reduction);
664 template <
class L,
class S,
class P>
666 real_type reduction,
int maxit,
int verbose) :
667 _op(op), _prec(prec), _sp(sp), _maxit(maxit), _verbose(verbose)
669 static_assert(static_cast<int>(L::category) == static_cast<int>(P::category),
670 "L and P must have the same category!");
671 static_assert(static_cast<int>(L::category) == static_cast<int>(S::category),
672 "L and S must have the same category!");
674 auto crit = std::make_shared<ResidReductionCriterion<X>>(_sp, reduction);
683 virtual void apply(X& x, X& b, Dune::InverseOperatorResult& res)
685 const real_type EPSILON = 1e-80;
687 field_type rho, rho_new, alpha, beta, h, omega;
707 _op.applyscaleadd(-1, x, r);
721 std::cout <<
"=== BiCGSTABSolver" << std::endl << std::flush;
727 res.converged =
true;
738 for (it = 0.5; it < _maxit; it += .5)
745 rho_new = _sp.dot(rt, r);
748 if (std::abs(rho) <= EPSILON)
749 DUNE_THROW(Dune::ISTLError,
"breakdown in BiCGSTAB - rho "
750 << rho <<
" <= EPSILON " << EPSILON
751 <<
" after " << it <<
" iterations");
752 if (std::abs(omega) <= EPSILON)
753 DUNE_THROW(Dune::ISTLError,
"breakdown in BiCGSTAB - omega "
754 << omega <<
" <= EPSILON " << EPSILON
755 <<
" after " << it <<
" iterations");
760 beta = (rho_new / rho) * (alpha / omega);
776 if (std::abs(h) < EPSILON)
777 DUNE_THROW(Dune::ISTLError,
"h=0 in BiCGSTAB");
795 res.converged =
true;
808 omega = _sp.dot(t, r) / _sp.dot(t, t);
826 res.converged =
true;
832 it = std::min(static_cast<double>(_maxit), it);
838 res.iterations =
static_cast<int>(std::ceil(it));
840 res.conv_rate = std::pow(res.reduction, 1.0/res.iterations);
841 res.elapsed = watch.elapsed();
845 Dune::SeqScalarProduct<X> ssp;
846 Dune::LinearOperator<X, X> &_op;
847 Dune::Preconditioner<X, X> &_prec;
848 Dune::ScalarProduct<X> &_sp;
849 std::shared_ptr<ConvergenceCriterion> convergenceCriterion_;
872 #if DUNE_VERSION_NEWER(DUNE_COMMON, 2,4)
873 typedef typename Dune::FieldTraits<field_type>::real_type
real_type;
884 template <
class L,
class P>
885 MINRESSolver(L& op, P& prec, real_type reduction,
int maxit,
int verbose) :
886 ssp(), _op(op), _prec(prec), _sp(ssp), _maxit(maxit), _verbose(verbose)
888 static_assert(static_cast<int>(L::category) == static_cast<int>(P::category),
889 "L and P must have the same category!");
890 static_assert(static_cast<int>(L::category) == static_cast<int>(Dune::SolverCategory::sequential),
891 "L must be sequential!");
893 auto crit = std::make_shared<ResidReductionCriterion<X>>(_sp, reduction);
901 template <
class L,
class S,
class P>
902 MINRESSolver(L& op, S& sp, P& prec, real_type reduction,
int maxit,
int verbose) :
903 _op(op), _prec(prec), _sp(sp), _maxit(maxit), _verbose(verbose)
905 static_assert(static_cast<int>(L::category) == static_cast<int>(P::category),
906 "L and P must have the same category!");
907 static_assert(static_cast<int>(L::category) == static_cast<int>(S::category),
908 "L and S must have the same category!");
910 auto crit = std::make_shared<ResidReductionCriterion<X>>(_sp, reduction);
919 virtual void apply(X& x, X& b, Dune::InverseOperatorResult& res)
929 _op.applyscaleadd(-1, x, b);
934 std::cout <<
"=== MINRESSolver" << std::endl << std::flush;
942 res.converged =
true;
951 field_type alpha, beta;
953 std::array<real_type, 2> c{{0.0, 0.0}};
955 std::array<field_type, 2> s{{0.0, 0.0}};
958 std::array<field_type, 3> T{{0.0, 0.0, 0.0}};
961 std::array<field_type, 2> xi{{1.0, 0.0}};
972 beta = std::sqrt(_sp.dot(z, b));
973 field_type beta0 = beta;
976 std::array<X, 3> p{{b, b, b}};
982 std::array<X, 3> q{{b, b, b}};
991 for (; i <= _maxit; i++) {
1000 q[i2].axpy(-beta, q[i0]);
1004 alpha = _sp.dot(z, q[i2]);
1005 q[i2].axpy(-alpha, q[i1]);
1008 _prec.apply(z, q[i2]);
1012 beta = std::sqrt(_sp.dot(q[i2], z));
1025 T[2] = c[(i+1)%2]*alpha - s[(i+1)%2]*T[1];
1026 T[1] = c[(i+1)%2]*T[1] + s[(i+1)%2]*alpha;
1032 generateGivensRotation(T[2], beta, c[i%2], s[i%2]);
1034 T[2] = c[i%2]*T[2] + s[i%2]*beta;
1036 xi[i%2] = -s[i%2]*xi[(i+1)%2];
1037 xi[(i+1)%2] *= c[i%2];
1041 p[i2].axpy(-T[1], p[i1]);
1042 p[i2].axpy(-T[0], p[i0]);
1046 x.axpy(beta0*xi[(i+1)%2], p[i2]);
1056 res.converged =
true;
1069 res.conv_rate = std::pow(res.reduction, 1.0/res.iterations);
1070 res.elapsed = watch.elapsed();
1075 void generateGivensRotation(field_type& dx, field_type& dy, real_type& cs, field_type& sn)
1077 real_type norm_dx = std::abs(dx);
1078 real_type norm_dy = std::abs(dy);
1079 if (norm_dy < 1e-15) {
1082 }
else if (norm_dx < 1e-15) {
1085 }
else if (norm_dy > norm_dx) {
1086 real_type temp = norm_dx/norm_dy;
1087 cs = 1.0/std::sqrt(1.0 + temp*temp);
1095 real_type temp = norm_dy/norm_dx;
1096 cs = 1.0/std::sqrt(1.0 + temp*temp);
1104 Dune::SeqScalarProduct<X> ssp;
1105 Dune::LinearOperator<X, X>& _op;
1106 Dune::Preconditioner<X, X>& _prec;
1107 Dune::ScalarProduct<X>& _sp;
1108 std::shared_ptr<ConvergenceCriterion> convergenceCriterion_;
1126 template <
class X,
class Y=X,
class F = Y>
1138 #if DUNE_VERSION_NEWER(DUNE_COMMON, 2,4)
1139 typedef typename Dune::FieldTraits<field_type>::real_type
real_type;
1153 template <
class L,
class P>
1156 ssp(), _sp(ssp), _restart(restart),
1157 _maxit(maxit), _verbose(verbose)
1159 static_assert(static_cast<int>(P::category) == static_cast<int>(L::category),
1160 "P and L must be the same category!");
1161 static_assert(static_cast<int>(L::category) == static_cast<int>(Dune::SolverCategory::sequential),
1162 "L must be sequential!");
1164 auto crit = std::make_shared<ResidReductionCriterion<X>>(_sp, reduction);
1174 template <
class L,
class S,
class P>
1177 _sp(sp), _restart(restart),
1178 _maxit(maxit), _verbose(verbose)
1180 static_assert(static_cast<int>(P::category) == static_cast<int>(L::category),
1181 "P and L must have the same category!");
1182 static_assert(static_cast<int>(P::category) == static_cast<int>(S::category),
1183 "P and S must have the same category!");
1185 auto crit = std::make_shared<ResidReductionCriterion<X>>(_sp, reduction);
1194 virtual void apply(X& x, X& b, Dune::InverseOperatorResult& res)
1196 const real_type EPSILON = 1e-80;
1197 const int m = _restart;
1198 real_type norm, norm_0;
1200 std::vector<field_type> s(m+1), sn(m);
1201 std::vector<real_type> cs(m);
1206 std::vector< std::vector<field_type> > H(m+1, s);
1207 std::vector<F> v(m+1, b);
1218 _A.applyscaleadd(-1.0, x, b);
1220 v[0] = 0.0; _W.apply(v[0], b);
1221 norm_0 = _sp.norm(v[0]);
1226 std::cout <<
"=== RestartedGMResSolver" << std::endl << std::flush;
1231 res.converged =
true;
1233 while (j <= _maxit && res.converged !=
true) {
1238 for (i=1; i < m+1; i++)
1241 for (i=0; i < m && j <= _maxit && res.converged !=
true; i++, j++) {
1246 _A.apply(v[i], v[i+1]);
1247 _W.apply(w, v[i+1]);
1248 for (
int k=0; k < i+1; k++) {
1253 H[k][i] = _sp.dot(v[k], w);
1255 w.axpy(-H[k][i], v[k]);
1257 H[i+1][i] = _sp.norm(w);
1258 if (std::abs(H[i+1][i]) < EPSILON)
1259 DUNE_THROW(Dune::ISTLError,
1260 "breakdown in GMRes - |w| == 0.0 after " << j <<
" iterations");
1263 v[i+1] = w; v[i+1] *= 1.0/H[i+1][i];
1266 for (
int k=0; k < i; k++)
1267 applyPlaneRotation(H[k][i], H[k+1][i], cs[k], sn[k]);
1270 generatePlaneRotation(H[i][i], H[i+1][i], cs[i], sn[i]);
1272 applyPlaneRotation(H[i][i], H[i+1][i], cs[i], sn[i]);
1273 applyPlaneRotation(s[i], s[i+1], cs[i], sn[i]);
1276 norm = std::abs(s[i+1]);
1282 res.converged =
true;
1288 update(w, i, H, s, v);
1295 if (res.converged !=
true && j <= _maxit ) {
1298 std::cout <<
"=== GMRes::restart" << std::endl;
1302 _A.applyscaleadd(-1.0, x, b);
1306 norm = _sp.norm(v[0]);
1315 res.iterations = j-1;
1317 res.conv_rate = std::pow(res.reduction, 1.0/res.iterations);
1318 res.elapsed = watch.elapsed();
1325 void update(X& w,
int i,
1326 const std::vector<std::vector<field_type> >& H,
1327 const std::vector<field_type>& s,
1328 const std::vector<X>& v) {
1330 std::vector<field_type> y(s);
1333 for (
int a=i-1; a >=0; a--) {
1334 field_type rhs(s[a]);
1335 for (
int b=a+1; b < i; b++)
1336 rhs -= H[a][b]*y[b];
1345 template <
typename T>
1346 typename std::enable_if<std::is_same<field_type, real_type>::value, T>::type conjugate(
const T& t) {
1350 template <
typename T>
1351 typename std::enable_if<!std::is_same<field_type, real_type>::value, T>::type conjugate(
const T& t) {
1356 generatePlaneRotation(field_type& dx, field_type& dy, real_type& cs, field_type& sn)
1358 real_type norm_dx = std::abs(dx);
1359 real_type norm_dy = std::abs(dy);
1360 if (norm_dy < 1e-15) {
1363 }
else if (norm_dx < 1e-15) {
1366 }
else if (norm_dy > norm_dx) {
1367 real_type temp = norm_dx/norm_dy;
1368 cs = 1.0/std::sqrt(1.0 + temp*temp);
1372 sn *= conjugate(dy)/norm_dy;
1374 real_type temp = norm_dy/norm_dx;
1375 cs = 1.0/std::sqrt(1.0 + temp*temp);
1377 sn *= conjugate(dy/dx);
1383 applyPlaneRotation(field_type& dx, field_type& dy, real_type& cs, field_type& sn)
1385 field_type temp = cs * dx + sn * dy;
1386 dy = -conjugate(sn) * dx + cs * dy;
1390 Dune::LinearOperator<X, X> &_A;
1391 Dune::Preconditioner<X, X> &_W;
1392 Dune::SeqScalarProduct<X> ssp;
1393 Dune::ScalarProduct<X> &_sp;
1425 #if DUNE_VERSION_NEWER(DUNE_COMMON, 2,4)
1426 typedef typename Dune::FieldTraits<field_type>::real_type
real_type;
1439 template <
class L,
class P>
1442 ssp(), _op(op), _prec(prec), _sp(ssp), _maxit(maxit),
1443 _verbose(verbose), _restart(
std::min(maxit, restart))
1445 static_assert(static_cast<int>(L::category) == static_cast<int>(P::category),
1446 "L and P have to have the same category!");
1447 static_assert(static_cast<int>(L::category) ==
1448 static_cast<int>(Dune::SolverCategory::sequential),
1449 "L has to be sequential!");
1451 auto crit = std::make_shared<ResidReductionCriterion<X>>(_sp, reduction);
1462 template <
class L,
class P,
class S>
1464 real_type reduction,
int maxit,
int verbose,
int restart=10) :
1465 _op(op), _prec(prec), _sp(sp), _maxit(maxit), _verbose(verbose),
1466 _restart(
std::min(maxit, restart))
1468 static_assert(static_cast<int>(L::category) == static_cast<int>(P::category),
1469 "L and P must have the same category!");
1470 static_assert(static_cast<int>(L::category) == static_cast<int>(S::category),
1471 "L and S must have the same category!");
1473 auto crit = std::make_shared<ResidReductionCriterion<X>>(_sp, reduction);
1481 virtual void apply(X& x, X& b, Dune::InverseOperatorResult& res)
1486 _op.applyscaleadd(-1, x, b);
1488 std::vector<std::shared_ptr<X> > p(_restart);
1489 std::vector<typename X::field_type> pp(_restart);
1493 p[0].reset(
new X(x));
1497 std::cout <<
"=== GeneralizedPCGSolver" << std::endl << std::flush;
1502 res.converged =
true;
1511 field_type rho, lambda;
1517 _prec.apply(*(p[0]), b);
1518 rho = _sp.dot(*(p[0]), b);
1519 _op.apply(*(p[0]), q);
1520 pp[0] = _sp.dot(*(p[0]), q);
1522 x.axpy(lambda, *(p[0]));
1531 res.converged =
true;
1534 while (i < _maxit) {
1536 int end = std::min(_restart, _maxit-i+1);
1537 for (ii=1; ii < end;++ii )
1542 _prec.apply(prec_res, b);
1544 p[ii].reset(
new X(prec_res));
1545 _op.apply(prec_res, q);
1547 for (
int j=0; j < ii;++j) {
1548 rho = _sp.dot(q, *(p[j]))/pp[j];
1549 p[ii]->axpy(-rho, *(p[j]));
1553 _op.apply(*(p[ii]), q);
1554 pp[ii] = _sp.dot(*(p[ii]), q);
1555 rho = _sp.dot(*(p[ii]), b);
1556 lambda = rho/pp[ii];
1557 x.axpy(lambda, *(p[ii]));
1565 res.converged =
true;
1571 if (end == _restart) {
1572 *(p[0])=*(p[_restart-1]);
1573 pp[0]=pp[_restart-1];
1583 res.conv_rate = std::pow(res.reduction, 1.0/res.iterations);
1584 res.elapsed = watch.elapsed();
1591 Dune::SeqScalarProduct<X> ssp;
1592 Dune::LinearOperator<X, X> &_op;
1593 Dune::Preconditioner<X, X> &_prec;
1594 Dune::ScalarProduct<X> &_sp;
X domain_type
Type of the domain of the operator to be inverted.
Definition: solvers.hh:106
X::field_type field_type
The field type of the operator to be inverted.
Definition: solvers.hh:1424
conjugate gradient method
Definition: solvers.hh:478
virtual ~InverseOperator()
Destructor.
Definition: solvers.hh:147
field_type real_type
Definition: solvers.hh:1429
X::field_type field_type
The field type of the operator to be inverted.
Definition: solvers.hh:1137
F basis_type
The field type of the basis vectors.
Definition: solvers.hh:1145
Abstract base class for all solvers.
Definition: solvers.hh:102
RestartedGMResSolver(L &op, S &sp, P &prec, real_type reduction, int restart, int maxit, int verbose)
Set up solver.
Definition: solvers.hh:1175
Base class for all convergence criteria which only defines an virtual API.
Definition: convergencecriterion.hh:51
LoopSolver(L &op, P &prec, real_type reduction, int maxit, int verbose)
Set up Loop solver.
Definition: solvers.hh:205
Provides a convergence criterion which looks at the reduction of the two-norm of the residual for the...
X range_type
The range type of the operator to be inverted.
Definition: solvers.hh:869
X::field_type field_type
The field type of the operator to be inverted.
Definition: solvers.hh:633
virtual void apply(X &x, X &b, Dune::InverseOperatorResult &res)
Apply inverse operator.
Definition: solvers.hh:1194
X::field_type field_type
The field type of the operator.
Definition: solvers.hh:112
X domain_type
The domain type of the operator to be inverted.
Definition: solvers.hh:1133
X range_type
The range type of the operator that we do the inverse for.
Definition: solvers.hh:175
X domain_type
The domain type of the operator that we do the inverse for.
Definition: solvers.hh:173
Provides a convergence criterion for the linear solvers which looks at the weighted maximum of the di...
X domain_type
The domain type of the operator to be inverted.
Definition: solvers.hh:484
X::field_type field_type
The field type of the operator that we do the inverse for.
Definition: solvers.hh:177
gradient method
Definition: solvers.hh:344
LoopSolver(L &op, S &sp, P &prec, real_type reduction, int maxit, int verbose)
Set up loop solver.
Definition: solvers.hh:239
field_type real_type
Definition: solvers.hh:359
BiCGSTABSolver(L &op, S &sp, P &prec, real_type reduction, int maxit, int verbose)
Set up solver.
Definition: solvers.hh:665
X range_type
The range type of the operator to be inverted.
Definition: solvers.hh:631
virtual void apply(X &x, X &b, Dune::InverseOperatorResult &res)
Apply inverse operator.
Definition: solvers.hh:405
virtual const Ewoms::ConvergenceCriterion< X > & convergenceCriterion() const
Return the criterion to be used to check for convergence of the linear solver.
Definition: solvers.hh:118
field_type real_type
Definition: solvers.hh:1142
CGSolver(L &op, S &sp, P &prec, real_type reduction, int maxit, int verbose)
Set up conjugate gradient solver.
Definition: solvers.hh:515
virtual void apply(X &x, X &b, Dune::InverseOperatorResult &res)
Apply inverse operator.
Definition: solvers.hh:919
X domain_type
The domain type of the operator to be inverted.
Definition: solvers.hh:629
RestartedGMResSolver(L &op, P &prec, real_type reduction, int restart, int maxit, int verbose)
Set up solver.
Definition: solvers.hh:1154
virtual void setConvergenceCriterion(std::shared_ptr< Ewoms::ConvergenceCriterion< X > > convCrit)
Set the criterion to be used to check for convergence of the linear solver.
Definition: solvers.hh:132
X range_type
The range type of the operator to be inverted.
Definition: solvers.hh:1422
virtual void apply(X &x, X &b, Dune::InverseOperatorResult &res)
Definition: solvers.hh:254
Definition: baseauxiliarymodule.hh:35
X::field_type field_type
The field type of the operator to be inverted.
Definition: solvers.hh:871
Preconditioned loop solver.
Definition: solvers.hh:167
GeneralizedPCGSolver(L &op, P &prec, real_type reduction, int maxit, int verbose, int restart=10)
Set up nonlinear preconditioned conjugate gradient solver.
Definition: solvers.hh:1440
MINRESSolver(L &op, P &prec, real_type reduction, int maxit, int verbose)
Set up MINRES solver.
Definition: solvers.hh:885
field_type real_type
Definition: solvers.hh:638
Dune::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: solvers.hh:490
virtual void apply(X &x, X &b, Dune::InverseOperatorResult &res)
Apply inverse operator.
Definition: solvers.hh:532
virtual Ewoms::ConvergenceCriterion< X > & convergenceCriterion()
Return the criterion to be used to check for convergence of the linear solver.
Definition: solvers.hh:125
MINRESSolver(L &op, S &sp, P &prec, real_type reduction, int maxit, int verbose)
Set up MINRES solver.
Definition: solvers.hh:902
Bi-conjugate Gradient Stabilized (BiCG-STAB)
Definition: solvers.hh:623
GradientSolver(L &op, P &prec, real_type reduction, int maxit, int verbose)
Set up solver.
Definition: solvers.hh:369
X range_type
The range type of the operator to be inverted.
Definition: solvers.hh:486
field_type real_type
Definition: solvers.hh:182
BiCGSTABSolver(L &op, P &prec, real_type reduction, int maxit, int verbose)
Set up solver.
Definition: solvers.hh:647
Minimal Residual Method (MINRES)
Definition: solvers.hh:861
Convergence criterion which looks at the weighted absolute value of the residual. ...
X::field_type field_type
The field type of the operator that we do the inverse for.
Definition: solvers.hh:354
Generalized preconditioned conjugate gradient solver.
Definition: solvers.hh:1414
virtual void apply(X &x, X &b, Dune::InverseOperatorResult &res)
Apply inverse operator.
Definition: solvers.hh:1481
X range_type
The range type of the operator that we do the inverse for.
Definition: solvers.hh:352
GeneralizedPCGSolver(L &op, S &sp, P &prec, real_type reduction, int maxit, int verbose, int restart=10)
Set up nonlinear preconditioned conjugate gradient solver.
Definition: solvers.hh:1463
GradientSolver(L &op, S &sp, P &prec, real_type reduction, int maxit, int verbose)
Set up solver.
Definition: solvers.hh:387
implements the Generalized Minimal Residual (GMRes) method
Definition: solvers.hh:1127
Y range_type
Type of the range of the operator to be inverted.
Definition: solvers.hh:109
X domain_type
The domain type of the operator to be inverted.
Definition: solvers.hh:1420
virtual void apply(X &x, Y &b, Dune::InverseOperatorResult &res)=0
Apply inverse operator,.
Y range_type
The range type of the operator to be inverted.
Definition: solvers.hh:1135
X domain_type
The domain type of the operator to be inverted.
Definition: solvers.hh:867
X domain_type
The domain type of the operator that we do the inverse for.
Definition: solvers.hh:350
virtual void apply(X &x, X &b, Dune::InverseOperatorResult &res)
Apply inverse operator.
Definition: solvers.hh:683
CGSolver(L &op, P &prec, real_type reduction, int maxit, int verbose)
Set up conjugate gradient solver.
Definition: solvers.hh:498
field_type real_type
Definition: solvers.hh:876
X::field_type field_type
The field type of the operator to be inverted.
Definition: solvers.hh:488