3.3.0
DUNE for Multi-{Phase, Component, Scale, Physics, ...} flow and transport in porous media
2p2c/sequential/variableclassadaptive.hh
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 * See the file COPYING for full copying permissions. *
5 * *
6 * This program 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 * This program 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 this program. If not, see <http://www.gnu.org/licenses/>. *
18 *****************************************************************************/
24#ifndef DUMUX_VARIABLECLASS2P2C_ADAPTIVE_HH
25#define DUMUX_VARIABLECLASS2P2C_ADAPTIVE_HH
26
27#include <dune/common/power.hh>
28
29// for parallelization
31
32namespace Dumux {
45template<class TypeTag>
47{
48private:
51
56 using AdaptedValues = typename CellData::AdaptedValues;
57
58 using Grid = typename GridView::Grid;
59 using IdType = typename Grid::LocalIdSet::IdType;
60 using IntersectionIterator = typename GridView::IntersectionIterator;
61 using Intersection = typename GridView::Intersection;
62 enum
63 {
64 dim = GridView::dimension
65 };
66 enum
67 {
68 first = 0,
69 second = 1,
70 diagonal1 = 2,
71 diagonal2 = 3
72 };
73 // convenience shortcuts for Vectors/Matrices
74 using GlobalPosition = Dune::FieldVector<Scalar, GridView::dimensionworld>;
75 using TransmissivityMatrix = Dune::FieldVector<Scalar,dim+1>;
76
77protected:
82 const static int storageRequirement = Dune::StaticPower<dim-1, dim>::power;
83
91 struct mpfaData
92 {
93 TransmissivityMatrix T1_[(dim-1)*2];
96 std::vector<IntersectionIterator> secondHalfEdgeIntersection_;
98
101 {
103 }
104
112 void setIntersection(IntersectionIterator& is23)
113 {
114 secondHalfEdgeIntersection_.push_back(is23);
116 }
118 const IntersectionIterator& getIntersection()
119 {
121 }
122 };
123 std::map<IdType, mpfaData> irregularInterfaceMap_;
124 const Grid& grid_;
125
126public:
134 {}
135
144 void adaptVariableSize(int size)
145 {
147 // clear mapper
149 }
158 void reconstructPrimVars(Problem& problem)
159 {
161
162 problem.pressureModel().adaptPressure();
163 }
164
179 void storeMpfaData(const typename GridView::Intersection & irregularIs,
180 const TransmissivityMatrix& T1,
181 const GlobalPosition& globalPos3,
182 const int& globalIdx3)
183 {
184 IdType intersectionID = grid_.localIdSet().subId(
185 irregularIs.inside(), irregularIs.indexInInside(), 1);
186 // mapping is only unique from smaller cell (if *inside and not *outside)
187 if (irregularIs.inside().level() < irregularIs.outside().level())
188 {
189 // IS is regarded from larger cell: get the unique number as seen from smaller
190 intersectionID = grid_.localIdSet().subId(
191 irregularIs.outside(), irregularIs.indexInOutside(), 1);
192
193 // store as if it was seen from smaller: change i & j
194 irregularInterfaceMap_[intersectionID].T1_[first][2] = - T1[0];
195 irregularInterfaceMap_[intersectionID].T1_[first][1] = - T1[1];
196 irregularInterfaceMap_[intersectionID].T1_[first][0] = - T1[2];
197 }
198 else // we look from smaller cell = unique interface
199 // proceed with numbering according to Aavatsmark, seen from cell i
200 irregularInterfaceMap_[intersectionID].T1_[first] = T1;
201
202 irregularInterfaceMap_[intersectionID].globalPos3_[0] = globalPos3;
203 irregularInterfaceMap_[intersectionID].globalIdx3_[0] = globalIdx3;
204 }
205
222 void storeMpfaData(const typename GridView::Intersection & irregularIs,
223 IntersectionIterator& secondHalfEdgeIntersectionIt,
224 const TransmissivityMatrix& T1,
225 const TransmissivityMatrix& T1_secondHalfEdge,
226 const GlobalPosition& globalPos3,
227 const int& globalIdx3)
228 {
229 IdType intersectionID
230 = grid_.localIdSet().subId(irregularIs.inside(),
231 irregularIs.indexInInside(), 1);
232
233 // mapping is only unique from smaller cell (if *inside and not *outside)
234 if (irregularIs.inside().level() < irregularIs.outside().level())
235 {
236 // IS is regarded from larger cell: get the unique number as seen from smaller
237 intersectionID
238 = grid_.localIdSet().subId(irregularIs.outside(),
239 irregularIs.indexInOutside(), 1);
240
241 // store as if it was seen from smaller: change i & j
242 irregularInterfaceMap_[intersectionID].T1_[first][2] = - T1[0];
243 irregularInterfaceMap_[intersectionID].T1_[first][1] = - T1[1];
244 irregularInterfaceMap_[intersectionID].T1_[first][0] = - T1[2];
245
246 irregularInterfaceMap_[intersectionID].T1_[second][2] = - T1_secondHalfEdge[0];
247 irregularInterfaceMap_[intersectionID].T1_[second][1] = - T1_secondHalfEdge[1];
248 irregularInterfaceMap_[intersectionID].T1_[second][0] = - T1_secondHalfEdge[2];
249
250 }
251 else // we look from smaller cell = unique interface
252 {
253 // proceed with numbering according to Aavatsmark, seen from cell i
254 irregularInterfaceMap_[intersectionID].T1_[first] = T1;
255 irregularInterfaceMap_[intersectionID].T1_[second] = T1_secondHalfEdge;
256 }
257
258 irregularInterfaceMap_[intersectionID].globalPos3_[0] = globalPos3;
259 irregularInterfaceMap_[intersectionID].globalIdx3_[0] = globalIdx3;
260 // second half edge
261 irregularInterfaceMap_[intersectionID].setIntersection(secondHalfEdgeIntersectionIt);
262 }
263
281 void storeMpfaData3D(const typename GridView::Intersection & irregularIs,
282 const TransmissivityMatrix& T1,
283 const GlobalPosition& globalPos3,
284 const int& globalIdx3,
285 const GlobalPosition& globalPos4,
286 const int& globalIdx4,
287 int subFaceIdx = 0)
288 {
289 // if the second interaction Region (subfaceIdx=1) needs to be stored,
290 // we need an offset for the second globalPos and globalIdxes
291 const int offset = subFaceIdx * 2;
292
293 IdType intersectionID
294 = grid_.localIdSet().subId(irregularIs.inside(),
295 irregularIs.indexInInside(), 1);
296
297 // mapping is only unique from smaller cell (if *inside and not *outside)
298 if (irregularIs.inside().level() < irregularIs.outside().level())
299 {
300 // IS is regarded from larger cell: get the unique ID as seen from smaller
301 intersectionID
302 = grid_.localIdSet().subId(irregularIs.outside(),
303 irregularIs.indexInOutside(), 1);
304
305 // store as if it was seen from smaller: change i & j
306 irregularInterfaceMap_[intersectionID].T1_[subFaceIdx][0] = -T1[1];
307 irregularInterfaceMap_[intersectionID].T1_[subFaceIdx][1] = -T1[0];
308 irregularInterfaceMap_[intersectionID].T1_[subFaceIdx][2] = -T1[2];
309 irregularInterfaceMap_[intersectionID].T1_[subFaceIdx][3] = -T1[3];
310 }
311 else
312 {
313 // proceed with numbering according to case2
314 irregularInterfaceMap_[intersectionID].T1_[subFaceIdx][0] = T1[0];
315 irregularInterfaceMap_[intersectionID].T1_[subFaceIdx][1] = T1[1];
316 irregularInterfaceMap_[intersectionID].T1_[subFaceIdx][2] = T1[2];
317 irregularInterfaceMap_[intersectionID].T1_[subFaceIdx][3] = T1[3];
318 }
319 irregularInterfaceMap_[intersectionID].globalPos3_[offset+0] = globalPos3;
320 irregularInterfaceMap_[intersectionID].globalIdx3_[offset+0] = globalIdx3;
321 irregularInterfaceMap_[intersectionID].globalPos3_[offset+1] = globalPos4;
322 irregularInterfaceMap_[intersectionID].globalIdx3_[offset+1] = globalIdx4;
323
324 using std::max;
325 irregularInterfaceMap_[intersectionID].interactionRegionsStored
326 = max(irregularInterfaceMap_[intersectionID].interactionRegionsStored, subFaceIdx+1);
327 }
328
342 void performTransmissitivityWeighting(const typename GridView::Intersection & irregularIs,
343 Scalar weight, int subFaceIdx = -1)
344 {
345 IdType intersectionID
346 = grid_.localIdSet().subId(irregularIs.inside(),
347 irregularIs.indexInInside(), 1);
348
349 // mapping is only unique from smaller cell (if *inside and not *outside)
350 if (irregularIs.inside().level() < irregularIs.outside().level())
351 {
352 // IS is regarded from larger cell: get the unique ID as seen from smaller
353 intersectionID
354 = grid_.localIdSet().subId(irregularIs.outside(),
355 irregularIs.indexInOutside(), 1);
356 }
357
358 // for subFaceIdx == -1, we weight all subfaces equally
359 if(subFaceIdx == -1)
360 {
361 for(int i = 0; i < irregularInterfaceMap_[intersectionID].interactionRegionsStored; i++)
362 irregularInterfaceMap_[intersectionID].T1_[i] *= weight;
363 }
364 else //weight specifically
365 irregularInterfaceMap_[intersectionID].T1_[subFaceIdx] *= weight;
366 }
367
384 int getMpfaData(const Intersection& irregularIs,
385 IntersectionIterator& secondHalfEdgeIntersectionIt,
386 TransmissivityMatrix& T1,
387 TransmissivityMatrix& T1_secondHalfEdge,
388 GlobalPosition& globalPos3,
389 int& globalIdx3)
390 {
391 IdType intersectionID
392 = grid_.localIdSet().subId(irregularIs.inside(),
393 irregularIs.indexInInside(), 1);
394 // mapping is only unique from smaller cell (if *inside and not *outside)
395 if (irregularIs.inside().level() < irregularIs.outside().level())
396 {
397 // IS is regarded from larger cell: get the unique number as seen from smaller
398 intersectionID
399 = grid_.localIdSet().subId(irregularIs.outside(),
400 irregularIs.indexInOutside(), 1);
401
402 // check if T1ransmissibility matrix was stored for that IF
403 if (irregularInterfaceMap_.find(intersectionID) == irregularInterfaceMap_.end())
404 return 0; // no stored data!
405
406 // If data is stored, it is so as if the IF is regarded from smaller cell.
407 // since we are looking from larger cell, cells i & j have to be changed
408 // Additionally, flux points in opposite direction: - sign
409 T1[0] = -irregularInterfaceMap_[intersectionID].T1_[first][2];
410 T1[1] = -irregularInterfaceMap_[intersectionID].T1_[first][1];
411 T1[2] = -irregularInterfaceMap_[intersectionID].T1_[first][0];
412 globalPos3 = irregularInterfaceMap_[intersectionID].globalPos3_[0];
413 globalIdx3 = irregularInterfaceMap_[intersectionID].globalIdx3_[0];
414 //second half edge
415 if(irregularInterfaceMap_[intersectionID].interactionRegionsStored == 2)
416 {
417 secondHalfEdgeIntersectionIt = irregularInterfaceMap_[intersectionID].getIntersection();
418 T1_secondHalfEdge[0] = -irregularInterfaceMap_[intersectionID].T1_[second][2];
419 T1_secondHalfEdge[1] = -irregularInterfaceMap_[intersectionID].T1_[second][1];
420 T1_secondHalfEdge[2] = -irregularInterfaceMap_[intersectionID].T1_[second][0];
421 return 2;
422 }
423 return 1;
424 }
425 // check if T1ransmissibility matrix was stored for that IF
426 if (irregularInterfaceMap_.find(intersectionID) == irregularInterfaceMap_.end())
427 return 0; // no stored data!
428
429 T1[0] = irregularInterfaceMap_[intersectionID].T1_[first][0];
430 T1[1] = irregularInterfaceMap_[intersectionID].T1_[first][1];
431 T1[2] = irregularInterfaceMap_[intersectionID].T1_[first][2];
432 globalPos3 = irregularInterfaceMap_[intersectionID].globalPos3_[0];
433 globalIdx3 = irregularInterfaceMap_[intersectionID].globalIdx3_[0];
434 //second half edge
435 if(irregularInterfaceMap_[intersectionID].interactionRegionsStored == 2)
436 {
437 secondHalfEdgeIntersectionIt = irregularInterfaceMap_[intersectionID].getIntersection();
438 T1_secondHalfEdge = irregularInterfaceMap_[intersectionID].T1_[second];
439 return 2;
440 }
441 return 1;
442 }
443
461 int getMpfaData3D(const Intersection& irregularIs,
462 TransmissivityMatrix& T1,
463 GlobalPosition& globalPos3,
464 int& globalIdx3,
465 GlobalPosition& globalPos4,
466 int& globalIdx4,
467 int subFaceIdx = 0)
468 {
469 // if the second interaction Region (subfaceIdx=1) needs to be stored,
470 // we need an offset for the second globalPos and globalIdxes
471 const int offset = subFaceIdx * 2;
472
473 IdType intersectionID
474 = grid_.localIdSet().subId(irregularIs.inside(),
475 irregularIs.indexInInside(), 1);
476 // mapping is only unique from smaller cell (if *inside and not *outside)
477 if (irregularIs.inside().level() < irregularIs.outside().level())
478 {
479 // IS is regarded from larger cell: get the unique number as seen from smaller
480 intersectionID
481 = grid_.localIdSet().subId(irregularIs.outside(),
482 irregularIs.indexInOutside(), 1);
483
484 // check if T1ransmissibility matrix was stored for that IF
485 if (irregularInterfaceMap_.find(intersectionID) == irregularInterfaceMap_.end())
486 return 0; // no stored data!
487
488 // If data is stored, it is so as if the IF is regarded from smaller cell.
489 // since we are looking from larger cell, cells i & j have to be changed
490 // Additionally, flux points in opposite direction: - sign
491 T1[0] = -irregularInterfaceMap_[intersectionID].T1_[subFaceIdx][1];
492 T1[1] = -irregularInterfaceMap_[intersectionID].T1_[subFaceIdx][0];
493 T1[2] = -irregularInterfaceMap_[intersectionID].T1_[subFaceIdx][2];
494 T1[3] = -irregularInterfaceMap_[intersectionID].T1_[subFaceIdx][3];
495
496 }
497 else
498 {
499 // check if T1ransmissibility matrix was stored for that IF
500 if (irregularInterfaceMap_.find(intersectionID) == irregularInterfaceMap_.end())
501 return 0; // no stored data!
502
503 T1[0] = irregularInterfaceMap_[intersectionID].T1_[subFaceIdx][0];
504 T1[1] = irregularInterfaceMap_[intersectionID].T1_[subFaceIdx][1];
505 T1[2] = irregularInterfaceMap_[intersectionID].T1_[subFaceIdx][2];
506 T1[3] = irregularInterfaceMap_[intersectionID].T1_[subFaceIdx][3];
507 }
508
509 // return what does not depend on direction: additional cells
510 globalPos3 = irregularInterfaceMap_[intersectionID].globalPos3_[offset+0];
511 globalIdx3 = irregularInterfaceMap_[intersectionID].globalIdx3_[offset+0];
512 globalPos4 = irregularInterfaceMap_[intersectionID].globalPos3_[offset+1];
513 globalIdx4 = irregularInterfaceMap_[intersectionID].globalIdx3_[offset+1];
514
515 return irregularInterfaceMap_[intersectionID].interactionRegionsStored;
516 }
518
519};
520} // end namespace Dumux
521#endif
Definition: adapt.hh:29
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:149
Base class holding the variables and discretized data for sequential models.
Definition: variableclass.hh:49
const GridView & gridView() const
Return gridView.
Definition: variableclass.hh:156
void adaptVariableSize(const int size)
Resizes sequential variable vectors.
Definition: variableclass.hh:102
Class holding additionally mpfa data of adaptive compositional models.
Definition: 2p2c/sequential/variableclassadaptive.hh:47
VariableClass2P2CAdaptive(const GridView &gridView)
Constructs a VariableClass object.
Definition: 2p2c/sequential/variableclassadaptive.hh:132
const Grid & grid_
The grid.
Definition: 2p2c/sequential/variableclassadaptive.hh:124
int getMpfaData3D(const Intersection &irregularIs, TransmissivityMatrix &T1, GlobalPosition &globalPos3, int &globalIdx3, GlobalPosition &globalPos4, int &globalIdx4, int subFaceIdx=0)
Provides access to stored 3D Mpfa Data on an intersection for up to 4 subfaces.
Definition: 2p2c/sequential/variableclassadaptive.hh:461
int getMpfaData(const Intersection &irregularIs, IntersectionIterator &secondHalfEdgeIntersectionIt, TransmissivityMatrix &T1, TransmissivityMatrix &T1_secondHalfEdge, GlobalPosition &globalPos3, int &globalIdx3)
Provides access to stored Mpfa Data on an intersection for both half-edges.
Definition: 2p2c/sequential/variableclassadaptive.hh:384
std::map< IdType, mpfaData > irregularInterfaceMap_
Storage container for mpfa information.
Definition: 2p2c/sequential/variableclassadaptive.hh:123
void reconstructPrimVars(Problem &problem)
Reconstruct missing primary variables (where elements are created/deleted)
Definition: 2p2c/sequential/variableclassadaptive.hh:158
void storeMpfaData(const typename GridView::Intersection &irregularIs, IntersectionIterator &secondHalfEdgeIntersectionIt, const TransmissivityMatrix &T1, const TransmissivityMatrix &T1_secondHalfEdge, const GlobalPosition &globalPos3, const int &globalIdx3)
Stores Mpfa Data on an intersection for both half-edges.
Definition: 2p2c/sequential/variableclassadaptive.hh:222
void adaptVariableSize(int size)
Resizes sequential variable vectors.
Definition: 2p2c/sequential/variableclassadaptive.hh:144
void storeMpfaData(const typename GridView::Intersection &irregularIs, const TransmissivityMatrix &T1, const GlobalPosition &globalPos3, const int &globalIdx3)
Stores Mpfa Data of one interaction region on an intersection.
Definition: 2p2c/sequential/variableclassadaptive.hh:179
void performTransmissitivityWeighting(const typename GridView::Intersection &irregularIs, Scalar weight, int subFaceIdx=-1)
Weigths the transmissivity coefficient by the flux area (3D)
Definition: 2p2c/sequential/variableclassadaptive.hh:342
static const int storageRequirement
Definition: 2p2c/sequential/variableclassadaptive.hh:82
void storeMpfaData3D(const typename GridView::Intersection &irregularIs, const TransmissivityMatrix &T1, const GlobalPosition &globalPos3, const int &globalIdx3, const GlobalPosition &globalPos4, const int &globalIdx4, int subFaceIdx=0)
Stores 3D Mpfa Data on an intersection.
Definition: 2p2c/sequential/variableclassadaptive.hh:281
Storage object for data related to the MPFA method.
Definition: 2p2c/sequential/variableclassadaptive.hh:92
mpfaData()
Constructor for the local storage object of mpfa data.
Definition: 2p2c/sequential/variableclassadaptive.hh:100
void setIntersection(IntersectionIterator &is23)
Stores an intersection for the 2D implementation.
Definition: 2p2c/sequential/variableclassadaptive.hh:112
int interactionRegionsStored
Definition: 2p2c/sequential/variableclassadaptive.hh:97
const IntersectionIterator & getIntersection()
Acess method to the stored intersection (only 2D)
Definition: 2p2c/sequential/variableclassadaptive.hh:118
TransmissivityMatrix T1_[(dim-1) *2]
Definition: 2p2c/sequential/variableclassadaptive.hh:93
GlobalPosition globalPos3_[storageRequirement]
Definition: 2p2c/sequential/variableclassadaptive.hh:95
int globalIdx3_[storageRequirement]
Definition: 2p2c/sequential/variableclassadaptive.hh:94
std::vector< IntersectionIterator > secondHalfEdgeIntersection_
Definition: 2p2c/sequential/variableclassadaptive.hh:96
Base class holding the variables and discretized data for sequential models.
Definition: sequential/variableclassadaptive.hh:47
void reconstructPrimVars(const Problem &problem)
Definition: sequential/variableclassadaptive.hh:140
Base class holding the variables for sequential models.