Fix various errors and clean up the ISTLSolvers.

Much of this went undiscovered because the ISTLSolverEbosBda class
was not really used before due to the typo disabling it.
This commit is contained in:
Atgeirr Flø Rasmussen 2023-09-22 15:48:52 +02:00
parent 2ffc7f8e6d
commit e4e32defc4
2 changed files with 45 additions and 52 deletions

View File

@ -184,7 +184,6 @@ std::unique_ptr<Matrix> blockJacobiAdjacency(const Grid& grid,
ISTLSolverEbos(const Simulator& simulator, const FlowLinearSolverParameters& parameters, bool forceSerial = false)
: simulator_(simulator),
iterations_( 0 ),
calls_( 0 ),
converged_(false),
matrix_(nullptr),
parameters_{parameters},
@ -198,7 +197,6 @@ std::unique_ptr<Matrix> blockJacobiAdjacency(const Grid& grid,
explicit ISTLSolverEbos(const Simulator& simulator)
: simulator_(simulator),
iterations_( 0 ),
calls_( 0 ),
solveCount_(0),
converged_(false),
matrix_(nullptr)
@ -208,7 +206,7 @@ std::unique_ptr<Matrix> blockJacobiAdjacency(const Grid& grid,
initialize();
}
void initialize(bool have_gpu = false)
void initialize()
{
OPM_TIMEBLOCK(IstlSolverEbos);
@ -284,12 +282,6 @@ std::unique_ptr<Matrix> blockJacobiAdjacency(const Grid& grid,
}
OpmLog::note(os.str());
}
if(have_gpu){
if (EWOMS_GET_PARAM(TypeTag, std::string, AcceleratorMode) != "none") {
OPM_THROW(std::logic_error,"Cannot use accelerated solver since CUDA, OpenCL and amgcl were not found by cmake");
}
}
}
// nothing to clean here
@ -384,7 +376,6 @@ std::unique_ptr<Matrix> blockJacobiAdjacency(const Grid& grid,
{
OPM_TIMEBLOCK(istlSolverEbosSolve);
++solveCount_;
calls_ += 1;
// Write linear system if asked for.
const int verbosity = prm_[activeSolverNum_].get("verbosity", 0);
const bool write_matrix = verbosity > 10;
@ -522,7 +513,7 @@ std::unique_ptr<Matrix> blockJacobiAdjacency(const Grid& grid,
if (this->parameters_[activeSolverNum_].cpr_reuse_setup_ == 4) {
// Recreate solver every 'step' solve calls.
const int step = this->parameters_[activeSolverNum_].cpr_reuse_interval_;
const bool create = ((calls_ % step) == 0);
const bool create = ((solveCount_ % step) == 0);
return create;
}
@ -567,7 +558,6 @@ std::unique_ptr<Matrix> blockJacobiAdjacency(const Grid& grid,
const Simulator& simulator_;
mutable int iterations_;
mutable int calls_;
mutable int solveCount_;
mutable bool converged_;
std::any parallelInformation_;

View File

@ -131,8 +131,7 @@ public:
ISTLSolverEbosBda(const Simulator& simulator, const FlowLinearSolverParameters& parameters)
: ParentType(simulator, parameters)
{
bool have_gpu = true;
this->initialize(have_gpu);
initializeBda();
}
/// Construct a system solver.
@ -140,59 +139,65 @@ public:
explicit ISTLSolverEbosBda(const Simulator& simulator)
: ParentType(simulator)
{
initializeBda();
}
void initialize()
void initializeBda()
{
OPM_TIMEBLOCK(initialize);
ParentType::initialize(false);
const bool on_io_rank = (this->simulator_.gridView().comm().rank() == 0);
{
std::string accelerator_mode = EWOMS_GET_PARAM(TypeTag, std::string, AcceleratorMode);
if ((this->simulator_.vanguard().grid().comm().size() > 1) && (accelerator_mode != "none")) {
if (on_io_rank) {
OpmLog::warning("Cannot use GPU with MPI, GPU is disabled");
}
accelerator_mode = "none";
OPM_TIMEBLOCK(initializeBda);
std::string accelerator_mode = EWOMS_GET_PARAM(TypeTag, std::string, AcceleratorMode);
// 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);
if (on_io_rank) {
OpmLog::warning("Cannot use AcceleratorMode feature with MPI, setting AcceleratorMode to 'none'.");
}
const int platformID = EWOMS_GET_PARAM(TypeTag, int, OpenclPlatformId);
const int deviceID = EWOMS_GET_PARAM(TypeTag, int, BdaDeviceId);
const int maxit = EWOMS_GET_PARAM(TypeTag, int, LinearSolverMaxIter);
const double tolerance = EWOMS_GET_PARAM(TypeTag, double, LinearSolverReduction);
const bool opencl_ilu_parallel = EWOMS_GET_PARAM(TypeTag, bool, OpenclIluParallel);
const int linear_solver_verbosity = this->parameters_.linear_solver_verbosity_;
std::string linsolver = EWOMS_GET_PARAM(TypeTag, std::string, LinearSolver);
bdaBridge_ = std::make_unique<detail::BdaSolverInfo<Matrix,Vector>>(accelerator_mode,
linear_solver_verbosity,
maxit,
tolerance,
platformID,
deviceID,
opencl_ilu_parallel,
linsolver);
accelerator_mode = "none";
}
if (accelerator_mode == "none") {
return;
}
// Initialize the BdaBridge
const int platformID = EWOMS_GET_PARAM(TypeTag, int, OpenclPlatformId);
const int deviceID = EWOMS_GET_PARAM(TypeTag, int, BdaDeviceId);
const int maxit = EWOMS_GET_PARAM(TypeTag, int, LinearSolverMaxIter);
const double tolerance = EWOMS_GET_PARAM(TypeTag, double, LinearSolverReduction);
const bool opencl_ilu_parallel = EWOMS_GET_PARAM(TypeTag, bool, OpenclIluParallel);
const int linear_solver_verbosity = this->parameters_[0].linear_solver_verbosity_;
std::string linsolver = EWOMS_GET_PARAM(TypeTag, std::string, LinearSolver);
bdaBridge_ = std::make_unique<detail::BdaSolverInfo<Matrix,Vector>>(accelerator_mode,
linear_solver_verbosity,
maxit,
tolerance,
platformID,
deviceID,
opencl_ilu_parallel,
linsolver);
}
void prepare(const Matrix& M, Vector& b)
{
OPM_TIMEBLOCK(prepare);
ParentType::prepare(M,b);
const bool firstcall = (this->matrix_ == nullptr);
ParentType::prepare(M,b);
#if HAVE_OPENCL
// update matrix entries for solvers.
if (firstcall) {
if (firstcall && bdaBridge_) {
// ebos will not change the matrix object. Hence simply store a pointer
// to the original one with a deleter that does nothing.
// Outch! We need to be able to scale the linear system! Hence const_cast
// setup sparsity pattern for jacobi matrix for preconditioner (only used for openclSolver)
#if HAVE_OPENCL
bdaBridge_->numJacobiBlocks_ = EWOMS_GET_PARAM(TypeTag, int, NumJacobiBlocks);
bdaBridge_->prepare(this->simulator_.vanguard().grid(),
this->simulator_.vanguard().cartesianIndexMapper(),
this->simulator_.vanguard().schedule().getWellsatEnd(),
this->simulator_.vanguard().cellPartition(),
this->getMatrix().nonzeroes(), this->useWellConn_);
#endif
}
#endif
}
@ -213,8 +218,12 @@ public:
bool solve(Vector& x)
{
OPM_TIMEBLOCK(solve);
this->calls_ += 1;
if (!bdaBridge_) {
return ParentType::solve(x);
}
OPM_TIMEBLOCK(istlSolverEbosBdaSolve);
this->solveCount_ += 1;
// Write linear system if asked for.
const int verbosity = this->prm_[this->activeSolverNum_].template get<int>("verbosity", 0);
const bool write_matrix = verbosity > 10;
@ -253,12 +262,6 @@ public:
}
protected:
void prepareFlexibleSolver()
{
if(bdaBridge_->gpuActive()){
ParentType::prepareFlexibleSolver();
}
}
std::unique_ptr<detail::BdaSolverInfo<Matrix, Vector>> bdaBridge_;
}; // end ISTLSolver