My Project
BdaSolver.hpp
1 /*
2  Copyright 2019 Equinor ASA
3 
4  This file is part of the Open Porous Media project (OPM).
5 
6  OPM is free software: you can redistribute it and/or modify
7  it under the terms of the GNU General Public License as published by
8  the Free Software Foundation, either version 3 of the License, or
9  (at your option) any later version.
10 
11  OPM is distributed in the hope that it will be useful,
12  but WITHOUT ANY WARRANTY; without even the implied warranty of
13  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
14  GNU General Public License for more details.
15 
16  You should have received a copy of the GNU General Public License
17  along with OPM. If not, see <http://www.gnu.org/licenses/>.
18 */
19 
20 #ifndef OPM_BDASOLVER_BACKEND_HEADER_INCLUDED
21 #define OPM_BDASOLVER_BACKEND_HEADER_INCLUDED
22 
23 
24 #include <opm/simulators/linalg/bda/BdaResult.hpp>
25 #include <opm/simulators/linalg/bda/BlockedMatrix.hpp>
26 
27 #include <string>
28 
29 namespace Opm {
30 
31 class WellContributions;
32 
33 namespace Accelerator {
34  enum class SolverStatus {
35  BDA_SOLVER_SUCCESS,
36  BDA_SOLVER_ANALYSIS_FAILED,
37  BDA_SOLVER_CREATE_PRECONDITIONER_FAILED,
38  BDA_SOLVER_UNKNOWN_ERROR
39  };
40 
43  template <unsigned int block_size>
44  class BdaSolver
45  {
46 
47  protected:
48 
49  // verbosity
50  // 0: print nothing during solves, only when initializing
51  // 1: print number of iterations and final norm
52  // 2: also print norm each iteration
53  // 3: also print timings of different backend functions
54 
55  int verbosity = 0;
56 
57  int maxit = 200;
58  double tolerance = 1e-2;
59 
60  std::string bitstream;
61 
62  int N; // number of rows
63  int Nb; // number of blocked rows (Nb*block_size == N)
64  int nnz; // number of nonzeroes (scalars)
65  int nnzb; // number of nonzero blocks (nnzb*block_size*block_size == nnz)
66 
67  unsigned int platformID = 0; // ID of OpenCL platform to be used, only used by openclSolver now
68  unsigned int deviceID = 0; // ID of the device to be used
69 
70  bool initialized = false;
71 
72  public:
80  BdaSolver(int linear_solver_verbosity, int max_it, double tolerance_) : verbosity(linear_solver_verbosity), maxit(max_it), tolerance(tolerance_) {};
81  BdaSolver(int linear_solver_verbosity, int max_it, double tolerance_, unsigned int deviceID_) : verbosity(linear_solver_verbosity), maxit(max_it), tolerance(tolerance_), deviceID(deviceID_) {};
82  BdaSolver(int linear_solver_verbosity, int max_it, double tolerance_, unsigned int platformID_, unsigned int deviceID_) : verbosity(linear_solver_verbosity), maxit(max_it), tolerance(tolerance_), platformID(platformID_), deviceID(deviceID_) {};
83  BdaSolver(std::string fpga_bitstream, int linear_solver_verbosity, int max_it, double tolerance_) : verbosity(linear_solver_verbosity), maxit(max_it), tolerance(tolerance_), bitstream(fpga_bitstream) {};
84 
86  virtual ~BdaSolver() {};
87 
89  virtual SolverStatus solve_system(std::shared_ptr<BlockedMatrix> matrix, double *b,
90  std::shared_ptr<BlockedMatrix> jacMatrix, WellContributions& wellContribs, BdaResult &res) = 0;
91 
92  virtual void get_result(double *x) = 0;
93 
94  }; // end class BdaSolver
95 
96 } // namespace Accelerator
97 } // namespace Opm
98 
99 #endif
This class is based on InverseOperatorResult struct from dune/istl/solver.hh It is needed to prevent ...
Definition: BdaResult.hpp:31
This class serves to simplify choosing between different backend solvers, such as cusparseSolver and ...
Definition: BdaSolver.hpp:45
virtual SolverStatus solve_system(std::shared_ptr< BlockedMatrix > matrix, double *b, std::shared_ptr< BlockedMatrix > jacMatrix, WellContributions &wellContribs, BdaResult &res)=0
Define as pure virtual functions, so derivedclass must implement them.
virtual ~BdaSolver()
Define virtual destructor, so that the derivedclass destructor will be called.
Definition: BdaSolver.hpp:86
BdaSolver(int linear_solver_verbosity, int max_it, double tolerance_)
Construct a BdaSolver, can be cusparseSolver, openclSolver, fpgaSolver.
Definition: BdaSolver.hpp:80
This class serves to eliminate the need to include the WellContributions into the matrix (with –matri...
Definition: WellContributions.hpp:53
This file contains a set of helper functions used by VFPProd / VFPInj.
Definition: BlackoilPhases.hpp:27