OilfieldToolsNet/opm-simulators
4
0
Fork 0
mirror of https://github.com/OPM/opm-simulators.git synced 2025-02-25 18:55:30 -06:00
Code Issues Packages Projects Releases Wiki Activity

Add paralell preconditioner, enable parallel case.

Browse Source
This commit is contained in:
Atgeirr Flø Rasmussen 2015-06-16 13:28:11 +02:00
parent a3d115ff22
commit 5e513642d7
2 changed files with 49 additions and 14 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

2
opm/autodiff/NewtonIterationBlackoilInterleaved.cpp
View File

@ -164,7 +164,7 @@ namespace Opm
Dune::InverseOperatorResult result;
// Parallel version is deactivated until we figure out how to do it properly.
#if 0 // HAVE_MPI
#if HAVE_MPI
if (parallelInformation_.type() == typeid(ParallelISTLInformation))
{
typedef Dune::OwnerOverlapCopyCommunication<int,int> Comm;

61
opm/autodiff/NewtonIterationBlackoilInterleaved.hpp
View File

@ -27,6 +27,7 @@
#include <opm/autodiff/NewtonIterationBlackoilInterface.hpp>
#include <opm/core/utility/parameters/ParameterGroup.hpp>
#include <opm/core/linalg/ParallelISTLInformation.hpp>
#include <opm/core/utility/platform_dependent/disable_warnings.h>
@ -34,6 +35,7 @@
#include <dune/istl/operators.hh>
#include <dune/istl/preconditioners.hh>
#include <dune/istl/solvers.hh>
#include <dune/istl/owneroverlapcopy.hh>
#include <dune/istl/paamg/amg.hh>
#include <opm/core/utility/platform_dependent/reenable_warnings.h>
@ -87,33 +89,66 @@ namespace Opm
/// \brief construct the CPR preconditioner and the solver.
/// \tparam P The type of the parallel information.
/// \param parallelInformation the information about the parallelization.
template<int category=Dune::SolverCategory::sequential, class O, class P>
template<int category=Dune::SolverCategory::sequential, class O, class POrComm>
void constructPreconditionerAndSolve(O& opA,
Vector& x, Vector& istlb,
const P& parallelInformation,
const POrComm& parallelInformation,
Dune::InverseOperatorResult& result) const
{
typedef Dune::ScalarProductChooser<Vector,P,category> ScalarProductChooser;
std::unique_ptr<typename ScalarProductChooser::ScalarProduct>
sp(ScalarProductChooser::construct(parallelInformation));
// Construct preconditioner.
typedef Dune::SeqILU0<Mat,Vector,Vector> Preconditioner;
// typedef Opm::CPRPreconditioner<Mat,Vector,Vector,P> Preconditioner;
parallelInformation.copyOwnerToAll(istlb, istlb);
const double relax = 1.0;
Preconditioner precond(opA.getmat(), relax);
// Construct scalar product.
typedef Dune::ScalarProductChooser<Vector, POrComm, category> ScalarProductChooser;
auto sp = ScalarProductChooser::construct(parallelInformation);
// Construct preconditioner.
auto precond = constructPrecond(opA, parallelInformation);
// Communicate if parallel.
parallelInformation.copyOwnerToAll(istlb, istlb);
// Solve.
solve(opA, x, istlb, *sp, precond, result);
}
typedef Dune::SeqILU0<Mat, Vector, Vector> SeqPreconditioner;
typedef Dune::OwnerOverlapCopyCommunication<int, int> Comm;
typedef Dune::BlockPreconditioner<Vector, Vector, Comm, SeqPreconditioner> ParPreconditioner;
template <class Operator>
SeqPreconditioner constructPrecond(Operator& opA, const Dune::Amg::SequentialInformation&) const
{
const double relax = 1.0;
SeqPreconditioner precond(opA.getmat(), relax);
return precond;
}
template <class Operator>
ParPreconditioner constructPrecond(Operator& opA, const Comm& comm) const
{
const double relax = 1.0;
SeqPreconditioner seq_precond(opA.getmat(), relax);
ParPreconditioner precond(seq_precond, comm);
return precond;
}
/// \brief Solve the system using the given preconditioner and scalar product.
template <class Operator, class ScalarProd, class Precond>
void solve(Operator& opA, Vector& x, Vector& istlb, ScalarProd& sp, Precond& precond, Dune::InverseOperatorResult& result) const
{
// TODO: Revise when linear solvers interface opm-core is done
// Construct linear solver.
// GMRes solver
if ( newton_use_gmres_ ) {
Dune::RestartedGMResSolver<Vector> linsolve(opA, *sp, precond,
Dune::RestartedGMResSolver<Vector> linsolve(opA, sp, precond,
linear_solver_reduction_, linear_solver_restart_, linear_solver_maxiter_, linear_solver_verbosity_);
// Solve system.
linsolve.apply(x, istlb, result);
}
else { // BiCGstab solver
Dune::BiCGSTABSolver<Vector> linsolve(opA, *sp, precond,
Dune::BiCGSTABSolver<Vector> linsolve(opA, sp, precond,
linear_solver_reduction_, linear_solver_maxiter_, linear_solver_verbosity_);
// Solve system.
linsolve.apply(x, istlb, result);