My Project
Loading...
Searching...
No Matches
GpuBridge.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 GPUBRIDGE_HEADER_INCLUDED
21#define GPUBRIDGE_HEADER_INCLUDED
22
23#include "dune/istl/solver.hh" // for struct InverseOperatorResult
24
25#include <opm/simulators/linalg/gpubridge/GpuSolver.hpp>
26
27namespace Opm
28{
29
30template<class Scalar> class WellContributions;
31
32typedef Dune::InverseOperatorResult InverseOperatorResult;
33
35template <class BridgeMatrix, class BridgeVector, int block_size>
37{
38private:
39 using Scalar = typename BridgeVector::field_type;
40 int verbosity = 0;
41 bool use_gpu = false;
42 std::string accelerator_mode;
43 std::unique_ptr<Accelerator::GpuSolver<Scalar,block_size>> backend;
44 std::shared_ptr<Accelerator::BlockedMatrix<Scalar>> matrix; // 'stores' matrix, actually points to h_rows, h_cols and the received BridgeMatrix for the nonzeroes
45 std::shared_ptr<Accelerator::BlockedMatrix<Scalar>> jacMatrix; // 'stores' preconditioner matrix, actually points to h_rows, h_cols and the received BridgeMatrix for the nonzeroes
46 std::vector<int> h_rows, h_cols; // store the sparsity pattern of the matrix
47 std::vector<int> h_jacRows, h_jacCols; // store the sparsity pattern of the jacMatrix
48 std::vector<typename BridgeMatrix::size_type> diagIndices; // contains offsets of the diagonal blocks wrt start of the row, used for replaceZeroDiagonal()
49 std::vector<typename BridgeMatrix::size_type> jacDiagIndices; // same but for jacMatrix
50
51public:
61 GpuBridge(std::string accelerator_mode,
63 int maxit,
64 Scalar tolerance,
65 unsigned int platformID,
66 unsigned int deviceID,
67 bool opencl_ilu_parallel,
68 std::string linsolver);
69
70
80 BridgeMatrix* jacMat,
81 int numJacobiBlocks,
84 InverseOperatorResult &result);
85
88 void get_result(BridgeVector &x);
89
92 bool getUseGpu()
93 {
94 return use_gpu;
95 }
96
101 static void copySparsityPatternFromISTL(const BridgeMatrix& mat,
102 std::vector<int>& h_rows,
103 std::vector<int>& h_cols);
104
110
112 std::string getAccleratorName()
113 {
114 return accelerator_mode;
115 }
116}; // end class GpuBridge
117
118}
119
120#endif
GpuBridge acts as interface between opm-simulators with the GpuSolvers.
Definition GpuBridge.hpp:37
std::string getAccleratorName()
Return the selected accelerator mode, this is input via the command-line.
Definition GpuBridge.hpp:112
void initWellContributions(WellContributions< Scalar > &wellContribs, unsigned N)
Initialize the WellContributions object with opencl context and queue those must be set before callin...
Definition GpuBridge.cpp:341
void get_result(BridgeVector &x)
Get the resulting x vector.
Definition GpuBridge.cpp:332
void solve_system(BridgeMatrix *bridgeMat, BridgeMatrix *jacMat, int numJacobiBlocks, BridgeVector &b, WellContributions< Scalar > &wellContribs, InverseOperatorResult &result)
Solve linear system, A*x = b.
Definition GpuBridge.cpp:229
static void copySparsityPatternFromISTL(const BridgeMatrix &mat, std::vector< int > &h_rows, std::vector< int > &h_cols)
Store sparsity pattern into vectors.
Definition GpuBridge.cpp:180
bool getUseGpu()
Return whether the GpuBridge will use the GPU or not return whether the GpuBridge will use the GPU or...
Definition GpuBridge.hpp:92
This class serves to eliminate the need to include the WellContributions into the matrix (with –matri...
Definition WellContributions.hpp:51
This file contains a set of helper functions used by VFPProd / VFPInj.
Definition blackoilboundaryratevector.hh:37
constexpr auto getPropValue()
get the value data member of a property
Definition propertysystem.hh:242