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 "dune/istl/bcrsmatrix.hh"
26 #include <opm/simulators/linalg/matrixblock.hh>
27 
28 #include <opm/simulators/linalg/bda/BdaSolver.hpp>
29 #include <opm/simulators/linalg/bda/ILUReorder.hpp>
30 #include <opm/simulators/linalg/bda/WellContributions.hpp>
31 
32 
33 #if HAVE_FPGA
34 #include <opm/simulators/linalg/bda/FPGASolverBackend.hpp>
35 #endif
36 
37 namespace Opm
38 {
39 
40 typedef Dune::InverseOperatorResult InverseOperatorResult;
41 using bda::ILUReorder;
42 
44 template <class BridgeMatrix, class BridgeVector, int block_size>
45 class BdaBridge
46 {
47 private:
48  int verbosity = 0;
49  bool use_gpu = false;
50  bool use_fpga = false;
51  std::string accelerator_mode;
52  std::unique_ptr<bda::BdaSolver<block_size> > backend;
53 
54 public:
64  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);
65 
66 
73  void solve_system(BridgeMatrix *mat, BridgeVector &b, WellContributions& wellContribs, InverseOperatorResult &result);
74 
77  void get_result(BridgeVector &x);
78 
81  bool getUseGpu(){
82  return use_gpu;
83  }
84 
88  void initWellContributions(WellContributions& wellContribs);
89 
92  bool getUseFpga(){
93  return use_fpga;
94  }
95 
97  std::string getAccleratorName(){
98  return accelerator_mode;
99  }
100 
101 }; // end class BdaBridge
102 
103 }
104 
105 #endif
BdaBridge acts as interface between opm-simulators with the BdaSolvers.
Definition: BdaBridge.hpp:46
void solve_system(BridgeMatrix *mat, BridgeVector &b, WellContributions &wellContribs, InverseOperatorResult &result)
Solve linear system, A*x = b.
Definition: BdaBridge.cpp:193
void get_result(BridgeVector &x)
Get the resulting x vector.
Definition: BdaBridge.cpp:262
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)
Construct a BdaBridge.
Definition: BdaBridge.cpp:56
bool getUseGpu()
Return whether the BdaBridge will use the GPU or not return whether the BdaBridge will use the GPU or...
Definition: BdaBridge.hpp:81
bool getUseFpga()
Return whether the BdaBridge will use the FPGA or not return whether the BdaBridge will use the FPGA ...
Definition: BdaBridge.hpp:92
std::string getAccleratorName()
Return the selected accelerator mode, this is input via the command-line.
Definition: BdaBridge.hpp:97
void initWellContributions(WellContributions &wellContribs)
Initialize the WellContributions object with opencl context and queue those must be set before callin...
Definition: BdaBridge.cpp:269
This class serves to eliminate the need to include the WellContributions into the matrix (with –matri...
Definition: WellContributions.hpp:61
This file contains a set of helper functions used by VFPProd / VFPInj.
Definition: BlackoilPhases.hpp:26