3.3.0
DUNE for Multi-{Phase, Component, Scale, Physics, ...} flow and transport in porous media
2p/sequential/celldata.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_ELEMENTDATA2P_HH
25#define DUMUX_ELEMENTDATA2P_HH
26
27#include "properties.hh"
28#include "fluxdata.hh"
29
30namespace Dumux {
31template<class TypeTag>
32class FluxData2P;
33
46template<class TypeTag, bool enableCompressibility>
48
61template<class TypeTag>
62class CellData2P<TypeTag, false>
63{
64private:
69
70 enum
71 {
72 dim = GridView::dimension, dimWorld = GridView::dimensionworld
73 };
74
76
77 enum
78 {
79 wPhaseIdx = Indices::wPhaseIdx, nPhaseIdx = Indices::nPhaseIdx
80 };
81 enum
82 {
83 numPhases = getPropValue<TypeTag, Properties::NumPhases>()
84 };
85private:
86 Scalar saturation_[numPhases];
87 Scalar pressure_[numPhases];
88 Scalar potential_[numPhases];
89 Scalar mobility_[numPhases];
90 Scalar fracFlowFunc_[numPhases];
91
92 Scalar capillaryPressure_;
93 Scalar update_;
94
95 FluxData fluxData_;
96public:
97
100 {
101 for (int i = 0; i < numPhases;i++)
102 {
103 saturation_[i] = 0.0;
104 pressure_[i] = 0.0;
105 potential_[i] = 0.0;
106 mobility_[i] = 0.0;
107 fracFlowFunc_[i] = 0.0;
108 }
109 capillaryPressure_ = 0.0;
110 update_ = 0.0;
111 }
112
115 {
116 return fluxData_;
117 }
118
120 const FluxData& fluxData() const
121 {
122 return fluxData_;
123 }
124
127 FluidState& fluidState()
128 {
129 DUNE_THROW(Dune::NotImplemented,"fluid states not stored in cell data of incompressible models!");
130 }
131
133 const FluidState& fluidState() const
134 {
135 DUNE_THROW(Dune::NotImplemented,"fluid states not stored in cell data of incompressible models!");
136 }
138
140 // functions returning primary variables
142
148 Scalar pressure(int phaseIdx)
149 {
150 return pressure_[phaseIdx];
151 }
152
158 Scalar pressure(int phaseIdx) const
159 {
160 return pressure_[phaseIdx];
161 }
162
169 void setPressure(int phaseIdx, Scalar press)
170 {
171 pressure_[phaseIdx] = press;
172 }
173
176 {
177 return pressure_[wPhaseIdx];
178 }
179
181 Scalar globalPressure() const
182 {
183 return pressure_[wPhaseIdx];
184 }
185
191 void setGlobalPressure(Scalar press)
192 {
193 pressure_[wPhaseIdx] = press;
194 }
195
201 Scalar potential(int phaseIdx)
202 {
203 return potential_[phaseIdx];
204 }
205
211 Scalar potential(int phaseIdx) const
212 {
213 return potential_[phaseIdx];
214 }
215
222 void setPotential(int phaseIdx, Scalar pot)
223 {
224 potential_[phaseIdx] = pot;
225 }
226
232 Scalar saturation(int phaseIdx)
233 {
234 return saturation_[phaseIdx];
235 }
236
242 Scalar saturation(int phaseIdx) const
243 {
244 return saturation_[phaseIdx];
245 }
246
253 void setSaturation(int phaseIdx, Scalar sat)
254 {
255 saturation_[phaseIdx] = sat;
256 }
257
259 // functions returning the vectors of secondary variables
261
267 Scalar mobility(int phaseIdx)
268 {
269 return mobility_[phaseIdx];
270 }
271
277 Scalar mobility(int phaseIdx) const
278 {
279 return mobility_[phaseIdx];
280 }
281
288 void setMobility(int phaseIdx, Scalar mobility)
289 {
290 mobility_[phaseIdx] = mobility;
291 }
292
298 Scalar fracFlowFunc(int phaseIdx)
299 {
300 return fracFlowFunc_[phaseIdx];
301 }
302
308 Scalar fracFlowFunc(int phaseIdx) const
309 {
310 return fracFlowFunc_[phaseIdx];
311 }
312
319 void setFracFlowFunc(int phaseIdx, Scalar fracFlowFunc)
320 {
321 fracFlowFunc_[phaseIdx] = fracFlowFunc;
322 }
323
326 {
327 return capillaryPressure_;
328 }
329
331 Scalar capillaryPressure() const
332 {
333 return capillaryPressure_;
334 }
335
341 void setCapillaryPressure(Scalar pc)
342 {
343 capillaryPressure_ = pc;
344 }
345
351 void setUpdate(Scalar update)
352 {
353 update_ = update;
354 }
355
358 {
359 return update_;
360 }
362 Scalar volumeCorrection() const
363 {
364 return update_;
365 }
366
368 //densities are not stored for the incompressible model, however the function is needed to avoid compiler errors
369 Scalar density(int phaseIdx)
370 {
371 DUNE_THROW(Dune::NotImplemented,"density function not implemented in cell data of incompressible models!");
372 }
373 Scalar density(int phaseIdx) const
374 {
375 DUNE_THROW(Dune::NotImplemented,"density function not implemented in cell data of incompressible models!");
376 }
377
378 //viscosities are not stored for the incompressible model, however the function is needed to avoid compiler errors
379 Scalar viscosity(int phaseIdx)
380 {
381 DUNE_THROW(Dune::NotImplemented,"viscosity function not implemented in cell data of incompressible models!");
382 }
383
384 Scalar viscosity(int phaseIdx) const
385 {
386 DUNE_THROW(Dune::NotImplemented,"viscosity function not implemented in cell data of incompressible models!");
387 }
389};
390
403template<class TypeTag>
404class CellData2P<TypeTag, true>
405{
406private:
411
412 enum
413 {
414 dim = GridView::dimension, dimWorld = GridView::dimensionworld
415 };
416
418
419 enum
420 {
421 wPhaseIdx = Indices::wPhaseIdx, nPhaseIdx = Indices::nPhaseIdx
422 };
423 enum
424 {
425 numPhases = getPropValue<TypeTag, Properties::NumPhases>()
426 };
427private:
428 Scalar potential_[numPhases];
429 Scalar mobility_[numPhases];
430 Scalar fracFlowFunc_[numPhases];
431 Scalar update_;
432
433 FluxData fluxData_;
434 FluidState fluidState_;
435public:
436
439 {
440 for (int i = 0; i < numPhases;i++)
441 {
442 mobility_[i] = 0.0;
443 potential_[i] = 0.0;
444 fracFlowFunc_[i] = 0.0;
445 }
446 update_ = 0.0;
447 }
448
451 {
452 return fluxData_;
453 }
454
456 const FluxData& fluxData() const
457 {
458 return fluxData_;
459 }
460
462 FluidState& fluidState()
463 {
464 return fluidState_;
465 }
467 const FluidState& fluidState() const
468 {
469 return fluidState_;
470 }
471
473 // functions returning primary variables
475
481 Scalar pressure(int phaseIdx)
482 {
483 return fluidState_.pressure(phaseIdx);
484 }
485
491 Scalar pressure(int phaseIdx) const
492 {
493 return fluidState_.pressure(phaseIdx);
494 }
495
502 void setPressure(int phaseIdx, Scalar press)
503 {
504 fluidState_.setPressure(phaseIdx, press);
505 }
506
512 Scalar potential(int phaseIdx)
513 {
514 return potential_[phaseIdx];
515 }
516
522 Scalar potential(int phaseIdx) const
523 {
524 return potential_[phaseIdx];
525 }
526
533 void setPotential(int phaseIdx, Scalar pot)
534 {
535 potential_[phaseIdx] = pot;
536 }
537
539 // global pressure is not supported in the compressible case, however the function hast to be defined to avoid compiler errors!
540 Scalar globalPressure()
541 {
542 DUNE_THROW(Dune::NotImplemented,"no global pressure defined for compressible models!");
543 }
544
545 // global pressure is not supported in the compressible case, however the function hast to be defined to avoid compiler errors!
546 Scalar globalPressure() const
547 {
548 DUNE_THROW(Dune::NotImplemented,"no global pressure defined for compressible models!");
549 }
550
551 // global pressure is not supported in the compressible case, however the function hast to be defined to avoid compiler errors!
552 void setGlobalPressure(Scalar press)
553 {
554 DUNE_THROW(Dune::NotImplemented,"no global pressure defined for compressible models!");
555 }
557
563 Scalar saturation(int phaseIdx)
564 {
565 return fluidState_.saturation(phaseIdx);
566 }
567
573 Scalar saturation(int phaseIdx) const
574 {
575 return fluidState_.saturation(phaseIdx);
576 }
577
584 void setSaturation(int phaseIdx, Scalar sat)
585 {
586 fluidState_.setSaturation(phaseIdx, sat);
587 }
588
590 // functions returning the vectors of secondary variables
592
598 Scalar mobility(int phaseIdx)
599 {
600 return mobility_[phaseIdx];
601 }
602
608 Scalar mobility(int phaseIdx) const
609 {
610 return mobility_[phaseIdx];
611 }
612
619 void setMobility(int phaseIdx, Scalar mobility)
620 {
621 mobility_[phaseIdx] = mobility;
622 }
623
629 Scalar fracFlowFunc(int phaseIdx)
630 {
631 return fracFlowFunc_[phaseIdx];
632 }
633
639 Scalar fracFlowFunc(int phaseIdx) const
640 {
641 return fracFlowFunc_[phaseIdx];
642 }
643
650 void setFracFlowFunc(int phaseIdx, Scalar fracFlowFunc)
651 {
652 fracFlowFunc_[phaseIdx] = fracFlowFunc;
653 }
654
657 {
658 return fluidState_.pressure(nPhaseIdx) - fluidState_.pressure(wPhaseIdx);
659 }
660
662 Scalar capillaryPressure() const
663 {
664 return fluidState_.pressure(nPhaseIdx) - fluidState_.pressure(wPhaseIdx);
665 }
666
668 //in the compressible model capillary pressure is not stored explicitly, however the function is needed to avoid compiler errors
669 void setCapillaryPressure(Scalar pc)
670 {
671 DUNE_THROW(Dune::NotImplemented,"no capillary pressure stored for compressible models!");
672 }
674
680 Scalar density(int phaseIdx)
681 {
682 return fluidState_.density(phaseIdx);
683 }
684
690 Scalar density(int phaseIdx) const
691 {
692 return fluidState_.density(phaseIdx);
693 }
694
701 void setDensity(int phaseIdx, Scalar density)
702 {
703 fluidState_.setDensity(phaseIdx, density);
704 }
705
711 Scalar viscosity(int phaseIdx)
712 {
713 return fluidState_.viscosity(phaseIdx);
714 }
715
721 Scalar viscosity(int phaseIdx) const
722 {
723 return fluidState_.viscosity(phaseIdx);
724 }
725
732 void setViscosity(int phaseIdx, Scalar viscosity)
733 {
734 fluidState_.setViscosity(phaseIdx, viscosity);
735 }
736
742 void setUpdate(Scalar update)
743 {
744 update_ = update;
745 }
746
749 {
750 return update_;
751 }
752
754 Scalar volumeCorrection() const
755 {
756 return update_;
757 }
758
759};
760} // end namespace Dumux
761#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
std::string mobility(int phaseIdx) noexcept
I/O name of mobility for multiphase systems.
Definition: name.hh:101
std::string viscosity(int phaseIdx) noexcept
I/O name of viscosity for multiphase systems.
Definition: name.hh:74
std::string density(int phaseIdx) noexcept
I/O name of density for multiphase systems.
Definition: name.hh:65
Class storing data assigned to a cell-cell interfaces, so-called flux-data.
Definition: 2p/sequential/fluxdata.hh:41
Class including data of one grid cell.
Definition: 2p/sequential/celldata.hh:47
Scalar globalPressure() const
Returns the global pressure of the cell.
Definition: 2p/sequential/celldata.hh:181
Scalar mobility(int phaseIdx)
Returns the cell phase mobility.
Definition: 2p/sequential/celldata.hh:267
Scalar pressure(int phaseIdx)
Returns the cell phase pressure.
Definition: 2p/sequential/celldata.hh:148
Scalar pressure(int phaseIdx) const
Returns the cell phase pressure.
Definition: 2p/sequential/celldata.hh:158
void setCapillaryPressure(Scalar pc)
Sets the cell capillary pressure.
Definition: 2p/sequential/celldata.hh:341
Scalar mobility(int phaseIdx) const
Returns the cell phase mobility.
Definition: 2p/sequential/celldata.hh:277
Scalar saturation(int phaseIdx) const
Returns the cell phase saturation.
Definition: 2p/sequential/celldata.hh:242
FluxData & fluxData()
Returns the flux data of the cell.
Definition: 2p/sequential/celldata.hh:114
void setMobility(int phaseIdx, Scalar mobility)
Sets the cell phase mobility.
Definition: 2p/sequential/celldata.hh:288
Scalar saturation(int phaseIdx)
Returns the cell phase saturation.
Definition: 2p/sequential/celldata.hh:232
const FluxData & fluxData() const
Returns the flux data of the cell.
Definition: 2p/sequential/celldata.hh:120
Scalar fracFlowFunc(int phaseIdx) const
Returns the cell phase fractional flow function.
Definition: 2p/sequential/celldata.hh:308
void setPressure(int phaseIdx, Scalar press)
Sets the cell phase pressure.
Definition: 2p/sequential/celldata.hh:169
Scalar capillaryPressure() const
Returns the cell capillary pressure.
Definition: 2p/sequential/celldata.hh:331
void setSaturation(int phaseIdx, Scalar sat)
Sets the cell phase saturation.
Definition: 2p/sequential/celldata.hh:253
Scalar volumeCorrection() const
Returns the cell volume correction needed in the pressure equation.
Definition: 2p/sequential/celldata.hh:362
CellData2P()
Constructs a CellData2P object.
Definition: 2p/sequential/celldata.hh:99
Scalar potential(int phaseIdx) const
Returns the cell phase potential.
Definition: 2p/sequential/celldata.hh:211
void setGlobalPressure(Scalar press)
Sets the cell global pressure.
Definition: 2p/sequential/celldata.hh:191
void setUpdate(Scalar update)
Store transport update.
Definition: 2p/sequential/celldata.hh:351
Scalar potential(int phaseIdx)
Returns the cell phase potential.
Definition: 2p/sequential/celldata.hh:201
Scalar volumeCorrection()
Returns the cell volume correction needed in the pressure equation.
Definition: 2p/sequential/celldata.hh:357
void setFracFlowFunc(int phaseIdx, Scalar fracFlowFunc)
Sets the cell phase fractional flow function.
Definition: 2p/sequential/celldata.hh:319
Scalar fracFlowFunc(int phaseIdx)
Returns the cell phase fractional flow function.
Definition: 2p/sequential/celldata.hh:298
void setPotential(int phaseIdx, Scalar pot)
Sets the cell phase potential.
Definition: 2p/sequential/celldata.hh:222
Scalar globalPressure()
Returns the global pressure of the cell.
Definition: 2p/sequential/celldata.hh:175
Scalar capillaryPressure()
Returns the cell capillary pressure.
Definition: 2p/sequential/celldata.hh:325
Scalar capillaryPressure() const
Returns the cell capillary pressure.
Definition: 2p/sequential/celldata.hh:662
void setMobility(int phaseIdx, Scalar mobility)
Sets the cell phase mobility.
Definition: 2p/sequential/celldata.hh:619
void setUpdate(Scalar update)
Store transport update.
Definition: 2p/sequential/celldata.hh:742
Scalar density(int phaseIdx) const
Returns the cell phase density.
Definition: 2p/sequential/celldata.hh:690
Scalar potential(int phaseIdx)
Returns the cell phase potential.
Definition: 2p/sequential/celldata.hh:512
void setPotential(int phaseIdx, Scalar pot)
Sets the cell phase potential.
Definition: 2p/sequential/celldata.hh:533
Scalar potential(int phaseIdx) const
Returns the cell phase potential.
Definition: 2p/sequential/celldata.hh:522
Scalar mobility(int phaseIdx) const
Returns the cell phase mobility.
Definition: 2p/sequential/celldata.hh:608
Scalar capillaryPressure()
Returns the cell capillary pressure.
Definition: 2p/sequential/celldata.hh:656
void setSaturation(int phaseIdx, Scalar sat)
Sets the cell phase saturation.
Definition: 2p/sequential/celldata.hh:584
void setPressure(int phaseIdx, Scalar press)
Sets the cell phase pressure.
Definition: 2p/sequential/celldata.hh:502
Scalar viscosity(int phaseIdx)
Returns the cell phase viscosity.
Definition: 2p/sequential/celldata.hh:711
Scalar density(int phaseIdx)
Returns the cell phase density.
Definition: 2p/sequential/celldata.hh:680
Scalar saturation(int phaseIdx) const
Returns the cell phase saturation.
Definition: 2p/sequential/celldata.hh:573
Scalar fracFlowFunc(int phaseIdx)
Returns the cell phase fractional flow function.
Definition: 2p/sequential/celldata.hh:629
Scalar mobility(int phaseIdx)
Returns the cell phase mobility.
Definition: 2p/sequential/celldata.hh:598
Scalar volumeCorrection() const
Returns the cell volume correction needed in the pressure equation.
Definition: 2p/sequential/celldata.hh:754
Scalar viscosity(int phaseIdx) const
Returns the cell phase viscosity.
Definition: 2p/sequential/celldata.hh:721
Scalar pressure(int phaseIdx) const
Returns the cell phase pressure.
Definition: 2p/sequential/celldata.hh:491
void setDensity(int phaseIdx, Scalar density)
Sets the cell phase density.
Definition: 2p/sequential/celldata.hh:701
Scalar volumeCorrection()
Returns the cell volume correction needed in the pressure equation.
Definition: 2p/sequential/celldata.hh:748
Scalar pressure(int phaseIdx)
Returns the cell phase pressure.
Definition: 2p/sequential/celldata.hh:481
CellData2P()
Constructs a CellData2P object.
Definition: 2p/sequential/celldata.hh:438
const FluxData & fluxData() const
Returns the flux data of the cell.
Definition: 2p/sequential/celldata.hh:456
Scalar fracFlowFunc(int phaseIdx) const
Returns the cell phase fractional flow function.
Definition: 2p/sequential/celldata.hh:639
const FluidState & fluidState() const
Returns the FluidState object for this cell.
Definition: 2p/sequential/celldata.hh:467
void setFracFlowFunc(int phaseIdx, Scalar fracFlowFunc)
Sets the cell phase fractional flow function.
Definition: 2p/sequential/celldata.hh:650
FluxData & fluxData()
Returns the flux data of the cell.
Definition: 2p/sequential/celldata.hh:450
FluidState & fluidState()
Returns the FluidState object for this cell.
Definition: 2p/sequential/celldata.hh:462
void setViscosity(int phaseIdx, Scalar viscosity)
Sets the cell phase viscosity.
Definition: 2p/sequential/celldata.hh:732
Scalar saturation(int phaseIdx)
Returns the cell phase saturation.
Definition: 2p/sequential/celldata.hh:563
Base file for properties related to sequential models.
Class including the variables and data of discretized data of the constitutive relations.