From f2ed2b6dc3d8f3147184807946775cd6119d662e Mon Sep 17 00:00:00 2001 From: hnil Date: Tue, 20 Aug 2019 14:59:07 +0200 Subject: [PATCH] added writing of coarse level matrix in cpr and corrected coarsening criterio in flexible solver --- opm/simulators/linalg/BlackoilAmg.hpp | 8 +++ opm/simulators/linalg/FlexibleSolver.hpp | 8 +++ .../linalg/PreconditionerFactory.hpp | 51 +++++++++++++++++-- .../linalg/PressureSolverPolicy.hpp | 10 ++++ 4 files changed, 72 insertions(+), 5 deletions(-) diff --git a/opm/simulators/linalg/BlackoilAmg.hpp b/opm/simulators/linalg/BlackoilAmg.hpp index e0963cdf3..a11a9d3f4 100644 --- a/opm/simulators/linalg/BlackoilAmg.hpp +++ b/opm/simulators/linalg/BlackoilAmg.hpp @@ -35,6 +35,7 @@ #include #include #include +#include "matrixmarket_ewoms.hh" namespace Dune { namespace Amg @@ -553,6 +554,13 @@ public: CoarseLevelSolver* createCoarseLevelSolver(LTP& transferPolicy) { coarseOperator_=transferPolicy.getCoarseLevelOperator(); + if( param_->cpr_solver_verbose_ >10 ){ + std::string filename = "course_matrix_istl_blackoil.txt"; + std::ofstream filem(filename); + //const auto& matrix = *transferPolicy.getCoarseMatrix(); + const auto& matrix = coarseOperator_->getmat(); + Dune::writeMatrixMarket(matrix, filem); + } const LevelTransferPolicy& transfer = reinterpret_cast(transferPolicy); AMGInverseOperator* inv = new AMGInverseOperator(param_, diff --git a/opm/simulators/linalg/FlexibleSolver.hpp b/opm/simulators/linalg/FlexibleSolver.hpp index 2f3f63637..657e0b9e2 100644 --- a/opm/simulators/linalg/FlexibleSolver.hpp +++ b/opm/simulators/linalg/FlexibleSolver.hpp @@ -130,6 +130,14 @@ private: tol, // desired residual reduction factor maxiter, // maximum number of iterations verbosity)); + } else if (solver_type == "cg") { + linsolver_.reset(new Dune::CGSolver(*linearoperator_, + *scalarproduct_, + *preconditioner_, + tol, // desired residual reduction factor + maxiter, // maximum number of iterations + verbosity)); + } else if (solver_type == "gmres") { int restart = prm.get("restart"); linsolver_.reset(new Dune::RestartedGMResSolver(*linearoperator_, diff --git a/opm/simulators/linalg/PreconditionerFactory.hpp b/opm/simulators/linalg/PreconditionerFactory.hpp index 859f0a8c1..70abc601c 100644 --- a/opm/simulators/linalg/PreconditionerFactory.hpp +++ b/opm/simulators/linalg/PreconditionerFactory.hpp @@ -29,6 +29,7 @@ #include #include +#include #include #include @@ -104,8 +105,10 @@ public: } private: + //using CriterionBase + // = Dune::Amg::AggregationCriterion>; using CriterionBase - = Dune::Amg::AggregationCriterion>; + = Dune::Amg::AggregationCriterion>; using Criterion = Dune::Amg::CoarsenCriterion; // Helpers for creation of AMG preconditioner. @@ -116,7 +119,9 @@ private: criterion.setAlpha(prm.get("alpha")); criterion.setBeta(prm.get("beta")); criterion.setMaxLevel(prm.get("maxlevel")); - criterion.setSkipIsolated(false); + criterion.setSkipIsolated(prm.get("skip_isolated")); + criterion.setNoPreSmoothSteps(prm.get("pre_smooth")); + criterion.setNoPostSmoothSteps(prm.get("post_smooth")); criterion.setDebugLevel(prm.get("verbosity")); return criterion; } @@ -135,11 +140,21 @@ private: } template - static PrecPtr makeAmgPreconditioner(const Operator& op, const boost::property_tree::ptree& prm) + static PrecPtr makeAmgPreconditioner(const Operator& op, const boost::property_tree::ptree& prm,bool use_kamg = false) { auto crit = amgCriterion(prm); auto sargs = amgSmootherArgs(prm); - return std::make_shared>(op, crit, sargs); + if(use_kamg){ + return std::make_shared< + DummyUpdatePreconditioner< + Dune::Amg::KAMG< Operator, Vector, Smoother> + > + >(op, crit, sargs, + prm.get("max_krylov"), + prm.get("min_reduction") ); + }else{ + return std::make_shared>(op, crit, sargs); + } } // Add a useful default set of preconditioners to the factory. @@ -191,7 +206,7 @@ private: }); doAddCreator("amg", [](const O& op, const P& prm, const C& comm) { const std::string smoother = prm.get("smoother"); - if (smoother == "ILU0") { + if (smoother == "ParOverILU0") { using Smoother = Opm::ParallelOverlappingILU0; auto crit = amgCriterion(prm); auto sargs = amgSmootherArgs(prm); @@ -275,6 +290,32 @@ private: throw std::runtime_error(msg); } }); + doAddCreator("kamg", [](const O& op, const P& prm) { + const std::string smoother = prm.get("smoother"); + if (smoother == "ILU0") { + using Smoother = SeqILU0; + return makeAmgPreconditioner(op, prm,true); + } else if (smoother == "Jac") { + using Smoother = SeqJac; + return makeAmgPreconditioner(op, prm,true); + } else if (smoother == "SOR") { + using Smoother = SeqSOR; + return makeAmgPreconditioner(op, prm,true); + // } else if (smoother == "GS") { + // using Smoother = SeqGS; + // return makeAmgPreconditioner(op, prm,true); + } else if (smoother == "SSOR") { + using Smoother = SeqSSOR; + return makeAmgPreconditioner(op, prm,true); + } else if (smoother == "ILUn") { + using Smoother = SeqILUn; + return makeAmgPreconditioner(op, prm,true); + } else { + std::string msg("No such smoother: "); + msg += smoother; + throw std::runtime_error(msg); + } + }); doAddCreator("famg", [](const O& op, const P& prm) { auto crit = amgCriterion(prm); Dune::Amg::Parameters parms; diff --git a/opm/simulators/linalg/PressureSolverPolicy.hpp b/opm/simulators/linalg/PressureSolverPolicy.hpp index 5692031ca..d5d845dc9 100644 --- a/opm/simulators/linalg/PressureSolverPolicy.hpp +++ b/opm/simulators/linalg/PressureSolverPolicy.hpp @@ -113,6 +113,16 @@ namespace Amg auto& tp = dynamic_cast(transferPolicy); // TODO: make this unnecessary. PressureInverseOperator* inv = new PressureInverseOperator(*coarseOperator_, prm_, tp.getCoarseLevelCommunication()); + if (prm_.get("verbosity") > 10) { + std::string filename = "course_matrix_istl_pressure_flex.txt"; + std::ofstream filem(filename); + if (!filem) { + throw std::runtime_error("Could not write matrix"); + } + //const auto& matrix = *levelTransferPolicy_.getCoarseMatrix(); + const auto& matrix = coarseOperator_->getmat(); + Dune::writeMatrixMarket(matrix, filem); + } return inv; }