From 3751ca6be45d19baa931b4b843ea6b820b3c653f Mon Sep 17 00:00:00 2001 From: Arne Morten Kvarving Date: Sat, 6 Jul 2024 10:22:47 +0200 Subject: [PATCH] move FlowLinearSolver parameters to TypeTag-free parameter system --- opm/simulators/flow/BlackoilModelNldd.hpp | 2 +- .../linalg/FlowLinearSolverParameters.hpp | 254 +++++------------- opm/simulators/linalg/ISTLSolver.hpp | 20 +- opm/simulators/linalg/ISTLSolverBda.hpp | 14 +- 4 files changed, 80 insertions(+), 210 deletions(-) diff --git a/opm/simulators/flow/BlackoilModelNldd.hpp b/opm/simulators/flow/BlackoilModelNldd.hpp index 4fc9a0923..ef4786017 100644 --- a/opm/simulators/flow/BlackoilModelNldd.hpp +++ b/opm/simulators/flow/BlackoilModelNldd.hpp @@ -160,7 +160,7 @@ public: // only. This must be addressed before going parallel. const auto& eclState = model_.simulator().vanguard().eclState(); FlowLinearSolverParameters loc_param; - loc_param.template init(eclState.getSimulationConfig().useCPR()); + loc_param.init(eclState.getSimulationConfig().useCPR()); // Override solver type with umfpack if small domain. // Otherwise hardcode to ILU0 if (domains_[index].cells.size() < 200) { diff --git a/opm/simulators/linalg/FlowLinearSolverParameters.hpp b/opm/simulators/linalg/FlowLinearSolverParameters.hpp index 78414de08..8a56e9278 100644 --- a/opm/simulators/linalg/FlowLinearSolverParameters.hpp +++ b/opm/simulators/linalg/FlowLinearSolverParameters.hpp @@ -63,154 +63,26 @@ struct LinearSolverBackend namespace Opm::Parameters { -template -struct LinearSolverReduction { using type = Properties::UndefinedProperty; }; - -template -struct RelaxedLinearSolverReduction { using type = Properties::UndefinedProperty; }; - -template -struct LinearSolverMaxIter { using type = Properties::UndefinedProperty; }; - -template -struct LinearSolverRestart { using type = Properties::UndefinedProperty; }; - -template -struct IluRelaxation { using type = Properties::UndefinedProperty; }; - -template -struct IluFillinLevel { using type = Properties::UndefinedProperty; }; - -template -struct MiluVariant { using type = Properties::UndefinedProperty; }; - -template -struct IluRedblack { using type = Properties::UndefinedProperty; }; - -template -struct IluReorderSpheres { using type = Properties::UndefinedProperty; }; - -template -struct UseGmres { using type = Properties::UndefinedProperty; }; - -template -struct LinearSolverIgnoreConvergenceFailure { using type = Properties::UndefinedProperty; }; - -template -struct ScaleLinearSystem { using type = Properties::UndefinedProperty; }; - -template -struct LinearSolver { using type = Properties::UndefinedProperty; }; - -template -struct LinearSolverPrintJsonDefinition { using type = Properties::UndefinedProperty; }; - -template -struct CprReuseSetup { using type = Properties::UndefinedProperty; }; - -template -struct CprReuseInterval { using type = Properties::UndefinedProperty; }; - -template -struct AcceleratorMode { using type = Properties::UndefinedProperty; }; - -template -struct BdaDeviceId { using type = Properties::UndefinedProperty; }; - -template -struct OpenclPlatformId { using type = Properties::UndefinedProperty; }; - -template -struct OpenclIluParallel { using type = Properties::UndefinedProperty; }; - -template -struct LinearSolverReduction -{ - using type = GetPropType; - static constexpr type value = 1e-2; -}; - -template -struct RelaxedLinearSolverReduction -{ - using type = GetPropType; - static constexpr type value = 1e-2; -}; - -template -struct LinearSolverMaxIter -{ static constexpr int value = 200; }; - -template -struct LinearSolverRestart -{ static constexpr int value = 40; }; - -template -struct IluRelaxation -{ - using type = GetPropType; - static constexpr type value = 0.9; -}; - -template -struct IluFillinLevel -{ static constexpr int value = 0; }; - -template -struct MiluVariant -{ static constexpr auto value = "ILU"; }; - -template -struct IluRedblack -{ static constexpr bool value = false; }; - -template -struct IluReorderSpheres -{ static constexpr bool value = false; }; - -template -struct UseGmres -{ static constexpr bool value = false; }; - -template -struct LinearSolverIgnoreConvergenceFailure -{ static constexpr bool value = false; }; - -template -struct ScaleLinearSystem -{ static constexpr bool value = false; }; - -template -struct LinearSolver -{ static constexpr auto value = "cprw"; }; - -template -struct LinearSolverPrintJsonDefinition -{ static constexpr auto value = true; }; - -template -struct CprReuseSetup -{ static constexpr int value = 4; }; - -template -struct CprReuseInterval -{ static constexpr int value = 30; }; - -template -struct AcceleratorMode -{ static constexpr auto value = "none"; }; - -template -struct BdaDeviceId -{ static constexpr int value = 0; }; - -template -struct OpenclPlatformId -{ static constexpr int value = 0; }; - -template -struct OpenclIluParallel -{ static constexpr bool value = true; }; // note: false should only be used in debug +struct LinearSolverReduction { static constexpr double value = 1e-2; }; +struct RelaxedLinearSolverReduction { static constexpr double value = 1e-2; }; +struct IluRelaxation { static constexpr double value = 0.9; }; +struct LinearSolverMaxIter { static constexpr int value = 200; }; +struct LinearSolverRestart { static constexpr int value = 40; }; +struct IluFillinLevel { static constexpr int value = 0; }; +struct MiluVariant { static constexpr auto value = "ILU"; }; +struct IluRedblack { static constexpr bool value = false; }; +struct IluReorderSpheres { static constexpr bool value = false; }; +struct UseGmres { static constexpr bool value = false; }; +struct LinearSolverIgnoreConvergenceFailure { static constexpr bool value = false; }; +struct ScaleLinearSystem { static constexpr bool value = false; }; +struct LinearSolver { static constexpr auto value = "cprw"; }; +struct LinearSolverPrintJsonDefinition { static constexpr auto value = true; }; +struct CprReuseSetup { static constexpr int value = 4; }; +struct CprReuseInterval { static constexpr int value = 30; }; +struct AcceleratorMode { static constexpr auto value = "none"; }; +struct BdaDeviceId { static constexpr int value = 0; }; +struct OpenclPlatformId { static constexpr int value = 0; }; +struct OpenclIluParallel { static constexpr bool value = true; }; // note: false should only be used in debug } // namespace Opm::Parameters @@ -241,109 +113,107 @@ struct FlowLinearSolverParameters int opencl_platform_id_; bool opencl_ilu_parallel_; - template void init(bool cprRequestedInDataFile) { // TODO: these parameters have undocumented non-trivial dependencies - linear_solver_reduction_ = Parameters::get(); - relaxed_linear_solver_reduction_ = Parameters::get(); - linear_solver_maxiter_ = Parameters::get(); - linear_solver_restart_ = Parameters::get(); + linear_solver_reduction_ = Parameters::Get(); + relaxed_linear_solver_reduction_ = Parameters::Get(); + linear_solver_maxiter_ = Parameters::Get(); + linear_solver_restart_ = Parameters::Get(); linear_solver_verbosity_ = Parameters::Get(); - ilu_relaxation_ = Parameters::get(); - ilu_fillin_level_ = Parameters::get(); - ilu_milu_ = convertString2Milu(Parameters::get()); - ilu_redblack_ = Parameters::get(); - ilu_reorder_sphere_ = Parameters::get(); - newton_use_gmres_ = Parameters::get(); - ignoreConvergenceFailure_ = Parameters::get(); - scale_linear_system_ = Parameters::get(); - linsolver_ = Parameters::get(); - linear_solver_print_json_definition_ = Parameters::get(); - cpr_reuse_setup_ = Parameters::get(); - cpr_reuse_interval_ = Parameters::get(); + ilu_relaxation_ = Parameters::Get(); + ilu_fillin_level_ = Parameters::Get(); + ilu_milu_ = convertString2Milu(Parameters::Get()); + ilu_redblack_ = Parameters::Get(); + ilu_reorder_sphere_ = Parameters::Get(); + newton_use_gmres_ = Parameters::Get(); + ignoreConvergenceFailure_ = Parameters::Get(); + scale_linear_system_ = Parameters::Get(); + linsolver_ = Parameters::Get(); + linear_solver_print_json_definition_ = Parameters::Get(); + cpr_reuse_setup_ = Parameters::Get(); + cpr_reuse_interval_ = Parameters::Get(); - if (!Parameters::isSet() && cprRequestedInDataFile) { + if (!Parameters::IsSet() && cprRequestedInDataFile) { linsolver_ = "cpr"; } else { - linsolver_ = Parameters::get(); + linsolver_ = Parameters::Get(); } - accelerator_mode_ = Parameters::get(); - bda_device_id_ = Parameters::get(); - opencl_platform_id_ = Parameters::get(); - opencl_ilu_parallel_ = Parameters::get(); + accelerator_mode_ = Parameters::Get(); + bda_device_id_ = Parameters::Get(); + opencl_platform_id_ = Parameters::Get(); + opencl_ilu_parallel_ = Parameters::Get(); } - template static void registerParameters() { - Parameters::registerParam + Parameters::Register ("The minimum reduction of the residual which the linear solver must achieve"); - Parameters::registerParam + Parameters::Register ("The minimum reduction of the residual which the linear solver need to " "achieve for the solution to be accepted"); - Parameters::registerParam + Parameters::Register ("The maximum number of iterations of the linear solver"); - Parameters::registerParam + Parameters::Register ("The number of iterations after which GMRES is restarted"); Parameters::Register ("The verbosity level of the linear solver (0: off, 2: all)"); - Parameters::registerParam + Parameters::Register ("The relaxation factor of the linear solver's ILU preconditioner"); - Parameters::registerParam + Parameters::Register ("The fill-in level of the linear solver's ILU preconditioner"); - Parameters::registerParam + Parameters::Register ("Specify which variant of the modified-ILU preconditioner ought to be used. " "Possible variants are: ILU (default, plain ILU), " "MILU_1 (lump diagonal with dropped row entries), " "MILU_2 (lump diagonal with the sum of the absolute values of the dropped row entries), " "MILU_3 (if diagonal is positive add sum of dropped row entries, otherwise subtract them), " "MILU_4 (if diagonal is positive add sum of dropped row entries, otherwise do nothing"); - Parameters::registerParam + Parameters::Register ("Use red-black partitioning for the ILU preconditioner"); - Parameters::registerParam + Parameters::Register ("Whether to reorder the entries of the matrix in the red-black " "ILU preconditioner in spheres starting at an edge. " "If false the original ordering is preserved in each color. " "Otherwise why try to ensure D4 ordering (in a 2D structured grid, " "the diagonal elements are consecutive)."); - Parameters::registerParam + Parameters::Register ("Use GMRES as the linear solver"); - Parameters::registerParam + Parameters::Register ("Continue with the simulation like nothing happened " "after the linear solver did not converge"); - Parameters::registerParam + Parameters::Register ("Scale linear system according to equation scale and primary variable types"); - Parameters::registerParam + Parameters::Register ("Configuration of solver. Valid options are: ilu0 (default), " "dilu, cprw, cpr (an alias for cprw), cpr_quasiimpes, " "cpr_trueimpes, cpr_trueimpesanalytic, amg or hybrid (experimental). " "Alternatively, you can request a configuration to be read from a " "JSON file by giving the filename here, ending with '.json.'"); - Parameters::registerParam + Parameters::Register ("Write the JSON definition of the linear solver setup to the DBG file."); - Parameters::registerParam + Parameters::Register ("Reuse preconditioner setup. Valid options are " "0: recreate the preconditioner for every linear solve, " "1: recreate once every timestep, " "2: recreate if last linear solve took more than 10 iterations, " "3: never recreate, " "4: recreated every CprReuseInterval"); - Parameters::registerParam + Parameters::Register ("Reuse preconditioner interval. Used when CprReuseSetup is set to 4, " "then the preconditioner will be fully recreated instead of reused " "every N linear solve, where N is this parameter."); - Parameters::registerParam + Parameters::Register ("Choose a linear solver, usage: " "'--accelerator-mode=[none|cusparse|opencl|amgcl|rocalution|rocsparse]'"); - Parameters::registerParam + Parameters::Register ("Choose device ID for cusparseSolver or openclSolver, " "use 'nvidia-smi' or 'clinfo' to determine valid IDs"); - Parameters::registerParam + Parameters::Register ("Choose platform ID for openclSolver, use 'clinfo' " "to determine valid platform IDs"); - Parameters::registerParam + Parameters::Register ("Parallelize ILU decomposition and application on GPU"); Parameters::SetDefault(0); diff --git a/opm/simulators/linalg/ISTLSolver.hpp b/opm/simulators/linalg/ISTLSolver.hpp index b65607267..cecf0aec2 100644 --- a/opm/simulators/linalg/ISTLSolver.hpp +++ b/opm/simulators/linalg/ISTLSolver.hpp @@ -170,7 +170,7 @@ std::unique_ptr blockJacobiAdjacency(const Grid& grid, static void registerParameters() { - FlowLinearSolverParameters::registerParameters(); + FlowLinearSolverParameters::registerParameters(); } /// Construct a system solver. @@ -203,7 +203,7 @@ std::unique_ptr blockJacobiAdjacency(const Grid& grid, matrix_(nullptr) { parameters_.resize(1); - parameters_[0].template init(simulator_.vanguard().eclState().getSimulationConfig().useCPR()); + parameters_[0].init(simulator_.vanguard().eclState().getSimulationConfig().useCPR()); initialize(); } @@ -220,21 +220,21 @@ std::unique_ptr blockJacobiAdjacency(const Grid& grid, parameters_.clear(); { FlowLinearSolverParameters para; - para.init(false); + para.init(false); para.linsolver_ = "cprw"; parameters_.push_back(para); prm_.push_back(setupPropertyTree(parameters_[0], - Parameters::isSet(), - Parameters::isSet())); + Parameters::IsSet(), + Parameters::IsSet())); } { FlowLinearSolverParameters para; - para.init(false); + para.init(false); para.linsolver_ = "ilu0"; parameters_.push_back(para); prm_.push_back(setupPropertyTree(parameters_[1], - Parameters::isSet(), - Parameters::isSet())); + Parameters::IsSet(), + Parameters::IsSet())); } // ------------ } else { @@ -242,8 +242,8 @@ std::unique_ptr blockJacobiAdjacency(const Grid& grid, assert(parameters_.size() == 1); assert(prm_.empty()); prm_.push_back(setupPropertyTree(parameters_[0], - Parameters::isSet(), - Parameters::isSet())); + Parameters::IsSet(), + Parameters::IsSet())); } flexibleSolver_.resize(prm_.size()); diff --git a/opm/simulators/linalg/ISTLSolverBda.hpp b/opm/simulators/linalg/ISTLSolverBda.hpp index 91691c75a..61cd95d06 100644 --- a/opm/simulators/linalg/ISTLSolverBda.hpp +++ b/opm/simulators/linalg/ISTLSolverBda.hpp @@ -152,7 +152,7 @@ public: { OPM_TIMEBLOCK(initializeBda); - std::string accelerator_mode = Parameters::get(); + std::string accelerator_mode = Parameters::Get(); // Force accelerator mode to none if using MPI. if ((this->simulator_.vanguard().grid().comm().size() > 1) && (accelerator_mode != "none")) { const bool on_io_rank = (this->simulator_.gridView().comm().rank() == 0); @@ -167,13 +167,13 @@ public: } // Initialize the BdaBridge - const int platformID = Parameters::get(); - const int deviceID = Parameters::get(); - const int maxit = Parameters::get(); - const double tolerance = Parameters::get(); - const bool opencl_ilu_parallel = Parameters::get(); + const int platformID = Parameters::Get(); + const int deviceID = Parameters::Get(); + const int maxit = Parameters::Get(); + const double tolerance = Parameters::Get(); + const bool opencl_ilu_parallel = Parameters::Get(); const int linear_solver_verbosity = this->parameters_[0].linear_solver_verbosity_; - std::string linsolver = Parameters::get(); + std::string linsolver = Parameters::Get(); bdaBridge_ = std::make_unique>(accelerator_mode, linear_solver_verbosity, maxit,