3.1-git
DUNE for Multi-{Phase, Component, Scale, Physics, ...} flow and transport in porous media
lmethod/2dpressurevelocity.hh
Go to the documentation of this file.
1/*****************************************************************************
2 * See the file COPYING for full copying permissions. *
3 * *
4 * This program is free software: you can redistribute it and/or modify *
5 * it under the terms of the GNU General Public License as published by *
6 * the Free Software Foundation, either version 3 of the License, or *
7 * (at your option) any later version. *
8 * *
9 * This program is distributed in the hope that it will be useful, *
10 * but WITHOUT ANY WARRANTY; without even the implied warranty of *
11 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the *
12 * GNU General Public License for more details. *
13 * *
14 * You should have received a copy of the GNU General Public License *
15 * along with this program. If not, see <http://www.gnu.org/licenses/>. *
16 *****************************************************************************/
22#ifndef DUMUX_MPFAL2DPRESSUREVELOCITY2P_HH
23#define DUMUX_MPFAL2DPRESSUREVELOCITY2P_HH
24
25#include <dune/common/float_cmp.hh>
26
27#include "2dpressure.hh"
28#include "2dvelocity.hh"
29
30namespace Dumux {
31
45template<class TypeTag> class FvMpfaL2dPressureVelocity2p: public FvMpfaL2dPressure2p<TypeTag>
46{
48
49 using GridView = typename GET_PROP_TYPE(TypeTag, GridView);
50
51 enum
52 {
53 dim = GridView::dimension, dimWorld = GridView::dimensionworld
54 };
55
56 using Scalar = typename GET_PROP_TYPE(TypeTag, Scalar);
57 using Problem = typename GET_PROP_TYPE(TypeTag, Problem);
58 using CellData = typename GET_PROP_TYPE(TypeTag, CellData);
59 using BoundaryTypes = typename GET_PROP_TYPE(TypeTag, BoundaryTypes);
60 using SolutionTypes = typename GET_PROP(TypeTag, SolutionTypes);
61 using PrimaryVariables = typename SolutionTypes::PrimaryVariables;
62
63 using ReferenceElements = Dune::ReferenceElements<Scalar, dim>;
64
65 using SpatialParams = typename GET_PROP_TYPE(TypeTag, SpatialParams);
66 using MaterialLaw = typename SpatialParams::MaterialLaw;
67
68 using Indices = typename GET_PROP_TYPE(TypeTag, ModelTraits)::Indices;
69
70 using FluidSystem = typename GET_PROP_TYPE(TypeTag, FluidSystem);
71 using FluidState = typename GET_PROP_TYPE(TypeTag, FluidState);
72
73 using IndexSet = typename GridView::IndexSet;
74 using Intersection = typename GridView::Intersection;
75
77
78 enum
79 {
80 pw = Indices::pressureW,
81 pn = Indices::pressureNw,
82 sw = Indices::saturationW,
83 sn = Indices::saturationNw,
84 wPhaseIdx = Indices::wPhaseIdx,
85 nPhaseIdx = Indices::nPhaseIdx,
86 pressureIdx = Indices::pressureIdx,
87 saturationIdx = Indices::saturationIdx,
88 pressEqIdx = Indices::pressureEqIdx,
89 satEqIdx = Indices::satEqIdx,
90 numPhases = GET_PROP_VALUE(TypeTag, NumPhases)
91 };
92
93 using Element = typename GridView::template Codim<0>::Entity;
94
95 using GlobalPosition = typename Element::Geometry::GlobalCoordinate;
96 using DimMatrix = Dune::FieldMatrix<Scalar, dim, dim>;
97
98public:
104 FvMpfaL2dPressureVelocity2p(Problem& problem) :
105 ParentType(problem), problem_(problem), velocity_(problem)
106 {
107 density_[wPhaseIdx] = 0.;
108 density_[nPhaseIdx] = 0.;
109 viscosity_[wPhaseIdx] = 0.;
110 viscosity_[nPhaseIdx] = 0.;
111
112 calcVelocityInTransport_ = getParam<bool>("MPFA.CalcVelocityInTransport");
113 }
114
116 void calculateVelocity();
117
119 void calculateVelocity(const Intersection&, CellData&);
120 void calculateVelocityOnBoundary(const Intersection& intersection, CellData& cellData);
121
124 {
125 this->updateMaterialLaws();
126
127 this->storePressureSolution();
128
131 }
132
139 {
140 const auto element = *problem_.gridView().template begin<0>();
141 FluidState fluidState;
142 fluidState.setPressure(wPhaseIdx, problem_.referencePressure(element));
143 fluidState.setPressure(nPhaseIdx, problem_.referencePressure(element));
144 fluidState.setTemperature(problem_.temperature(element));
145 fluidState.setSaturation(wPhaseIdx, 1.);
146 fluidState.setSaturation(nPhaseIdx, 0.);
147 density_[wPhaseIdx] = FluidSystem::density(fluidState, wPhaseIdx);
148 density_[nPhaseIdx] = FluidSystem::density(fluidState, nPhaseIdx);
149 viscosity_[wPhaseIdx] = FluidSystem::viscosity(fluidState, wPhaseIdx);
150 viscosity_[nPhaseIdx] = FluidSystem::viscosity(fluidState, nPhaseIdx);
151
153 velocity_.initialize();
156
157 return;
158 }
159
165 void update()
166 {
170 }
171
178 {
179 return calcVelocityInTransport_;
180 }
181
192 template<class MultiWriter>
193 void addOutputVtkFields(MultiWriter &writer)
194 {
196 velocity_.addOutputVtkFields(writer);
197 }
198
199private:
200 Problem& problem_;
202
203 Scalar density_[numPhases];
204 Scalar viscosity_[numPhases];
205 bool calcVelocityInTransport_;
206
208 static const int pressureType_ = GET_PROP_VALUE(TypeTag, PressureFormulation);
210 static const int saturationType_ = GET_PROP_VALUE(TypeTag, SaturationFormulation);
211};
212// end of template
213
214
220template<class TypeTag>
222{
223 // run through all elements
224 for (const auto& vertex : vertices(problem_.gridView()))
225 {
226 int vIdxGlobal = problem_.variables().index(vertex);
227
228 InteractionVolume& interactionVolume = this->interactionVolumes_[vIdxGlobal];
229
230 if (interactionVolume.isInnerVolume())
231 {
232 auto element1 = interactionVolume.getSubVolumeElement(0);
233 auto element2 = interactionVolume.getSubVolumeElement(1);
234 auto element3 = interactionVolume.getSubVolumeElement(2);
235 auto element4 = interactionVolume.getSubVolumeElement(3);
236
237 // cell index
238 int eIdxGlobal1 = problem_.variables().index(element1);
239 int eIdxGlobal2 = problem_.variables().index(element2);
240 int eIdxGlobal3 = problem_.variables().index(element3);
241 int eIdxGlobal4 = problem_.variables().index(element4);
242
243 //get the cell Data
244 CellData& cellData1 = problem_.variables().cellData(eIdxGlobal1);
245 CellData& cellData2 = problem_.variables().cellData(eIdxGlobal2);
246 CellData& cellData3 = problem_.variables().cellData(eIdxGlobal3);
247 CellData& cellData4 = problem_.variables().cellData(eIdxGlobal4);
248
249 velocity_.calculateInnerInteractionVolumeVelocity(interactionVolume, cellData1, cellData2,
250 cellData3, cellData4, this->innerBoundaryVolumeFaces_);
251
252 }
253 // at least one face on boundary!
254 else
255 {
256 for (int elemIdx = 0; elemIdx < 2 * dim; elemIdx++)
257 {
258 bool isOutside = false;
259 for (int fIdx = 0; fIdx < dim; fIdx++)
260 {
261 int intVolFaceIdx = interactionVolume.getFaceIndexFromSubVolume(elemIdx, fIdx);
262 if (interactionVolume.isOutsideFace(intVolFaceIdx))
263 {
264 isOutside = true;
265 break;
266 }
267 }
268 if (isOutside)
269 {
270 continue;
271 }
272 // cell index
273 int eIdxGlobal = problem_.variables().index(interactionVolume.getSubVolumeElement(elemIdx));
274 //get the cell Data
275 CellData& cellData = problem_.variables().cellData(eIdxGlobal);
276
277 velocity_.calculateBoundaryInteractionVolumeVelocity(interactionVolume, cellData, elemIdx);
278 }
279
280 } // end boundaries
281
282 } // end vertex iterator
283
284// }
285 return;
286} // end method calcTotalVelocity
287
296template<class TypeTag>
297void FvMpfaL2dPressureVelocity2p<TypeTag>::calculateVelocity(const Intersection& intersection, CellData& cellData)
298{
299 int numVertices = intersection.geometry().corners();
300
301 auto elementI = intersection.inside();
302 auto elementJ = intersection.outside();
303
304 int eIdxGlobalI = problem_.variables().index(elementI);
305 int eIdxGlobalJ = problem_.variables().index(elementJ);
306
307 CellData& cellDataJ = problem_.variables().cellData(eIdxGlobalJ);
308
309 const auto referenceElement = ReferenceElements::general(elementI.type());
310
311 int indexInInside = intersection.indexInInside();
312 int indexInOutside = intersection.indexInOutside();
313
314 Dune::FieldVector<CellData, 4> cellDataTemp;
315
316 for (int vIdx = 0; vIdx < numVertices; vIdx++)
317 {
318 int localVertIdx = referenceElement.subEntity(indexInInside, dim - 1, vIdx, dim);
319
320 int vIdxGlobal = problem_.variables().index(elementI.template subEntity<dim>(localVertIdx));
321
322 InteractionVolume& interactionVolume = this->interactionVolumes_[vIdxGlobal];
323
324 if (interactionVolume.isInnerVolume())
325 {
326
327 auto element1 = interactionVolume.getSubVolumeElement(0);
328 auto element2 = interactionVolume.getSubVolumeElement(1);
329 auto element3 = interactionVolume.getSubVolumeElement(2);
330 auto element4 = interactionVolume.getSubVolumeElement(3);
331
332 // cell index
333 int eIdxGlobal[4];
334 eIdxGlobal[0] = problem_.variables().index(element1);
335 eIdxGlobal[1] = problem_.variables().index(element2);
336 eIdxGlobal[2] = problem_.variables().index(element3);
337 eIdxGlobal[3] = problem_.variables().index(element4);
338
339 //get the cell Data
340 cellDataTemp[0] = problem_.variables().cellData(eIdxGlobal[0]);
341 cellDataTemp[1] = problem_.variables().cellData(eIdxGlobal[1]);
342 cellDataTemp[2] = problem_.variables().cellData(eIdxGlobal[2]);
343 cellDataTemp[3] = problem_.variables().cellData(eIdxGlobal[3]);
344
345 velocity_.calculateInnerInteractionVolumeVelocity(interactionVolume, cellDataTemp[0], cellDataTemp[1],
346 cellDataTemp[2], cellDataTemp[3], this->innerBoundaryVolumeFaces_);
347
348 for (int i = 0; i < 4; i++)
349 {
350 if (eIdxGlobal[i] == eIdxGlobalI)
351 {
352 cellData.fluxData().setVelocity(wPhaseIdx, indexInInside,
353 cellDataTemp[i].fluxData().velocity(wPhaseIdx, indexInInside));
354 cellData.fluxData().setVelocity(nPhaseIdx, indexInInside,
355 cellDataTemp[i].fluxData().velocity(nPhaseIdx, indexInInside));
356 cellData.fluxData().setUpwindPotential(wPhaseIdx, indexInInside,
357 cellDataTemp[i].fluxData().upwindPotential(wPhaseIdx, indexInInside));
358 cellData.fluxData().setUpwindPotential(nPhaseIdx, indexInInside,
359 cellDataTemp[i].fluxData().upwindPotential(nPhaseIdx, indexInInside));
360 }
361 else if (eIdxGlobal[i] == eIdxGlobalJ)
362 {
363 cellDataJ.fluxData().setVelocity(wPhaseIdx, indexInOutside,
364 cellDataTemp[i].fluxData().velocity(wPhaseIdx, indexInOutside));
365 cellDataJ.fluxData().setVelocity(nPhaseIdx, indexInOutside,
366 cellDataTemp[i].fluxData().velocity(nPhaseIdx, indexInOutside));
367 cellDataJ.fluxData().setUpwindPotential(wPhaseIdx, indexInOutside,
368 cellDataTemp[i].fluxData().upwindPotential(wPhaseIdx, indexInOutside));
369 cellDataJ.fluxData().setUpwindPotential(nPhaseIdx, indexInOutside,
370 cellDataTemp[i].fluxData().upwindPotential(nPhaseIdx, indexInOutside));
371 }
372 }
373 }
374 }
375 cellData.fluxData().setVelocityMarker(indexInInside);
376 cellDataJ.fluxData().setVelocityMarker(indexInOutside);
377}
378
387template<class TypeTag>
388void FvMpfaL2dPressureVelocity2p<TypeTag>::calculateVelocityOnBoundary(const Intersection& intersection, CellData& cellData)
389{
390 auto element = intersection.inside();
391
392 //get face index
393 int isIndex = intersection.indexInInside();
394
395 //get face normal
396 const Dune::FieldVector<Scalar, dim>& unitOuterNormal = intersection.centerUnitOuterNormal();
397
398 BoundaryTypes bcType;
399 //get boundary type
400 problem_.boundaryTypes(bcType, intersection);
401 PrimaryVariables boundValues(0.0);
402
403 if (bcType.isDirichlet(pressEqIdx))
404 {
405 problem_.dirichlet(boundValues, intersection);
406
407 // get global coordinates of cell centers
408 const GlobalPosition& globalPosI = element.geometry().center();
409
410 // center of face in global coordinates
411 const GlobalPosition& globalPosJ = intersection.geometry().center();
412
413 // get mobilities and fractional flow factors
414 Scalar lambdaWI = cellData.mobility(wPhaseIdx);
415 Scalar lambdaNwI = cellData.mobility(nPhaseIdx);
416
417 // get capillary pressure
418 Scalar pcI = cellData.capillaryPressure();
419
420 // distance vector between barycenters
421 GlobalPosition distVec = globalPosJ - globalPosI;
422
423 // compute distance between cell centers
424 Scalar dist = distVec.two_norm();
425
426 //permeability vector at boundary
427 // compute vectorized permeabilities
428 DimMatrix meanPermeability(0);
429
430 problem_.spatialParams().meanK(meanPermeability, problem_.spatialParams().intrinsicPermeability(element));
431
432 Dune::FieldVector<Scalar, dim> permeability(0);
433 meanPermeability.mv(unitOuterNormal, permeability);
434
435 //determine saturation at the boundary -> if no saturation is known directly at the boundary use the cell saturation
436 Scalar satW = 0;
437 if (bcType.isDirichlet(satEqIdx))
438 {
439 switch (saturationType_)
440 {
441 case sw:
442 {
443 satW = boundValues[saturationIdx];
444 break;
445 }
446 case sn:
447 {
448 satW = 1 - boundValues[saturationIdx];
449 break;
450 }
451 }
452 }
453 else
454 {
455 satW = cellData.saturation(wPhaseIdx);
456 }
457
458 Scalar pressBound = boundValues[pressureIdx];
459 Scalar pcBound = MaterialLaw::pc(problem_.spatialParams().materialLawParams(element), satW);
460
461 //determine phase pressures from primary pressure variable
462 Scalar pressWBound = 0;
463 Scalar pressNwBound = 0;
464 if (pressureType_ == pw)
465 {
466 pressWBound = pressBound;
467 pressNwBound = pressBound + pcBound;
468 }
469 else if (pressureType_ == pn)
470 {
471 pressWBound = pressBound - pcBound;
472 pressNwBound = pressBound;
473 }
474
475 Scalar lambdaWBound = MaterialLaw::krw(problem_.spatialParams().materialLawParams(element), satW)
476 / viscosity_[wPhaseIdx];
477 Scalar lambdaNwBound = MaterialLaw::krn(problem_.spatialParams().materialLawParams(element), satW)
478 / viscosity_[nPhaseIdx];
479
480 Scalar potentialDiffW = cellData.fluxData().upwindPotential(wPhaseIdx, isIndex);
481 Scalar potentialDiffNw = cellData.fluxData().upwindPotential(nPhaseIdx, isIndex);
482
483 //calculate potential gradient
484 potentialDiffW = (cellData.pressure(wPhaseIdx) - pressWBound);
485 potentialDiffNw = (cellData.pressure(nPhaseIdx) - pressNwBound);
486
487 potentialDiffW += density_[wPhaseIdx] * (distVec * problem_.gravity());
488 potentialDiffNw += density_[nPhaseIdx] * (distVec * problem_.gravity());
489
490 //store potential gradients for further calculations
491 cellData.fluxData().setUpwindPotential(wPhaseIdx, isIndex, potentialDiffW);
492 cellData.fluxData().setUpwindPotential(nPhaseIdx, isIndex, potentialDiffNw);
493
494 //do the upwinding of the mobility depending on the phase potentials
495 Scalar lambdaW = (potentialDiffW > 0.) ? lambdaWI : lambdaWBound;
496 lambdaW = (Dune::FloatCmp::eq<Scalar, Dune::FloatCmp::absolute>(potentialDiffW, 0.0, 1.0e-30)) ? 0.5 * (lambdaWI + lambdaWBound) : lambdaW;
497 Scalar lambdaNw = (potentialDiffNw > 0.) ? lambdaNwI : lambdaNwBound;
498 lambdaNw = (Dune::FloatCmp::eq<Scalar, Dune::FloatCmp::absolute>(potentialDiffNw, 0.0, 1.0e-30)) ? 0.5 * (lambdaNwI + lambdaNwBound) : lambdaNw;
499
500
501 Scalar scalarPerm = permeability.two_norm();
502
503 //calculate the gravity term
504 Dune::FieldVector<Scalar, dimWorld> velocityW(unitOuterNormal);
505 Dune::FieldVector<Scalar, dimWorld> velocityNw(unitOuterNormal);
506
507 //calculate unit distVec
508 distVec /= dist;
509 Scalar areaScaling = (unitOuterNormal * distVec);
510 //this treatment of g allows to account for gravity flux through faces where the face normal has no z component (e.g. parallelepiped grids)
511 Scalar gravityTermW = (problem_.gravity() * distVec) * density_[wPhaseIdx] * areaScaling;
512 Scalar gravityTermNw = (problem_.gravity() * distVec) * density_[nPhaseIdx] * areaScaling;
513
514 //calculate velocity depending on the pressure used -> use pc = pn - pw
515 switch (pressureType_)
516 {
517 case pw:
518 {
519 velocityW *= lambdaW * scalarPerm * ((cellData.pressure(wPhaseIdx) - pressBound) / dist + gravityTermW);
520 velocityNw *= lambdaNw * scalarPerm * ((cellData.pressure(wPhaseIdx) - pressBound) / dist + gravityTermNw)
521 + 0.5 * (lambdaNwI + lambdaNwBound) * scalarPerm * (pcI - pcBound) / dist;
522 break;
523 }
524 case pn:
525 {
526 velocityW *= lambdaW * scalarPerm * ((cellData.pressure(nPhaseIdx) - pressBound) / dist + gravityTermW)
527 - 0.5 * (lambdaWI + lambdaWBound) * scalarPerm * (pcI - pcBound) / dist;
528 velocityNw *= lambdaNw * scalarPerm * ((cellData.pressure(nPhaseIdx) - pressBound) / dist + gravityTermNw);
529 break;
530 }
531 }
532
533 //store velocities
534 cellData.fluxData().setVelocity(wPhaseIdx, isIndex, velocityW);
535 cellData.fluxData().setVelocity(nPhaseIdx, isIndex, velocityNw);
536 cellData.fluxData().setVelocityMarker(isIndex);
537
538 } //end dirichlet boundary
539
540 else if (bcType.isNeumann(pressEqIdx))
541 {
542 problem_.neumann(boundValues, intersection);
543
544 Dune::FieldVector<Scalar, dimWorld> velocityW(unitOuterNormal);
545 Dune::FieldVector<Scalar, dimWorld> velocityNw(unitOuterNormal);
546
547 velocityW *= boundValues[wPhaseIdx];
548 velocityNw *= boundValues[nPhaseIdx];
549
550 velocityW /= density_[wPhaseIdx];
551 velocityNw /= density_[nPhaseIdx];
552
553 //store potential gradients for further calculations
554 cellData.fluxData().setUpwindPotential(wPhaseIdx, isIndex, boundValues[wPhaseIdx]);
555 cellData.fluxData().setUpwindPotential(nPhaseIdx, isIndex, boundValues[nPhaseIdx]);
556
557 cellData.fluxData().setVelocity(wPhaseIdx, isIndex, velocityW);
558 cellData.fluxData().setVelocity(nPhaseIdx, isIndex, velocityNw);
559 cellData.fluxData().setVelocityMarker(isIndex);
560 } //end neumann boundary
561 else
562 {
563 DUNE_THROW(Dune::NotImplemented, "No valid boundary condition type defined for pressure equation!");
564 }
565}
566
567} // end namespace Dumux
568#endif
#define GET_PROP_VALUE(TypeTag, PropTagName)
Definition: propertysystemmacros.hh:282
#define GET_PROP(TypeTag, PropTagName)
Definition: propertysystemmacros.hh:281
#define GET_PROP_TYPE(TypeTag, PropTagName)
Definition: propertysystemmacros.hh:283
make the local view function available whenever we use the grid geometry
Definition: adapt.hh:29
Property tag SaturationFormulation
The formulation of the saturation model.
Definition: porousmediumflow/2p/sequential/properties.hh:53
Property tag NumPhases
Number of phases in the system.
Definition: porousmediumflow/sequential/properties.hh:69
std::string viscosity(int phaseIdx) noexcept
I/O name of viscosity for multiphase systems.
Definition: name.hh:74
std::string permeability() noexcept
I/O name of permeability.
Definition: name.hh:143
std::string density(int phaseIdx) noexcept
I/O name of density for multiphase systems.
Definition: name.hh:65
Finite volume MPFA L-method discretization of a two-phase flow pressure equation of the sequential IM...
Definition: lmethod/2dpressure.hh:72
void update()
Pressure update.
Definition: lmethod/2dpressure.hh:291
void initialize()
Initializes the pressure model.
Definition: lmethod/2dpressure.hh:194
FVMPFALInteractionVolume< TypeTag > InteractionVolume
Type of the interaction volume objects.
Definition: lmethod/2dpressure.hh:151
void storePressureSolution()
Globally stores the pressure solution.
Definition: lmethod/2dpressure.hh:228
void addOutputVtkFields(MultiWriter &writer)
Adds pressure output to the output file.
Definition: lmethod/2dpressure.hh:340
void updateMaterialLaws()
constitutive functions are initialized and stored in the variables object
Definition: lmethod/2dpressure.hh:1838
Class for the calculation of velocities from the pressure solution of an IMPES scheme using a MPFA L-...
Definition: lmethod/2dpressurevelocity.hh:46
void updateVelocity()
Function for updating the velocity field if iterations are necessary in the transport solution.
Definition: lmethod/2dpressurevelocity.hh:123
void calculateVelocityOnBoundary(const Intersection &intersection, CellData &cellData)
Calculates the velocity at a boundary.
Definition: lmethod/2dpressurevelocity.hh:388
FvMpfaL2dPressureVelocity2p(Problem &problem)
Constructs a FvMpfaL2dPressureVelocity2p object.
Definition: lmethod/2dpressurevelocity.hh:104
void initialize()
Initializes pressure and velocity.
Definition: lmethod/2dpressurevelocity.hh:138
void addOutputVtkFields(MultiWriter &writer)
Adds velocity output to the output file.
Definition: lmethod/2dpressurevelocity.hh:193
bool calculateVelocityInTransport()
Indicates if velocity is reconstructed in the pressure step or in the transport step.
Definition: lmethod/2dpressurevelocity.hh:177
void calculateVelocity()
Calculates the velocities at all cell-cell interfaces.
Definition: lmethod/2dpressurevelocity.hh:221
void update()
Pressure and velocity update.
Definition: lmethod/2dpressurevelocity.hh:165
Class for calculating 2-d velocities from cell-wise constant pressure values.
Definition: lmethod/2dvelocity.hh:58
Class including the information of an interaction volume of a MPFA L-method that does not change with...
Definition: linteractionvolume.hh:41
bool isInnerVolume()
Returns true if the interaction volume is completely inside the model domain.
Definition: linteractionvolume.hh:316
Element getSubVolumeElement(int subVolumeIdx)
Get an element of the interaction volume.
Definition: linteractionvolume.hh:281
int getFaceIndexFromSubVolume(int subVolumeIdx, int subVolumeFaceIdx)
Map from local interaction volume numbering on element to numbering on interaction volume.
Definition: linteractionvolume.hh:270
bool isOutsideFace(int subVolumeFaceIdx)
Returns true if an interaction volume flux face is outside the model domain.
Definition: linteractionvolume.hh:301
The finite volume base class for the solution of a pressure equation.
Definition: sequential/cellcentered/pressure.hh:48
Finite volume MPFA O-method discretization of a two-phase pressure equation of the sequential IMPES m...
Velocity calculation using a 2-d MPFA O-method.