My Project
WellContributions.hpp
1 /*
2  Copyright 2020 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 WELLCONTRIBUTIONS_HEADER_INCLUDED
21 #define WELLCONTRIBUTIONS_HEADER_INCLUDED
22 
23 #if HAVE_CUDA
24 #include <cuda_runtime.h>
25 #endif
26 
27 #if HAVE_OPENCL
28 #include <opm/simulators/linalg/bda/opencl.hpp>
29 #include <opm/simulators/linalg/bda/openclKernels.hpp>
30 #endif
31 
32 #include <vector>
33 
34 #include <opm/simulators/linalg/bda/MultisegmentWellContribution.hpp>
35 #if HAVE_SUITESPARSE_UMFPACK
36 #include<umfpack.h>
37 #endif
38 #include <dune/common/version.hh>
39 
40 namespace Opm
41 {
42 
49 
61 {
62 public:
63 #if DUNE_VERSION_NEWER(DUNE_ISTL, 2, 7)
64  using UMFPackIndex = SuiteSparse_long;
65 #else
66  using UMFPackIndex = int;
67 #endif
69  enum class MatrixType {
70  C,
71  D,
72  B
73  };
74 
75 private:
76  bool opencl_gpu = false;
77  bool cuda_gpu = false;
78  bool allocated = false;
79 
80  unsigned int N; // number of rows (not blockrows) in vectors x and y
81  unsigned int dim; // number of columns in blocks in B and C, equal to StandardWell::numEq
82  unsigned int dim_wells; // number of rows in blocks in B and C, equal to StandardWell::numStaticWellEq
83  unsigned int num_blocks = 0; // total number of blocks in all wells
84  unsigned int num_std_wells = 0; // number of StandardWells in this object
85  unsigned int num_ms_wells = 0; // number of MultisegmentWells in this object, must equal multisegments.size()
86  unsigned int num_blocks_so_far = 0; // keep track of where next data is written
87  unsigned int num_std_wells_so_far = 0; // keep track of where next data is written
88  unsigned int *val_pointers = nullptr; // val_pointers[wellID] == index of first block for this well in Ccols and Bcols
89 
90  double *h_x = nullptr;
91  double *h_y = nullptr;
92  std::vector<MultisegmentWellContribution*> multisegments;
93 
94 #if HAVE_OPENCL
95  cl::Context *context;
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;
100 
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;
104 
105  bool reorder = false;
106  int *h_toOrder = nullptr;
107 #endif
108 
109 #if HAVE_CUDA
110  cudaStream_t stream;
111 
112  // data for StandardWells, could remain nullptrs if not used
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;
121 
123  void allocStandardWells();
124 
126  void freeCudaMemory();
127 
133  void addMatrixGpu(MatrixType type, int *colIndices, double *values, unsigned int val_size);
134 #endif
135 
136 public:
137 #if HAVE_CUDA
140  void setCudaStream(cudaStream_t stream);
141 
146  void apply(double *d_x, double *d_y);
147 #endif
148 
149 #if HAVE_OPENCL
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_);
153 
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);
162 #endif
163 
164  unsigned int getNumWells(){
165  return num_std_wells + num_ms_wells;
166  }
167 
170  void addNumBlocks(unsigned int numBlocks);
171 
173  void alloc();
174 
178  WellContributions(std::string accelerator_mode, bool useWellConn);
179 
182 
186  void setBlockSize(unsigned int dim, unsigned int dim_wells);
187 
193  void addMatrix(MatrixType type, int *colIndices, double *values, unsigned int val_size);
194 
208  void addMultisegmentWellContribution(unsigned int dim, unsigned int dim_wells,
209  unsigned int Mb,
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);
214 };
215 } //namespace Opm
216 
217 #endif
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