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::EnableMICP>() == false) {
159 throw std::runtime_error("Input specifies MICP while simulator has it disabled");
160 }
161 }
162
163 if (phases.active(Phase::BRINE)) {
164 if (getPropValue<TypeTag, Properties::EnableBrine>() == false) {
165 throw std::runtime_error("Input specifies Brine while simulator has it disabled");
166 }
167 }
168
169 if (phases.active(Phase::POLYMER)) {
170 if (getPropValue<TypeTag, Properties::EnablePolymer>() == false) {
171 throw std::runtime_error("Input specifies Polymer while simulator has it disabled");
172 }
173 }
174
175 // checking for correct phases is more difficult TODO!
176 if (phases.active(Phase::ZFRACTION)) {
177 if (getPropValue<TypeTag, Properties::EnableExtbo>() == false) {
178 throw std::runtime_error("Input specifies ExBo while simulator has it disabled");
179 }
180 }
181 if (phases.active(Phase::FOAM)) {
182 if (getPropValue<TypeTag, Properties::EnableFoam>() == false) {
183 throw std::runtime_error("Input specifies Foam while simulator has it disabled");
184 }
185 }
186
187 if (phases.active(Phase::SOLVENT)) {
188 if (getPropValue<TypeTag, Properties::EnableSolvent>() == false) {
189 throw std::runtime_error("Input specifies Solvent while simulator has it disabled");
190 }
191 }
192 if(phases.active(Phase::WATER)){
193 if(waterEnabled == false){
194 throw std::runtime_error("Input specifies water while simulator has it disabled");
195 }
196 }
197 if(phases.active(Phase::GAS)){
198 if(gasEnabled == false){
199 throw std::runtime_error("Input specifies gas while simulator has it disabled");
200 }
201 }
202 if(phases.active(Phase::OIL)){
203 if(oilEnabled == false){
204 throw std::runtime_error("Input specifies oil while simulator has it disabled");
205 }
206 }
207
208 }
209
216 {
217 globalTrans_.reset();
218 }
219
221 {
222 assert( globalTrans_ != nullptr );
223 return *globalTrans_;
224 }
225
232 {
233#if HAVE_MPI
234 if (const auto& extPFile = this->externalPartitionFile();
235 !extPFile.empty() && (extPFile != "none"))
236 {
238 }
239
240 this->doLoadBalance_(this->edgeWeightsMethod(), this->ownersFirst(),
241 this->addCorners(), this->numOverlap(),
242 this->partitionMethod(), this->serialPartitioning(),
245 this->imbalanceTol(),
246 this->gridView(), this->schedule(),
247 this->eclState(), this->parallelWells_,
248 this->numJacobiBlocks(), this->enableEclOutput());
249#endif
250
251 this->updateGridView_();
253 this->updateCellDepths_();
254 this->updateCellThickness_();
255
256#if HAVE_MPI
257 this->distributeFieldProps_(this->eclState());
258#endif
259 }
260
264 void addLgrs()
265 {
266 // Check if input file contains Lgrs. Add them, if any.
267 // In a parallel run, this adds the LGRs on the distributed simulation grid.
268 if (const auto& lgrs = this->eclState().getLgrs(); lgrs.size() > 0) {
269 OpmLog::info("\nAdding LGRs to the grid and updating its leaf grid view");
270 this->addLgrsUpdateLeafView(lgrs, lgrs.size(), *this->grid_);
271
272 this->updateGridView_();
273 this->updateCellDepths_();
274 this->updateCellThickness_();
275
276 if (this->grid_->comm().size()>1) {
277 // Add LGRs and update the leaf grid view in the global (undistributed) simulation grid.
278 // Purpose: To enable synchronization of cell ids in 'serial mode',
279 // we rely on the "parent-to-children" cell id mapping.
280 OpmLog::info("\nAdding LGRs to the global view and updating its leaf grid view");
281 this->grid_->switchToGlobalView();
282 this->addLgrsUpdateLeafView(lgrs, lgrs.size(), *this->grid_);
283 this->grid_->switchToDistributedView();
284 this->grid_->syncDistributedGlobalCellIds();
285 }
286 }
287 }
288
289 unsigned int gridEquilIdxToGridIdx(unsigned int elemIndex) const {
290 return elemIndex;
291 }
292
293 unsigned int gridIdxToEquilGridIdx(unsigned int elemIndex) const {
294 return elemIndex;
295 }
303 std::function<std::array<double,dimensionworld>(int)>
305 {
306 return this->cellCentroids_(this->cartesianIndexMapper(), true);
307 }
308
309 const std::vector<int>& globalCell()
310 {
311 return this->grid().globalCell();
312 }
313
314protected:
316 {
317 this->doCreateGrids_(this->eclState());
318 }
319
320 void allocTrans() override
321 {
322 OPM_TIMEBLOCK(allocateTrans);
324 this->gridView(),
325 this->cartesianIndexMapper(),
326 this->grid(),
327 this->cellCentroids(),
328 getPropValue<TypeTag, Properties::EnableEnergy>(),
329 getPropValue<TypeTag, Properties::EnableDiffusion>(),
330 getPropValue<TypeTag, Properties::EnableDispersion>()));
332 }
333
334 double getTransmissibility(unsigned I, unsigned J) const override
335 {
336 return globalTrans_->transmissibility(I,J);
337 }
338
339#if HAVE_MPI
340 const std::string& zoltanParams() const override
341 {
342 return this->zoltanParams_;
343 }
344
345 double zoltanPhgEdgeSizeThreshold() const override
346 {
347 return this->zoltanPhgEdgeSizeThreshold_;
348 }
349
350 const std::string& metisParams() const override
351 {
352 return this->metisParams_;
353 }
354#endif
355
356 // removing some connection located in inactive grid cells
358 {
359 this->doFilterConnections_(this->schedule());
360 }
361
362 // \Note: this globalTrans_ is used for domain decomposition and INIT file output.
363 // It only contains trans_ due to permeability and does not contain thermalHalfTrans_,
364 // diffusivity_ abd dispersivity_. The main reason is to reduce the memory usage for rank 0
365 // during parallel running.
366 std::unique_ptr<TransmissibilityType> globalTrans_;
367};
368
369} // namespace Opm
370
371#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:215
void addLgrs()
Add LGRs and update Leaf Grid View in the simulation grid.
Definition: CpGridVanguard.hpp:264
const std::string & zoltanParams() const override
Definition: CpGridVanguard.hpp:340
double zoltanPhgEdgeSizeThreshold() const override
Definition: CpGridVanguard.hpp:345
void filterConnections_()
Definition: CpGridVanguard.hpp:357
const TransmissibilityType & globalTransmissibility() const
Definition: CpGridVanguard.hpp:220
void allocTrans() override
Definition: CpGridVanguard.hpp:320
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:315
unsigned int gridIdxToEquilGridIdx(unsigned int elemIndex) const
Definition: CpGridVanguard.hpp:293
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:289
const std::string & metisParams() const override
Definition: CpGridVanguard.hpp:350
std::function< std::array< double, dimensionworld >(int)> cellCentroids() const
Get function to query cell centroids for a distributed grid.
Definition: CpGridVanguard.hpp:304
std::unique_ptr< TransmissibilityType > globalTrans_
Definition: CpGridVanguard.hpp:366
static constexpr bool oilEnabled
Definition: CpGridVanguard.hpp:108
void loadBalance()
Distribute the simulation grid over multiple processes.
Definition: CpGridVanguard.hpp:231
Transmissibility< Grid, GridView, ElementMapper, CartesianIndexMapper, Scalar > TransmissibilityType
Definition: CpGridVanguard.hpp:103
const std::vector< int > & globalCell()
Definition: CpGridVanguard.hpp:309
double getTransmissibility(unsigned I, unsigned J) const override
Definition: CpGridVanguard.hpp:334
Helper class for grid instantiation of ECL file-format using problems.
Definition: FlowBaseVanguard.hpp:84
void updateCartesianToCompressedMapping_()
Definition: FlowBaseVanguard.hpp:333
void updateCellThickness_()
Definition: FlowBaseVanguard.hpp:378
void updateCellDepths_()
Definition: FlowBaseVanguard.hpp:355
void callImplementationInit()
Definition: FlowBaseVanguard.hpp:321
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:447
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:301
int numJacobiBlocks() const
Number of blocks in the Block-Jacobi preconditioner.
Definition: FlowGenericVanguard.hpp:247
const Schedule & schedule() const
Return a reference to the object that managages the ECL schedule.
Definition: FlowGenericVanguard.hpp:176
std::string zoltanParams_
Definition: FlowGenericVanguard.hpp:386
bool addCorners() const
Definition: FlowGenericVanguard.hpp:263
bool serialPartitioning() const
Parameter that decides if partitioning for parallel runs should be performed on a single process only...
Definition: FlowGenericVanguard.hpp:279
ParallelWellStruct parallelWells_
Information about wells in parallel.
Definition: FlowGenericVanguard.hpp:425
bool enableDistributedWells() const
Whether perforations of a well might be distributed.
Definition: FlowGenericVanguard.hpp:307
const EclipseState & eclState() const
Return a reference to the internalized ECL deck.
Definition: FlowGenericVanguard.hpp:167
double imbalanceTol() const
Parameter that sets the imbalance tolarance, depending on the chosen partition method.
Definition: FlowGenericVanguard.hpp:286
double zoltanPhgEdgeSizeThreshold_
Definition: FlowGenericVanguard.hpp:385
const std::string & externalPartitionFile() const
Definition: FlowGenericVanguard.hpp:298
bool enableEclOutput() const
Whether or not to emit result files that are compatible with a commercial reservoir simulator.
Definition: FlowGenericVanguard.hpp:314
bool ownersFirst() const
Parameter that decide if cells owned by rank are ordered before ghost cells.
Definition: FlowGenericVanguard.hpp:259
Dune::EdgeWeightMethod edgeWeightsMethod() const
Parameter deciding the edge-weight strategy of the load balancer.
Definition: FlowGenericVanguard.hpp:241
Dune::PartitionMethod partitionMethod() const
Parameter deciding which partition method to use.
Definition: FlowGenericVanguard.hpp:272
bool allow_splitting_inactive_wells_
Definition: FlowGenericVanguard.hpp:395
std::string metisParams_
Definition: FlowGenericVanguard.hpp:388
int numOverlap() const
Definition: FlowGenericVanguard.hpp:266
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: blackoilboundaryratevector.hh:39
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