My Project
Loading...
Searching...
No Matches
ISTLSolver.hpp
1/*
2 Copyright 2016 IRIS AS
3 Copyright 2019, 2020 Equinor ASA
4 Copyright 2020 SINTEF Digital, Mathematics and Cybernetics
5
6 This file is part of the Open Porous Media project (OPM).
7
8 OPM is free software: you can redistribute it and/or modify
9 it under the terms of the GNU General Public License as published by
10 the Free Software Foundation, either version 3 of the License, or
11 (at your option) any later version.
12
13 OPM is distributed in the hope that it will be useful,
14 but WITHOUT ANY WARRANTY; without even the implied warranty of
15 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
16 GNU General Public License for more details.
17
18 You should have received a copy of the GNU General Public License
19 along with OPM. If not, see <http://www.gnu.org/licenses/>.
20*/
21
22#ifndef OPM_ISTLSOLVER_HEADER_INCLUDED
23#define OPM_ISTLSOLVER_HEADER_INCLUDED
24
25#include <dune/istl/owneroverlapcopy.hh>
26#include <dune/istl/solver.hh>
27
28#include <opm/common/ErrorMacros.hpp>
29#include <opm/common/Exceptions.hpp>
30#include <opm/common/TimingMacros.hpp>
31
36#include <opm/simulators/flow/BlackoilModelParameters.hpp>
39#include <opm/simulators/linalg/ExtractParallelGridInformationToISTL.hpp>
40#include <opm/simulators/linalg/FlowLinearSolverParameters.hpp>
41#include <opm/simulators/linalg/matrixblock.hh>
43#include <opm/simulators/linalg/PreconditionerWithUpdate.hpp>
44#include <opm/simulators/linalg/WellOperators.hpp>
45#include <opm/simulators/linalg/WriteSystemMatrixHelper.hpp>
46#include <opm/simulators/linalg/findOverlapRowsAndColumns.hpp>
47#include <opm/simulators/linalg/getQuasiImpesWeights.hpp>
48#include <opm/simulators/linalg/setupPropertyTree.hpp>
49
50#include <any>
51#include <cstddef>
52#include <functional>
53#include <memory>
54#include <set>
55#include <sstream>
56#include <string>
57#include <tuple>
58#include <vector>
59
60namespace Opm::Properties {
61
62namespace TTag {
64 using InheritsFrom = std::tuple<FlowIstlSolverParams>;
65};
66}
67
68template <class TypeTag, class MyTypeTag>
69struct WellModel;
70
73template<class TypeTag>
74struct SparseMatrixAdapter<TypeTag, TTag::FlowIstlSolver>
75{
76private:
80
81public:
82 using type = typename Linear::IstlSparseMatrixAdapter<Block>;
83};
84
85} // namespace Opm::Properties
86
87namespace Opm
88{
89
90
91namespace detail
92{
93
94template<class Matrix, class Vector, class Comm>
96{
97 using AbstractSolverType = Dune::InverseOperator<Vector, Vector>;
98 using AbstractOperatorType = Dune::AssembledLinearOperator<Matrix, Vector, Vector>;
100
101 void create(const Matrix& matrix,
102 bool parallel,
103 const PropertyTree& prm,
104 std::size_t pressureIndex,
105 std::function<Vector()> weightCalculator,
106 const bool forceSerial,
107 Comm& comm);
108
109 std::unique_ptr<AbstractSolverType> solver_;
110 std::unique_ptr<AbstractOperatorType> op_;
111 std::unique_ptr<LinearOperatorExtra<Vector,Vector>> wellOperator_;
112 AbstractPreconditionerType* pre_ = nullptr;
113 std::size_t interiorCellNum_ = 0;
114};
115
116
117#ifdef HAVE_MPI
119void copyParValues(std::any& parallelInformation, std::size_t size,
120 Dune::OwnerOverlapCopyCommunication<int,int>& comm);
121#endif
122
125template<class Matrix>
126void makeOverlapRowsInvalid(Matrix& matrix,
127 const std::vector<int>& overlapRows);
128
131template<class Matrix, class Grid>
132std::unique_ptr<Matrix> blockJacobiAdjacency(const Grid& grid,
133 const std::vector<int>& cell_part,
134 std::size_t nonzeroes,
135 const std::vector<std::set<int>>& wellConnectionsGraph);
136}
137
142 template <class TypeTag>
144 {
145 protected:
153 using Matrix = typename SparseMatrixAdapter::IstlMatrix;
156 using AbstractSolverType = Dune::InverseOperator<Vector, Vector>;
157 using AbstractOperatorType = Dune::AssembledLinearOperator<Matrix, Vector, Vector>;
161 constexpr static std::size_t pressureIndex = GetPropType<TypeTag, Properties::Indices>::pressureSwitchIdx;
162
163 enum { enableMICP = getPropValue<TypeTag, Properties::EnableMICP>() };
164 enum { enablePolymerMolarWeight = getPropValue<TypeTag, Properties::EnablePolymerMW>() };
165 constexpr static bool isIncompatibleWithCprw = enableMICP || enablePolymerMolarWeight;
166
167#if HAVE_MPI
168 using CommunicationType = Dune::OwnerOverlapCopyCommunication<int,int>;
169#else
170 using CommunicationType = Dune::Communication<int>;
171#endif
172
173 public:
174 using AssembledLinearOperatorType = Dune::AssembledLinearOperator< Matrix, Vector, Vector >;
175
176 static void registerParameters()
177 {
178 FlowLinearSolverParameters::registerParameters();
179 }
180
188 ISTLSolver(const Simulator& simulator,
189 const FlowLinearSolverParameters& parameters,
190 bool forceSerial = false)
191 : simulator_(simulator),
192 iterations_( 0 ),
193 converged_(false),
194 matrix_(nullptr),
195 parameters_{parameters},
196 forceSerial_(forceSerial)
197 {
198 initialize();
199 }
200
203 explicit ISTLSolver(const Simulator& simulator)
204 : simulator_(simulator),
205 iterations_( 0 ),
206 solveCount_(0),
207 converged_(false),
208 matrix_(nullptr)
209 {
210 parameters_.resize(1);
211 parameters_[0].init(simulator_.vanguard().eclState().getSimulationConfig().useCPR());
212 initialize();
213 }
214
215 void initialize()
216 {
217 OPM_TIMEBLOCK(IstlSolver);
218
219 if (isIncompatibleWithCprw) {
220 // Some model variants are incompatible with the CPRW linear solver.
221 if (parameters_[0].linsolver_ == "cprw" || parameters_[0].linsolver_ == "hybrid") {
222 std::string incompatible_model = "Unknown";
223 if (enableMICP) {
224 incompatible_model = "MICP";
225 } else if (enablePolymerMolarWeight) {
226 incompatible_model = "Polymer injectivity";
227 }
228 OPM_THROW(std::runtime_error,
229 incompatible_model + " model is incompatible with the CPRW linear solver.\n"
230 "Choose a different option, for example --linear-solver=ilu0");
231 }
232 }
233
234 if (parameters_[0].linsolver_ == "hybrid") {
235 // Experimental hybrid configuration.
236 // When chosen, will set up two solvers, one with CPRW
237 // and the other with ILU0 preconditioner. More general
238 // options may be added later.
239 prm_.clear();
240 parameters_.clear();
241 {
242 FlowLinearSolverParameters para;
243 para.init(false);
244 para.linsolver_ = "cprw";
245 parameters_.push_back(para);
246 prm_.push_back(setupPropertyTree(parameters_[0],
247 Parameters::IsSet<Parameters::LinearSolverMaxIter>(),
248 Parameters::IsSet<Parameters::LinearSolverReduction>()));
249 }
250 {
251 FlowLinearSolverParameters para;
252 para.init(false);
253 para.linsolver_ = "ilu0";
254 parameters_.push_back(para);
255 prm_.push_back(setupPropertyTree(parameters_[1],
256 Parameters::IsSet<Parameters::LinearSolverMaxIter>(),
257 Parameters::IsSet<Parameters::LinearSolverReduction>()));
258 }
259 // ------------
260 } else {
261 assert(parameters_.size() == 1);
262 assert(prm_.empty());
263
264 // Do a normal linear solver setup.
265 if (parameters_[0].is_nldd_local_solver_) {
266 prm_.push_back(setupPropertyTree(parameters_[0],
267 Parameters::IsSet<Parameters::NlddLocalLinearSolverMaxIter>(),
268 Parameters::IsSet<Parameters::NlddLocalLinearSolverReduction>()));
269 }
270 else {
271 prm_.push_back(setupPropertyTree(parameters_[0],
272 Parameters::IsSet<Parameters::LinearSolverMaxIter>(),
273 Parameters::IsSet<Parameters::LinearSolverReduction>()));
274 }
275 }
276 flexibleSolver_.resize(prm_.size());
277
278 const bool on_io_rank = (simulator_.gridView().comm().rank() == 0);
279#if HAVE_MPI
280 comm_.reset( new CommunicationType( simulator_.vanguard().grid().comm() ) );
281#endif
282 extractParallelGridInformationToISTL(simulator_.vanguard().grid(), parallelInformation_);
283
284 // For some reason simulator_.model().elementMapper() is not initialized at this stage
285 //const auto& elemMapper = simulator_.model().elementMapper(); //does not work.
286 // Set it up manually
287 ElementMapper elemMapper(simulator_.vanguard().gridView(), Dune::mcmgElementLayout());
288 detail::findOverlapAndInterior(simulator_.vanguard().grid(), elemMapper, overlapRows_, interiorRows_);
289 useWellConn_ = Parameters::Get<Parameters::MatrixAddWellContributions>();
290 const bool ownersFirst = Parameters::Get<Parameters::OwnerCellsFirst>();
291 if (!ownersFirst) {
292 const std::string msg = "The linear solver no longer supports --owner-cells-first=false.";
293 if (on_io_rank) {
294 OpmLog::error(msg);
295 }
296 OPM_THROW_NOLOG(std::runtime_error, msg);
297 }
298
299 const int interiorCellNum_ = detail::numMatrixRowsToUseInSolver(simulator_.vanguard().grid(), true);
300 for (auto& f : flexibleSolver_) {
301 f.interiorCellNum_ = interiorCellNum_;
302 }
303
304#if HAVE_MPI
305 if (isParallel()) {
306 const std::size_t size = simulator_.vanguard().grid().leafGridView().size(0);
307 detail::copyParValues(parallelInformation_, size, *comm_);
308 }
309#endif
310
311 // Print parameters to PRT/DBG logs.
312 if (on_io_rank && parameters_[activeSolverNum_].linear_solver_print_json_definition_) {
313 std::ostringstream os;
314 os << "Property tree for linear solvers:\n";
315 for (std::size_t i = 0; i<prm_.size(); i++) {
316 prm_[i].write_json(os, true);
317 }
318 OpmLog::note(os.str());
319 }
320 }
321
322 // nothing to clean here
323 void eraseMatrix()
324 {
325 }
326
327 void setActiveSolver(const int num)
328 {
329 if (num > static_cast<int>(prm_.size()) - 1) {
330 OPM_THROW(std::logic_error, "Solver number " + std::to_string(num) + " not available.");
331 }
332 activeSolverNum_ = num;
333 if (simulator_.gridView().comm().rank() == 0) {
334 OpmLog::debug("Active solver = " + std::to_string(activeSolverNum_)
335 + " (" + parameters_[activeSolverNum_].linsolver_ + ")");
336 }
337 }
338
339 int numAvailableSolvers()
340 {
341 return flexibleSolver_.size();
342 }
343
344 void initPrepare(const Matrix& M, Vector& b)
345 {
346 const bool firstcall = (matrix_ == nullptr);
347
348 // update matrix entries for solvers.
349 if (firstcall) {
350 // model will not change the matrix object. Hence simply store a pointer
351 // to the original one with a deleter that does nothing.
352 // Outch! We need to be able to scale the linear system! Hence const_cast
353 matrix_ = const_cast<Matrix*>(&M);
354
355 useWellConn_ = Parameters::Get<Parameters::MatrixAddWellContributions>();
356 // setup sparsity pattern for jacobi matrix for preconditioner (only used for openclSolver)
357 } else {
358 // Pointers should not change
359 if ( &M != matrix_ ) {
360 OPM_THROW(std::logic_error,
361 "Matrix objects are expected to be reused when reassembling!");
362 }
363 }
364 rhs_ = &b;
365
366 // TODO: check all solvers, not just one.
367 if (isParallel() && prm_[activeSolverNum_].template get<std::string>("preconditioner.type") != "ParOverILU0") {
368 detail::makeOverlapRowsInvalid(getMatrix(), overlapRows_);
369 }
370 }
371
372 void prepare(const SparseMatrixAdapter& M, Vector& b)
373 {
374 prepare(M.istlMatrix(), b);
375 }
376
377 void prepare(const Matrix& M, Vector& b)
378 {
379 OPM_TIMEBLOCK(istlSolverPrepare);
380
381 initPrepare(M,b);
382
383 prepareFlexibleSolver();
384 }
385
386
387 void setResidual(Vector& /* b */)
388 {
389 // rhs_ = &b; // Must be handled in prepare() instead.
390 }
391
392 void getResidual(Vector& b) const
393 {
394 b = *rhs_;
395 }
396
397 void setMatrix(const SparseMatrixAdapter& /* M */)
398 {
399 // matrix_ = &M.istlMatrix(); // Must be handled in prepare() instead.
400 }
401
402 int getSolveCount() const {
403 return solveCount_;
404 }
405
406 void resetSolveCount() {
407 solveCount_ = 0;
408 }
409
410 bool solve(Vector& x)
411 {
412 OPM_TIMEBLOCK(istlSolverSolve);
413 ++solveCount_;
414 // Write linear system if asked for.
415 const int verbosity = prm_[activeSolverNum_].get("verbosity", 0);
416 const bool write_matrix = verbosity > 10;
417 if (write_matrix) {
418 Helper::writeSystem(simulator_, //simulator is only used to get names
419 getMatrix(),
420 *rhs_,
421 comm_.get());
422 }
423
424 // Solve system.
425 Dune::InverseOperatorResult result;
426 {
427 OPM_TIMEBLOCK(flexibleSolverApply);
428 assert(flexibleSolver_[activeSolverNum_].solver_);
429 flexibleSolver_[activeSolverNum_].solver_->apply(x, *rhs_, result);
430 }
431
432 // Check convergence, iterations etc.
433 checkConvergence(result);
434
435 return converged_;
436 }
437
438
444
446 int iterations () const { return iterations_; }
447
449 const std::any& parallelInformation() const { return parallelInformation_; }
450
451 const CommunicationType* comm() const { return comm_.get(); }
452
453 void setDomainIndex(const int index)
454 {
455 domainIndex_ = index;
456 }
457
458 bool isNlddLocalSolver() const
459 {
460 return parameters_[activeSolverNum_].is_nldd_local_solver_;
461 }
462
463 protected:
464#if HAVE_MPI
465 using Comm = Dune::OwnerOverlapCopyCommunication<int, int>;
466#endif
467
468 void checkConvergence( const Dune::InverseOperatorResult& result ) const
469 {
470 // store number of iterations
471 iterations_ = result.iterations;
472 converged_ = result.converged;
473 if(!converged_){
474 if(result.reduction < parameters_[activeSolverNum_].relaxed_linear_solver_reduction_){
475 std::stringstream ss;
476 ss<< "Full linear solver tolerance not achieved. The reduction is:" << result.reduction
477 << " after " << result.iterations << " iterations ";
478 OpmLog::warning(ss.str());
479 converged_ = true;
480 }
481 }
482 // Check for failure of linear solver.
483 if (!parameters_[activeSolverNum_].ignoreConvergenceFailure_ && !converged_) {
484 const std::string msg("Convergence failure for linear solver.");
485 OPM_THROW_NOLOG(NumericalProblem, msg);
486 }
487 }
488 protected:
489
490 bool isParallel() const {
491#if HAVE_MPI
492 return !forceSerial_ && comm_->communicator().size() > 1;
493#else
494 return false;
495#endif
496 }
497
498 void prepareFlexibleSolver()
499 {
500 OPM_TIMEBLOCK(flexibleSolverPrepare);
501 if (shouldCreateSolver()) {
502 if (!useWellConn_) {
503 if (isNlddLocalSolver()) {
504 auto wellOp = std::make_unique<DomainWellModelAsLinearOperator<WellModel, Vector, Vector>>(simulator_.problem().wellModel());
505 wellOp->setDomainIndex(domainIndex_);
506 flexibleSolver_[activeSolverNum_].wellOperator_ = std::move(wellOp);
507 }
508 else {
509 auto wellOp = std::make_unique<WellModelOperator>(simulator_.problem().wellModel());
510 flexibleSolver_[activeSolverNum_].wellOperator_ = std::move(wellOp);
511 }
512 }
513 std::function<Vector()> weightCalculator = this->getWeightsCalculator(prm_[activeSolverNum_], getMatrix(), pressureIndex);
514 OPM_TIMEBLOCK(flexibleSolverCreate);
515 flexibleSolver_[activeSolverNum_].create(getMatrix(),
516 isParallel(),
517 prm_[activeSolverNum_],
518 pressureIndex,
519 weightCalculator,
520 forceSerial_,
521 *comm_);
522 }
523 else
524 {
525 OPM_TIMEBLOCK(flexibleSolverUpdate);
526 flexibleSolver_[activeSolverNum_].pre_->update();
527 }
528 }
529
530
534 {
535 // Decide if we should recreate the solver or just do
536 // a minimal preconditioner update.
537 if (flexibleSolver_.empty()) {
538 return true;
539 }
540 if (!flexibleSolver_[activeSolverNum_].solver_) {
541 return true;
542 }
543
544 if (flexibleSolver_[activeSolverNum_].pre_->hasPerfectUpdate()) {
545 return false;
546 }
547
548 // For AMG based preconditioners, the hierarchy depends on the matrix values
549 // so it is recreated at certain intervals
550 if (this->parameters_[activeSolverNum_].cpr_reuse_setup_ == 0) {
551 // Always recreate solver.
552 return true;
553 }
554 if (this->parameters_[activeSolverNum_].cpr_reuse_setup_ == 1) {
555 // Recreate solver on the first iteration of every timestep.
556 const int newton_iteration = this->simulator_.model().newtonMethod().numIterations();
557 return newton_iteration == 0;
558 }
559 if (this->parameters_[activeSolverNum_].cpr_reuse_setup_ == 2) {
560 // Recreate solver if the last solve used more than 10 iterations.
561 return this->iterations() > 10;
562 }
563 if (this->parameters_[activeSolverNum_].cpr_reuse_setup_ == 3) {
564 // Never recreate the solver
565 return false;
566 }
567 if (this->parameters_[activeSolverNum_].cpr_reuse_setup_ == 4) {
568 // Recreate solver every 'step' solve calls.
569 const int step = this->parameters_[activeSolverNum_].cpr_reuse_interval_;
570 const bool create = ((solveCount_ % step) == 0);
571 return create;
572 }
573 // If here, we have an invalid parameter.
574 const bool on_io_rank = (simulator_.gridView().comm().rank() == 0);
575 std::string msg = "Invalid value: " + std::to_string(this->parameters_[activeSolverNum_].cpr_reuse_setup_)
576 + " for --cpr-reuse-setup parameter, run with --help to see allowed values.";
577 if (on_io_rank) {
578 OpmLog::error(msg);
579 }
580 throw std::runtime_error(msg);
581
582 return false;
583 }
584
585
586 // Weights to make approximate pressure equations.
587 // Calculated from the storage terms (only) of the
588 // conservation equations, ignoring all other terms.
589 std::function<Vector()> getWeightsCalculator(const PropertyTree& prm,
590 const Matrix& matrix,
591 std::size_t pressIndex) const
592 {
593 std::function<Vector()> weightsCalculator;
594
595 using namespace std::string_literals;
596
597 auto preconditionerType = prm.get("preconditioner.type"s, "cpr"s);
598 if (preconditionerType == "cpr" || preconditionerType == "cprt"
599 || preconditionerType == "cprw" || preconditionerType == "cprwt") {
600 const bool transpose = preconditionerType == "cprt" || preconditionerType == "cprwt";
601 const auto weightsType = prm.get("preconditioner.weight_type"s, "quasiimpes"s);
602 if (weightsType == "quasiimpes") {
603 // weights will be created as default in the solver
604 // assignment p = pressureIndex prevent compiler warning about
605 // capturing variable with non-automatic storage duration
606 weightsCalculator = [matrix, transpose, pressIndex]() {
607 return Amg::getQuasiImpesWeights<Matrix, Vector>(matrix,
608 pressIndex,
609 transpose);
610 };
611 } else if ( weightsType == "trueimpes" ) {
612 weightsCalculator =
613 [this, pressIndex]
614 {
615 Vector weights(rhs_->size());
616 ElementContext elemCtx(simulator_);
617 Amg::getTrueImpesWeights(pressIndex, weights,
618 simulator_.vanguard().gridView(),
619 elemCtx, simulator_.model(),
621 return weights;
622 };
623 } else if (weightsType == "trueimpesanalytic" ) {
624 weightsCalculator =
625 [this, pressIndex]
626 {
627 Vector weights(rhs_->size());
628 ElementContext elemCtx(simulator_);
629 Amg::getTrueImpesWeightsAnalytic(pressIndex, weights,
630 simulator_.vanguard().gridView(),
631 elemCtx, simulator_.model(),
633 return weights;
634 };
635 } else {
636 OPM_THROW(std::invalid_argument,
637 "Weights type " + weightsType +
638 "not implemented for cpr."
639 " Please use quasiimpes, trueimpes or trueimpesanalytic.");
640 }
641 }
642 return weightsCalculator;
643 }
644
645
646 Matrix& getMatrix()
647 {
648 return *matrix_;
649 }
650
651 const Matrix& getMatrix() const
652 {
653 return *matrix_;
654 }
655
656 const Simulator& simulator_;
657 mutable int iterations_;
658 mutable int solveCount_;
659 mutable bool converged_;
660 std::any parallelInformation_;
661
662 // non-const to be able to scale the linear system
663 Matrix* matrix_;
664 Vector *rhs_;
665
666 int activeSolverNum_ = 0;
667 std::vector<detail::FlexibleSolverInfo<Matrix,Vector,CommunicationType>> flexibleSolver_;
668 std::vector<int> overlapRows_;
669 std::vector<int> interiorRows_;
670
671 int domainIndex_ = -1;
672
673 bool useWellConn_;
674
675 std::vector<FlowLinearSolverParameters> parameters_;
676 bool forceSerial_ = false;
677 std::vector<PropertyTree> prm_;
678
679 std::shared_ptr< CommunicationType > comm_;
680 }; // end ISTLSolver
681
682} // namespace Opm
683
684#endif // OPM_ISTLSOLVER_HEADER_INCLUDED
Helper class for grid instantiation of ECL file-format using problems.
Interface class adding the update() method to the preconditioner interface.
Definition PreconditionerWithUpdate.hpp:32
This class solves the fully implicit black-oil system by solving the reduced system (after eliminatin...
Definition ISTLSolver.hpp:144
ISTLSolver(const Simulator &simulator, const FlowLinearSolverParameters &parameters, bool forceSerial=false)
Construct a system solver.
Definition ISTLSolver.hpp:188
int iterations() const
Solve the system of linear equations Ax = b, with A being the combined derivative matrix of the resid...
Definition ISTLSolver.hpp:446
bool shouldCreateSolver() const
Return true if we should (re)create the whole solver, instead of just calling update() on the precond...
Definition ISTLSolver.hpp:533
const std::any & parallelInformation() const
Definition ISTLSolver.hpp:449
ISTLSolver(const Simulator &simulator)
Construct a system solver.
Definition ISTLSolver.hpp:203
A sparse matrix interface backend for BCRSMatrix from dune-istl.
Definition istlsparsematrixadapter.hh:43
Definition MatrixMarketSpecializations.hpp:17
Definition PropertyTree.hpp:37
static unsigned threadId()
Return the index of the current OpenMP thread.
Definition threadmanager.cpp:84
Definition WellOperators.hpp:70
Declare the properties used by the infrastructure code of the finite volume discretizations.
A sparse matrix interface backend for BCRSMatrix from dune-istl.
Defines the common properties required by the porous medium multi-phase models.
This file contains a set of helper functions used by VFPProd / VFPInj.
Definition blackoilboundaryratevector.hh:37
typename Properties::Detail::GetPropImpl< TypeTag, Property >::type::type GetPropType
get the type alias defined in the property (equivalent to old macro GET_PROP_TYPE(....
Definition propertysystem.hh:235
constexpr auto getPropValue()
get the value data member of a property
Definition propertysystem.hh:242
PropertyTree setupPropertyTree(FlowLinearSolverParameters p, bool linearSolverMaxIterSet, bool linearSolverReductionSet)
Set up a property tree intended for FlexibleSolver by either reading the tree from a JSON file or cre...
Definition setupPropertyTree.cpp:40
This file provides the infrastructure to retrieve run-time parameters.
The Opm property system, traits with inheritance.
This class carries all parameters for the NewtonIterationBlackoilInterleaved class.
Definition FlowLinearSolverParameters.hpp:95
The class that allows to manipulate sparse matrices.
Definition linalgproperties.hh:50
Definition ISTLSolver.hpp:63
Definition ISTLSolver.hpp:69
Definition ISTLSolver.hpp:96