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 
26 #include <string>
27 
28 namespace Opm {
29 
30 class WellContributions;
31 
32 namespace Accelerator {
33  enum class SolverStatus {
34  BDA_SOLVER_SUCCESS,
35  BDA_SOLVER_ANALYSIS_FAILED,
36  BDA_SOLVER_CREATE_PRECONDITIONER_FAILED,
37  BDA_SOLVER_UNKNOWN_ERROR
38  };
39 
42  template <unsigned int block_size>
43  class BdaSolver
44  {
45 
46  protected:
47 
48  // verbosity
49  // 0: print nothing during solves, only when initializing
50  // 1: print number of iterations and final norm
51  // 2: also print norm each iteration
52  // 3: also print timings of different backend functions
53 
54  int verbosity = 0;
55 
56  int maxit = 200;
57  double tolerance = 1e-2;
58 
59  std::string bitstream;
60 
61  int N; // number of rows
62  int Nb; // number of blocked rows (Nb*block_size == N)
63  int nnz; // number of nonzeroes (scalars)
64  int nnzb; // number of nonzero blocks (nnzb*block_size*block_size == nnz)
65 
66  unsigned int platformID = 0; // ID of OpenCL platform to be used, only used by openclSolver now
67  unsigned int deviceID = 0; // ID of the device to be used
68 
69  bool initialized = false;
70 
71  public:
79  BdaSolver(int linear_solver_verbosity, int max_it, double tolerance_) : verbosity(linear_solver_verbosity), maxit(max_it), tolerance(tolerance_) {};
80  BdaSolver(int linear_solver_verbosity, int max_it, double tolerance_, unsigned int deviceID_) : verbosity(linear_solver_verbosity), maxit(max_it), tolerance(tolerance_), deviceID(deviceID_) {};
81  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_) {};
82  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) {};
83 
85  virtual ~BdaSolver() {};
86 
88  virtual SolverStatus solve_system(int N, int nnz, int dim,
89  double *vals, int *rows, int *cols,
90  double *b, 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:44
virtual SolverStatus solve_system(int N, int nnz, int dim, double *vals, int *rows, int *cols, double *b, 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:85
BdaSolver(int linear_solver_verbosity, int max_it, double tolerance_)
Construct a BdaSolver, can be cusparseSolver, openclSolver, fpgaSolver.
Definition: BdaSolver.hpp:79
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