CpGridVanguard.hpp
Go to the documentation of this file.
1// -*- mode: C++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*-
2// vi: set et ts=4 sw=4 sts=4:
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 2 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 Consult the COPYING file in the top-level source directory of this
20 module for the precise wording of the license and the list of
21 copyright holders.
22*/
27#ifndef OPM_CPGRID_VANGUARD_HPP
28#define OPM_CPGRID_VANGUARD_HPP
29
30#include <opm/common/TimingMacros.hpp>
31
34
39
40#include <array>
41#include <functional>
42#include <memory>
43#include <stdexcept>
44#include <tuple>
45#include <vector>
46
47namespace Opm {
48template <class TypeTag>
49class CpGridVanguard;
50}
51
52namespace Opm::Properties {
53
54namespace TTag {
56 using InheritsFrom = std::tuple<FlowBaseVanguard>;
57};
58}
59
60// declare the properties
61template<class TypeTag>
62struct Vanguard<TypeTag, TTag::CpGridVanguard> {
64};
65template<class TypeTag>
66struct Grid<TypeTag, TTag::CpGridVanguard> {
67 using type = Dune::CpGrid;
68};
69template<class TypeTag>
70struct EquilGrid<TypeTag, TTag::CpGridVanguard> {
72};
73
74} // namespace Opm::Properties
75
76namespace Opm {
77
85template <class TypeTag>
86class CpGridVanguard : public FlowBaseVanguard<TypeTag>
87 , public GenericCpGridVanguard<GetPropType<TypeTag, Properties::ElementMapper>,
88 GetPropType<TypeTag, Properties::GridView>,
89 GetPropType<TypeTag, Properties::Scalar>>
90{
91 friend class FlowBaseVanguard<TypeTag>;
93
97
98public:
104 static constexpr int dimensionworld = Grid::dimensionworld;
106 static constexpr bool waterEnabled = Indices::waterEnabled;
107 static constexpr bool gasEnabled = Indices::gasEnabled;
108 static constexpr bool oilEnabled = Indices::oilEnabled;
109private:
110 using Element = typename GridView::template Codim<0>::Entity;
111
112public:
113 explicit CpGridVanguard(Simulator& simulator)
114 : FlowBaseVanguard<TypeTag>(simulator)
115 {
116 this->checkConsistency();
118 }
119
120 int compressedIndexForInteriorLGR(const std::string& lgr_tag, const Connection& conn) const override
121 {
122 const std::array<int,3> lgr_ijk = {conn.getI(), conn.getJ(), conn.getK()};
123 const auto& lgr_level = this->grid().getLgrNameToLevel().at(lgr_tag);
124 if (ParentType::lgrMappers_.has_value() == false) {
125 ParentType::lgrMappers_.emplace(this->grid().mapLocalCartesianIndexSetsToLeafIndexSet());
126 }
127 const auto& lgr_dim = this->grid().currentData()[lgr_level]->logicalCartesianSize();
128 const auto lgr_cartesian_index = (lgr_ijk[2]*lgr_dim[0]*lgr_dim[1]) + (lgr_ijk[1]*lgr_dim[0]) + (lgr_ijk[0]);
129 return ParentType::lgrMappers_.value()[lgr_level].at(lgr_cartesian_index);
130 }
135 {
136 const auto& runspec = this->eclState().runspec();
137 const auto& config = this->eclState().getSimulationConfig();
138 const auto& phases = runspec.phases();
139
140 // check for correct module setup
141 if (config.isThermal()) {
142 if (getPropValue<TypeTag, Properties::EnableEnergy>() == false) {
143 throw std::runtime_error("Input specifies energy while simulator has disabled it, try xxx_energy");
144 }
145 } else {
146 if (getPropValue<TypeTag, Properties::EnableEnergy>() == true) {
147 throw std::runtime_error("Input specifies no energy while simulator has energy, try run without _energy");
148 }
149 }
150
151 if (config.isDiffusive()) {
152 if (getPropValue<TypeTag, Properties::EnableDiffusion>() == false) {
153 throw std::runtime_error("Input specifies diffusion while simulator has disabled it, try xxx_diffusion");
154 }
155 }
156
157 if (runspec.micp()) {
158 if (getPropValue<TypeTag, Properties::EnableBioeffects>() == false) {
159 throw std::runtime_error("Input specifies MICP while simulator has it disabled");
160 }
161 }
162
163 if (runspec.biof()) {
164 if (getPropValue<TypeTag, Properties::EnableBioeffects>() == false) {
165 throw std::runtime_error("Input specifies Biofilm while simulator has it disabled");
166 }
167 }
168
169 if (phases.active(Phase::BRINE)) {
170 if (getPropValue<TypeTag, Properties::EnableBrine>() == false) {
171 throw std::runtime_error("Input specifies Brine while simulator has it disabled");
172 }
173 }
174
175 if (phases.active(Phase::POLYMER)) {
176 if (getPropValue<TypeTag, Properties::EnablePolymer>() == false) {
177 throw std::runtime_error("Input specifies Polymer while simulator has it disabled");
178 }
179 }
180
181 // checking for correct phases is more difficult TODO!
182 if (phases.active(Phase::ZFRACTION)) {
183 if (getPropValue<TypeTag, Properties::EnableExtbo>() == false) {
184 throw std::runtime_error("Input specifies ExBo while simulator has it disabled");
185 }
186 }
187 if (phases.active(Phase::FOAM)) {
188 if (getPropValue<TypeTag, Properties::EnableFoam>() == false) {
189 throw std::runtime_error("Input specifies Foam while simulator has it disabled");
190 }
191 }
192
193 if (phases.active(Phase::SOLVENT)) {
194 if (getPropValue<TypeTag, Properties::EnableSolvent>() == false) {
195 throw std::runtime_error("Input specifies Solvent while simulator has it disabled");
196 }
197 }
198 if(phases.active(Phase::WATER)){
199 if(waterEnabled == false){
200 throw std::runtime_error("Input specifies water while simulator has it disabled");
201 }
202 }
203 if(phases.active(Phase::GAS)){
204 if(gasEnabled == false){
205 throw std::runtime_error("Input specifies gas while simulator has it disabled");
206 }
207 }
208 if(phases.active(Phase::OIL)){
209 if(oilEnabled == false){
210 throw std::runtime_error("Input specifies oil while simulator has it disabled");
211 }
212 }
213
214 }
215
222 {
223 globalTrans_.reset();
224 }
225
227 {
228 assert( globalTrans_ != nullptr );
229 return *globalTrans_;
230 }
231
238 {
239#if HAVE_MPI
240 if (const auto& extPFile = this->externalPartitionFile();
241 !extPFile.empty() && (extPFile != "none"))
242 {
244 }
245
246 this->doLoadBalance_(this->edgeWeightsMethod(), this->ownersFirst(),
247 this->addCorners(), this->numOverlap(),
248 this->partitionMethod(), this->serialPartitioning(),
251 this->imbalanceTol(),
252 this->gridView(), this->schedule(),
253 this->eclState(), this->parallelWells_,
254 this->numJacobiBlocks(), this->enableEclOutput());
255#endif
256
257 this->updateGridView_();
259 this->updateCellDepths_();
260 this->updateCellThickness_();
261
262#if HAVE_MPI
263 this->distributeFieldProps_(this->eclState());
264#endif
265 }
266
270 void addLgrs()
271 {
272 // Check if input file contains Lgrs. Add them, if any.
273 // In a parallel run, this adds the LGRs on the distributed simulation grid.
274 if (const auto& lgrs = this->eclState().getLgrs(); lgrs.size() > 0) {
275 OpmLog::info("\nAdding LGRs to the grid and updating its leaf grid view");
276 this->addLgrsUpdateLeafView(lgrs, lgrs.size(), *this->grid_);
277
278 this->updateGridView_();
279 this->updateCellDepths_();
280 this->updateCellThickness_();
281
282 if (this->grid_->comm().size()>1) {
283 // Add LGRs and update the leaf grid view in the global (undistributed) simulation grid.
284 // Purpose: To enable synchronization of cell ids in 'serial mode',
285 // we rely on the "parent-to-children" cell id mapping.
286 OpmLog::info("\nAdding LGRs to the global view and updating its leaf grid view");
287 this->grid_->switchToGlobalView();
288 this->addLgrsUpdateLeafView(lgrs, lgrs.size(), *this->grid_);
289 this->grid_->switchToDistributedView();
290 this->grid_->syncDistributedGlobalCellIds();
291 }
292 }
293 }
294
295 unsigned int gridEquilIdxToGridIdx(unsigned int elemIndex) const {
296 return elemIndex;
297 }
298
299 unsigned int gridIdxToEquilGridIdx(unsigned int elemIndex) const {
300 return elemIndex;
301 }
309 std::function<std::array<double,dimensionworld>(int)>
311 {
312 return this->cellCentroids_(this->cartesianIndexMapper(), true);
313 }
314
315 const std::vector<int>& globalCell()
316 {
317 return this->grid().globalCell();
318 }
319
320protected:
322 {
323 this->doCreateGrids_(this->edgeConformal(), this->eclState());
324 }
325
326 void allocTrans() override
327 {
328 OPM_TIMEBLOCK(allocateTrans);
330 this->gridView(),
331 this->cartesianIndexMapper(),
332 this->grid(),
333 this->cellCentroids(),
334 getPropValue<TypeTag, Properties::EnableEnergy>(),
335 getPropValue<TypeTag, Properties::EnableDiffusion>(),
336 getPropValue<TypeTag, Properties::EnableDispersion>()));
338 }
339
340 double getTransmissibility(unsigned I, unsigned J) const override
341 {
342 return globalTrans_->transmissibility(I,J);
343 }
344
345#if HAVE_MPI
346 const std::string& zoltanParams() const override
347 {
348 return this->zoltanParams_;
349 }
350
351 double zoltanPhgEdgeSizeThreshold() const override
352 {
353 return this->zoltanPhgEdgeSizeThreshold_;
354 }
355
356 const std::string& metisParams() const override
357 {
358 return this->metisParams_;
359 }
360#endif
361
362 // removing some connection located in inactive grid cells
364 {
365 this->doFilterConnections_(this->schedule());
366 }
367
368 // \Note: this globalTrans_ is used for domain decomposition and INIT file output.
369 // It only contains trans_ due to permeability and does not contain thermalHalfTrans_,
370 // diffusivity_ abd dispersivity_. The main reason is to reduce the memory usage for rank 0
371 // during parallel running.
372 std::unique_ptr<TransmissibilityType> globalTrans_;
373};
374
375} // namespace Opm
376
377#endif // OPM_CPGRID_VANGUARD_HPP
This file ensures that flow can be compiled in the presence of dune-fem.
Declares the properties required by the black oil model.
Definition: CollectDataOnIORank.hpp:49
void updateGridView_()
Definition: basevanguard.hh:132
const GridView & gridView() const
Returns a reference to the grid view to be used.
Definition: basevanguard.hh:70
Helper class for grid instantiation of ECL file-format using problems.
Definition: CpGridVanguard.hpp:90
GetPropType< TypeTag, Properties::Indices > Indices
Definition: CpGridVanguard.hpp:105
void releaseGlobalTransmissibilities()
Free the memory occupied by the global transmissibility object.
Definition: CpGridVanguard.hpp:221
void addLgrs()
Add LGRs and update Leaf Grid View in the simulation grid.
Definition: CpGridVanguard.hpp:270
const std::string & zoltanParams() const override
Definition: CpGridVanguard.hpp:346
double zoltanPhgEdgeSizeThreshold() const override
Definition: CpGridVanguard.hpp:351
void filterConnections_()
Definition: CpGridVanguard.hpp:363
const TransmissibilityType & globalTransmissibility() const
Definition: CpGridVanguard.hpp:226
void allocTrans() override
Definition: CpGridVanguard.hpp:326
static constexpr int dimensionworld
Definition: CpGridVanguard.hpp:104
static constexpr bool waterEnabled
Definition: CpGridVanguard.hpp:106
CpGridVanguard(Simulator &simulator)
Definition: CpGridVanguard.hpp:113
void createGrids_()
Definition: CpGridVanguard.hpp:321
unsigned int gridIdxToEquilGridIdx(unsigned int elemIndex) const
Definition: CpGridVanguard.hpp:299
GetPropType< TypeTag, Properties::EquilGrid > EquilGrid
Definition: CpGridVanguard.hpp:101
static constexpr bool gasEnabled
Definition: CpGridVanguard.hpp:107
int compressedIndexForInteriorLGR(const std::string &lgr_tag, const Connection &conn) const override
Definition: CpGridVanguard.hpp:120
void checkConsistency()
Definition: CpGridVanguard.hpp:134
unsigned int gridEquilIdxToGridIdx(unsigned int elemIndex) const
Definition: CpGridVanguard.hpp:295
const std::string & metisParams() const override
Definition: CpGridVanguard.hpp:356
std::function< std::array< double, dimensionworld >(int)> cellCentroids() const
Get function to query cell centroids for a distributed grid.
Definition: CpGridVanguard.hpp:310
std::unique_ptr< TransmissibilityType > globalTrans_
Definition: CpGridVanguard.hpp:372
static constexpr bool oilEnabled
Definition: CpGridVanguard.hpp:108
void loadBalance()
Distribute the simulation grid over multiple processes.
Definition: CpGridVanguard.hpp:237
Transmissibility< Grid, GridView, ElementMapper, CartesianIndexMapper, Scalar > TransmissibilityType
Definition: CpGridVanguard.hpp:103
const std::vector< int > & globalCell()
Definition: CpGridVanguard.hpp:315
double getTransmissibility(unsigned I, unsigned J) const override
Definition: CpGridVanguard.hpp:340
Helper class for grid instantiation of ECL file-format using problems.
Definition: FlowBaseVanguard.hpp:84
void updateCartesianToCompressedMapping_()
Definition: FlowBaseVanguard.hpp:334
void updateCellThickness_()
Definition: FlowBaseVanguard.hpp:379
void updateCellDepths_()
Definition: FlowBaseVanguard.hpp:356
void callImplementationInit()
Definition: FlowBaseVanguard.hpp:322
std::optional< std::vector< std::unordered_map< std::size_t, std::size_t > > > lgrMappers_
Mapping between LGR cartesian and compressed cells. It is initialized as it is called.
Definition: FlowBaseVanguard.hpp:448
std::function< std::array< double, dimensionworld >(int)> cellCentroids_(const CartMapper &cartMapper, const bool &isCpGrid) const
Get function to query cell centroids for a distributed grid.
Definition: FlowBaseVanguard.hpp:302
int numJacobiBlocks() const
Number of blocks in the Block-Jacobi preconditioner.
Definition: FlowGenericVanguard.hpp:248
const Schedule & schedule() const
Return a reference to the object that managages the ECL schedule.
Definition: FlowGenericVanguard.hpp:177
std::string zoltanParams_
Definition: FlowGenericVanguard.hpp:392
bool addCorners() const
Definition: FlowGenericVanguard.hpp:267
bool serialPartitioning() const
Parameter that decides if partitioning for parallel runs should be performed on a single process only...
Definition: FlowGenericVanguard.hpp:283
ParallelWellStruct parallelWells_
Information about wells in parallel.
Definition: FlowGenericVanguard.hpp:431
bool enableDistributedWells() const
Whether perforations of a well might be distributed.
Definition: FlowGenericVanguard.hpp:311
const EclipseState & eclState() const
Return a reference to the internalized ECL deck.
Definition: FlowGenericVanguard.hpp:168
double imbalanceTol() const
Parameter that sets the imbalance tolarance, depending on the chosen partition method.
Definition: FlowGenericVanguard.hpp:290
double zoltanPhgEdgeSizeThreshold_
Definition: FlowGenericVanguard.hpp:391
const std::string & externalPartitionFile() const
Definition: FlowGenericVanguard.hpp:302
bool enableEclOutput() const
Whether or not to emit result files that are compatible with a commercial reservoir simulator.
Definition: FlowGenericVanguard.hpp:318
bool ownersFirst() const
Parameter that decide if cells owned by rank are ordered before ghost cells.
Definition: FlowGenericVanguard.hpp:260
Dune::EdgeWeightMethod edgeWeightsMethod() const
Parameter deciding the edge-weight strategy of the load balancer.
Definition: FlowGenericVanguard.hpp:242
Dune::PartitionMethod partitionMethod() const
Parameter deciding which partition method to use.
Definition: FlowGenericVanguard.hpp:276
bool edgeConformal() const
Definition: FlowGenericVanguard.hpp:263
bool allow_splitting_inactive_wells_
Definition: FlowGenericVanguard.hpp:401
std::string metisParams_
Definition: FlowGenericVanguard.hpp:394
int numOverlap() const
Definition: FlowGenericVanguard.hpp:270
Definition: GenericCpGridVanguard.hpp:79
const CartesianIndexMapper & cartesianIndexMapper() const
Returns the object which maps a global element index of the simulation grid to the corresponding elem...
void doLoadBalance_(const Dune::EdgeWeightMethod edgeWeightsMethod, const bool ownersFirst, const bool addCorners, const int numOverlap, const Dune::PartitionMethod partitionMethod, const bool serialPartitioning, const bool enableDistributedWells, const bool allowSplittingInactiveWells, const double imbalanceTol, const GetPropType< TypeTag, Properties::GridView > &gridView, const Schedule &schedule, EclipseState &eclState, FlowGenericVanguard::ParallelWellStruct &parallelWells, const int numJacobiBlocks, const bool enableEclOutput)
Distribute the simulation grid over multiple processes.
static void setExternalLoadBalancer(const std::function< std::vector< int >(const Dune::CpGrid &)> &loadBalancer)
Sets a function that returns external load balancing information when passed the grid.
Definition: GenericCpGridVanguard.hpp:125
Definition: Transmissibility.hpp:54
Definition: GenericCpGridVanguard.hpp:57
Defines the common properties required by the porous medium multi-phase models.
Definition: blackoilmodel.hh:79
Definition: blackoilbioeffectsmodules.hh:43
typename Properties::Detail::GetPropImpl< TypeTag, Property >::type::type GetPropType
get the type alias defined in the property (equivalent to old macro GET_PROP_TYPE(....
Definition: propertysystem.hh:233
GetPropType< TypeTag, Properties::Grid > type
Definition: CpGridVanguard.hpp:71
Definition: FlowBaseVanguard.hpp:70
Dune::CpGrid type
Definition: CpGridVanguard.hpp:67
The type of the DUNE grid.
Definition: basicproperties.hh:100
Definition: CpGridVanguard.hpp:55
std::tuple< FlowBaseVanguard > InheritsFrom
Definition: CpGridVanguard.hpp:56
Property which provides a Vanguard (manages grids)
Definition: basicproperties.hh:96