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/BlockedMatrix.hpp>
27#include <opm/simulators/linalg/bda/ILUReorder.hpp>
28
29namespace Opm
30{
31
32class WellContributions;
33
34typedef Dune::InverseOperatorResult InverseOperatorResult;
35using Opm::Accelerator::ILUReorder;
36
38template <class BridgeMatrix, class BridgeVector, int block_size>
40{
41private:
42 int verbosity = 0;
43 bool use_gpu = false;
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; // 'stores' matrix, actually points to h_rows, h_cols and the received BridgeMatrix for the nonzeroes
48 std::shared_ptr<Opm::Accelerator::BlockedMatrix> jacMatrix; // 'stores' preconditioner matrix, actually points to h_rows, h_cols and the received BridgeMatrix for the nonzeroes
49 std::vector<int> h_rows, h_cols; // store the sparsity pattern of the matrix
50 std::vector<int> h_jacRows, h_jacCols; // store the sparsity pattern of the jacMatrix
51 std::vector<typename BridgeMatrix::size_type> diagIndices; // contains offsets of the diagonal blocks wrt start of the row, used for replaceZeroDiagonal()
52 std::vector<typename BridgeMatrix::size_type> jacDiagIndices; // same but for jacMatrix
53
54public:
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);
67
68
77 void solve_system(BridgeMatrix *bridgeMat, BridgeMatrix *jacMat, int numJacobiBlocks, BridgeVector &b, WellContributions& wellContribs, InverseOperatorResult &result);
78
81 void get_result(BridgeVector &x);
82
85 bool getUseGpu(){
86 return use_gpu;
87 }
88
93 static void copySparsityPatternFromISTL(const BridgeMatrix& mat, std::vector<int>& h_rows, std::vector<int>& h_cols);
94
99 void initWellContributions(WellContributions& wellContribs, unsigned N);
100
104 return use_fpga;
105 }
106
108 std::string getAccleratorName(){
109 return accelerator_mode;
110 }
111
112}; // end class BdaBridge
113
114}
115
116#endif
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