My Project
BdaBridge.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 BDABRIDGE_HEADER_INCLUDED
21 #define BDABRIDGE_HEADER_INCLUDED
22 
23 #include "dune/istl/solver.hh" // for struct InverseOperatorResult
24 
25 #include <opm/simulators/linalg/bda/BdaSolver.hpp>
26 #include <opm/simulators/linalg/bda/ILUReorder.hpp>
27 
28 namespace Opm
29 {
30 
31 class WellContributions;
32 
33 typedef Dune::InverseOperatorResult InverseOperatorResult;
34 using Opm::Accelerator::ILUReorder;
35 
37 template <class BridgeMatrix, class BridgeVector, int block_size>
38 class BdaBridge
39 {
40 private:
41  int verbosity = 0;
42  bool use_gpu = false;
43  bool use_fpga = false;
44  std::string accelerator_mode;
45  std::unique_ptr<Opm::Accelerator::BdaSolver<block_size> > backend;
46 
47 public:
58  BdaBridge(std::string accelerator_mode, std::string fpga_bitstream, int linear_solver_verbosity, int maxit, double tolerance,
59  unsigned int platformID, unsigned int deviceID, std::string opencl_ilu_reorder, std::string linsolver);
60 
61 
68  void solve_system(BridgeMatrix *mat, BridgeVector &b, WellContributions& wellContribs, InverseOperatorResult &result);
69 
72  void get_result(BridgeVector &x);
73 
76  bool getUseGpu(){
77  return use_gpu;
78  }
79 
84  static void copySparsityPatternFromISTL(const BridgeMatrix& mat, std::vector<int>& h_rows, std::vector<int>& h_cols);
85 
89  void initWellContributions(WellContributions& wellContribs);
90 
93  bool getUseFpga(){
94  return use_fpga;
95  }
96 
98  std::string getAccleratorName(){
99  return accelerator_mode;
100  }
101 
102 }; // end class BdaBridge
103 
104 }
105 
106 #endif
BdaBridge acts as interface between opm-simulators with the BdaSolvers.
Definition: BdaBridge.hpp:39
void solve_system(BridgeMatrix *mat, BridgeVector &b, WellContributions &wellContribs, InverseOperatorResult &result)
Solve linear system, A*x = b.
Definition: BdaBridge.cpp:195
void get_result(BridgeVector &x)
Get the resulting x vector.
Definition: BdaBridge.cpp:264
bool getUseGpu()
Return whether the BdaBridge will use the GPU or not return whether the BdaBridge will use the GPU or...
Definition: BdaBridge.hpp:76
BdaBridge(std::string accelerator_mode, std::string fpga_bitstream, int linear_solver_verbosity, int maxit, double tolerance, unsigned int platformID, unsigned int deviceID, std::string opencl_ilu_reorder, std::string linsolver)
Construct a BdaBridge.
Definition: BdaBridge.cpp:60
static void copySparsityPatternFromISTL(const BridgeMatrix &mat, std::vector< int > &h_rows, std::vector< int > &h_cols)
Store sparsity pattern into vectors.
Definition: BdaBridge.cpp:173
bool getUseFpga()
Return whether the BdaBridge will use the FPGA or not return whether the BdaBridge will use the FPGA ...
Definition: BdaBridge.hpp:93
std::string getAccleratorName()
Return the selected accelerator mode, this is input via the command-line.
Definition: BdaBridge.hpp:98
void initWellContributions(WellContributions &wellContribs)
Initialize the WellContributions object with opencl context and queue those must be set before callin...
Definition: BdaBridge.cpp:271
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