20 #ifndef WELLCONTRIBUTIONS_HEADER_INCLUDED
21 #define WELLCONTRIBUTIONS_HEADER_INCLUDED
24 #include <cuda_runtime.h>
28 #include <opm/simulators/linalg/bda/opencl.hpp>
29 #include <opm/simulators/linalg/bda/openclKernels.hpp>
34 #include <opm/simulators/linalg/bda/MultisegmentWellContribution.hpp>
35 #if HAVE_SUITESPARSE_UMFPACK
38 #include <dune/common/version.hh>
63 #if DUNE_VERSION_NEWER(DUNE_ISTL, 2, 7)
64 using UMFPackIndex = SuiteSparse_long;
66 using UMFPackIndex = int;
76 bool opencl_gpu =
false;
77 bool cuda_gpu =
false;
78 bool allocated =
false;
82 unsigned int dim_wells;
83 unsigned int num_blocks = 0;
84 unsigned int num_std_wells = 0;
85 unsigned int num_ms_wells = 0;
86 unsigned int num_blocks_so_far = 0;
87 unsigned int num_std_wells_so_far = 0;
88 unsigned int *val_pointers =
nullptr;
90 double *h_x =
nullptr;
91 double *h_y =
nullptr;
92 std::vector<MultisegmentWellContribution*> multisegments;
96 cl::CommandQueue *queue;
97 bda::stdwell_apply_kernel_type *kernel;
98 bda::stdwell_apply_no_reorder_kernel_type *kernel_no_reorder;
99 std::vector<cl::Event> events;
101 std::unique_ptr<cl::Buffer> d_Cnnzs_ocl, d_Dnnzs_ocl, d_Bnnzs_ocl;
102 std::unique_ptr<cl::Buffer> d_Ccols_ocl, d_Bcols_ocl;
103 std::unique_ptr<cl::Buffer> d_val_pointers_ocl;
105 bool reorder =
false;
106 int *h_toOrder =
nullptr;
113 double *d_Cnnzs =
nullptr;
114 double *d_Dnnzs =
nullptr;
115 double *d_Bnnzs =
nullptr;
116 int *d_Ccols =
nullptr;
117 int *d_Bcols =
nullptr;
118 double *d_z1 =
nullptr;
119 double *d_z2 =
nullptr;
120 unsigned int *d_val_pointers =
nullptr;
123 void allocStandardWells();
126 void freeCudaMemory();
133 void addMatrixGpu(
MatrixType type,
int *colIndices,
double *values,
unsigned int val_size);
140 void setCudaStream(cudaStream_t stream);
146 void apply(
double *d_x,
double *d_y);
150 void setKernel(bda::stdwell_apply_kernel_type *kernel_,
151 bda::stdwell_apply_no_reorder_kernel_type *kernel_no_reorder_);
152 void setOpenCLEnv(cl::Context *context_, cl::CommandQueue *queue_);
158 void setReordering(
int *toOrder,
bool reorder);
159 void apply_stdwells(cl::Buffer d_x, cl::Buffer d_y, cl::Buffer d_toOrder);
160 void apply_mswells(cl::Buffer d_x, cl::Buffer d_y);
161 void apply(cl::Buffer d_x, cl::Buffer d_y, cl::Buffer d_toOrder);
164 unsigned int getNumWells(){
165 return num_std_wells + num_ms_wells;
186 void setBlockSize(
unsigned int dim,
unsigned int dim_wells);
210 std::vector<double> &Bvalues, std::vector<unsigned int> &BcolIndices, std::vector<unsigned int> &BrowPointers,
211 unsigned int DnumBlocks,
double *Dvalues,
212 UMFPackIndex *DcolPointers, UMFPackIndex *DrowIndices,
213 std::vector<double> &Cvalues);
This class serves to eliminate the need to include the WellContributions into the matrix (with –matri...
Definition: WellContributions.hpp:61
WellContributions(std::string accelerator_mode, bool useWellConn)
Create a new WellContributions.
Definition: WellContributions.cpp:31
~WellContributions()
Destroy a WellContributions, and free memory.
Definition: WellContributions.cpp:51
void setBlockSize(unsigned int dim, unsigned int dim_wells)
Indicate how large the blocks of the StandardWell (C and B) are.
Definition: WellContributions.cpp:215
void addMatrix(MatrixType type, int *colIndices, double *values, unsigned int val_size)
Store a matrix in this object, in blocked csr format, can only be called after alloc() is called.
Definition: WellContributions.cpp:152
void addMultisegmentWellContribution(unsigned int dim, unsigned int dim_wells, unsigned int Mb, std::vector< double > &Bvalues, std::vector< unsigned int > &BcolIndices, std::vector< unsigned int > &BrowPointers, unsigned int DnumBlocks, double *Dvalues, UMFPackIndex *DcolPointers, UMFPackIndex *DrowIndices, std::vector< double > &Cvalues)
Add a MultisegmentWellContribution, actually creates an object on heap that is destroyed in the destr...
Definition: WellContributions.cpp:261
MatrixType
StandardWell has C, D and B matrices that need to be copied.
Definition: WellContributions.hpp:69
void alloc()
Allocate memory for the StandardWells.
Definition: WellContributions.cpp:236
void addNumBlocks(unsigned int numBlocks)
Indicate how large the next StandardWell is, this function cannot be called after alloc() is called.
Definition: WellContributions.cpp:227
This file contains a set of helper functions used by VFPProd / VFPInj.
Definition: BlackoilPhases.hpp:26