20#ifndef BDABRIDGE_HEADER_INCLUDED
21#define BDABRIDGE_HEADER_INCLUDED
23#include "dune/istl/solver.hh"
25#include <opm/simulators/linalg/bda/BdaSolver.hpp>
26#include <opm/simulators/linalg/bda/BlockedMatrix.hpp>
27#include <opm/simulators/linalg/bda/ILUReorder.hpp>
32class WellContributions;
34typedef Dune::InverseOperatorResult InverseOperatorResult;
35using Opm::Accelerator::ILUReorder;
38template <
class Br
idgeMatrix,
class Br
idgeVector,
int block_size>
44 bool use_fpga =
false;
45 std::string accelerator_mode;
46 std::unique_ptr<Opm::Accelerator::BdaSolver<block_size> > backend;
47 std::shared_ptr<Opm::Accelerator::BlockedMatrix> matrix;
48 std::shared_ptr<Opm::Accelerator::BlockedMatrix> jacMatrix;
49 std::vector<int> h_rows, h_cols;
50 std::vector<int> h_jacRows, h_jacCols;
51 std::vector<typename BridgeMatrix::size_type> diagIndices;
52 std::vector<typename BridgeMatrix::size_type> jacDiagIndices;
65 BdaBridge(std::string accelerator_mode, std::string fpga_bitstream,
int linear_solver_verbosity,
int maxit,
double tolerance,
66 unsigned int platformID,
unsigned int deviceID, std::string opencl_ilu_reorder, std::string linsolver);
77 void solve_system(BridgeMatrix *bridgeMat, BridgeMatrix *jacMat,
int numJacobiBlocks, BridgeVector &b,
WellContributions& wellContribs, InverseOperatorResult &result);
109 return accelerator_mode;
BdaBridge acts as interface between opm-simulators with the BdaSolvers.
Definition: BdaBridge.hpp:40
void get_result(BridgeVector &x)
Get the resulting x vector.
Definition: BdaBridge.cpp:303
bool getUseGpu()
Return whether the BdaBridge will use the GPU or not return whether the BdaBridge will use the GPU or...
Definition: BdaBridge.hpp:85
void initWellContributions(WellContributions &wellContribs, unsigned N)
Initialize the WellContributions object with opencl context and queue those must be set before callin...
Definition: BdaBridge.cpp:310
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:59
void solve_system(BridgeMatrix *bridgeMat, BridgeMatrix *jacMat, int numJacobiBlocks, BridgeVector &b, WellContributions &wellContribs, InverseOperatorResult &result)
Solve linear system, A*x = b.
Definition: BdaBridge.cpp:216
static void copySparsityPatternFromISTL(const BridgeMatrix &mat, std::vector< int > &h_rows, std::vector< int > &h_cols)
Store sparsity pattern into vectors.
Definition: BdaBridge.cpp:171
bool getUseFpga()
Return whether the BdaBridge will use the FPGA or not return whether the BdaBridge will use the FPGA ...
Definition: BdaBridge.hpp:103
std::string getAccleratorName()
Return the selected accelerator mode, this is input via the command-line.
Definition: BdaBridge.hpp:108
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