Opm::BdaBridge< BridgeMatrix, BridgeVector, block_size > Class Template Reference

BdaBridge acts as interface between opm-simulators with the BdaSolvers. More...

#include <BdaBridge.hpp>

Public Member Functions

 BdaBridge (std::string accelerator_mode, int linear_solver_verbosity, int maxit, double tolerance, unsigned int platformID, unsigned int deviceID, bool opencl_ilu_parallel, std::string linsolver)
 
void solve_system (BridgeMatrix *bridgeMat, BridgeMatrix *jacMat, int numJacobiBlocks, BridgeVector &b, WellContributions &wellContribs, InverseOperatorResult &result)
 
void get_result (BridgeVector &x)
 
bool getUseGpu ()
 
void initWellContributions (WellContributions &wellContribs, unsigned N)
 
std::string getAccleratorName ()
 Return the selected accelerator mode, this is input via the command-line. More...
 

Static Public Member Functions

static void copySparsityPatternFromISTL (const BridgeMatrix &mat, std::vector< int > &h_rows, std::vector< int > &h_cols)
 

Detailed Description

template<class BridgeMatrix, class BridgeVector, int block_size>
class Opm::BdaBridge< BridgeMatrix, BridgeVector, block_size >

BdaBridge acts as interface between opm-simulators with the BdaSolvers.

Constructor & Destructor Documentation

◆ BdaBridge()

template<class BridgeMatrix , class BridgeVector , int block_size>
Opm::BdaBridge< BridgeMatrix, BridgeVector, block_size >::BdaBridge ( std::string  accelerator_mode,
int  linear_solver_verbosity,
int  maxit,
double  tolerance,
unsigned int  platformID,
unsigned int  deviceID,
bool  opencl_ilu_parallel,
std::string  linsolver 
)

Construct a BdaBridge

Parameters
[in]accelerator_modeto select if an accelerated solver is used, is passed via command-line: '–accelerator-mode=[none|cusparse|opencl|amgcl|rocalution|rocsparse]'
[in]linear_solver_verbosityverbosity of BdaSolver
[in]maxitmaximum number of iterations for BdaSolver
[in]tolerancerequired relative tolerance for BdaSolver
[in]platformIDthe OpenCL platform ID to be used
[in]deviceIDthe device ID to be used by the cusparse- and openclSolvers, too high values could cause runtime errors
[in]opencl_ilu_parallelwhether to parallelize the ILU decomposition and application in OpenCL with level_scheduling
[in]linsolverindicating the preconditioner, equal to the –linear-solver cmdline argument

Member Function Documentation

◆ copySparsityPatternFromISTL()

template<class BridgeMatrix , class BridgeVector , int block_size>
static void Opm::BdaBridge< BridgeMatrix, BridgeVector, block_size >::copySparsityPatternFromISTL ( const BridgeMatrix &  mat,
std::vector< int > &  h_rows,
std::vector< int > &  h_cols 
)
static

Store sparsity pattern into vectors

Parameters
[in]matinput matrix, probably BCRSMatrix
[out]h_rowsrowpointers
[out]h_colscolumnindices

◆ get_result()

template<class BridgeMatrix , class BridgeVector , int block_size>
void Opm::BdaBridge< BridgeMatrix, BridgeVector, block_size >::get_result ( BridgeVector &  x)

Get the resulting x vector

Parameters
[in,out]xvector x, should be of type Dune::BlockVector

◆ getAccleratorName()

template<class BridgeMatrix , class BridgeVector , int block_size>
std::string Opm::BdaBridge< BridgeMatrix, BridgeVector, block_size >::getAccleratorName ( )
inline

Return the selected accelerator mode, this is input via the command-line.

◆ getUseGpu()

template<class BridgeMatrix , class BridgeVector , int block_size>
bool Opm::BdaBridge< BridgeMatrix, BridgeVector, block_size >::getUseGpu ( )
inline

Return whether the BdaBridge will use the GPU or not return whether the BdaBridge will use the GPU or not

◆ initWellContributions()

template<class BridgeMatrix , class BridgeVector , int block_size>
void Opm::BdaBridge< BridgeMatrix, BridgeVector, block_size >::initWellContributions ( WellContributions wellContribs,
unsigned  N 
)

Initialize the WellContributions object with opencl context and queue those must be set before calling BlackOilWellModel::getWellContributions() in ISTL

Parameters
[in]wellContribscontainer to hold all WellContributions
[in]Nnumber of rows in scalar vector that wellContribs will be applied on

◆ solve_system()

template<class BridgeMatrix , class BridgeVector , int block_size>
void Opm::BdaBridge< BridgeMatrix, BridgeVector, block_size >::solve_system ( BridgeMatrix *  bridgeMat,
BridgeMatrix *  jacMat,
int  numJacobiBlocks,
BridgeVector &  b,
WellContributions wellContribs,
InverseOperatorResult result 
)

Solve linear system, A*x = b

Warning
Values of A might get overwritten!
Parameters
[in]bridgeMatmatrix A, should be of type Dune::BCRSMatrix
[in]jacMatmatrix A, but modified for the preconditioner, should be of type Dune::BCRSMatrix
[in]numJacobiBlocksnumber of jacobi blocks in jacMat
[in]bvector b, should be of type Dune::BlockVector
[in]wellContribscontains all WellContributions, to apply them separately, instead of adding them to matrix A
[in,out]resultsummary of solver result

The documentation for this class was generated from the following file: