3.1-git
DUNE for Multi-{Phase, Component, Scale, Physics, ...} flow and transport in porous media
nonequilibrium.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 *****************************************************************************/
26#ifndef DUMUX_NONEQUILIBRIUM_FLUID_STATE_HH
27#define DUMUX_NONEQUILIBRIUM_FLUID_STATE_HH
28
29#include <cmath>
30#include <algorithm>
31#include <dune/common/exceptions.hh>
32
33namespace Dumux {
34
41template <class ScalarType, class FluidSystem>
43{
44public:
45 static constexpr int numPhases = FluidSystem::numPhases;
46 static constexpr int numComponents = FluidSystem::numComponents;
47
49 using Scalar = ScalarType;
50
51 /*****************************************************
52 * Generic access to fluid properties (No assumptions
53 * on thermodynamic equilibrium required)
54 *****************************************************/
59 int wettingPhase() const { return wPhaseIdx_; }
68 Scalar saturation(int phaseIdx) const
69 { return saturation_[phaseIdx]; }
70
80 Scalar moleFraction(int phaseIdx, int compIdx) const
81 { return moleFraction_[phaseIdx][compIdx]; }
82
83
99 Scalar massFraction(int phaseIdx, int compIdx) const
100 {
101 using std::abs;
102 using std::max;
103 return abs(sumMoleFractions_[phaseIdx])
104 * moleFraction_[phaseIdx][compIdx]
105 * FluidSystem::molarMass(compIdx)
106 / max(1e-40, abs(averageMolarMass_[phaseIdx]));
107 }
108
117 Scalar averageMolarMass(int phaseIdx) const
118 { return averageMolarMass_[phaseIdx]; }
119
129 Scalar molarity(int phaseIdx, int compIdx) const
130 { return molarDensity(phaseIdx)*moleFraction(phaseIdx, compIdx); }
131
135 Scalar fugacityCoefficient(int phaseIdx, int compIdx) const
136 { return fugacityCoefficient_[phaseIdx][compIdx]; }
137
155 Scalar fugacity(int phaseIdx, int compIdx) const
156 {
157 return pressure_[phaseIdx]
158 *fugacityCoefficient_[phaseIdx][compIdx]
159 *moleFraction_[phaseIdx][compIdx];
160 }
161
162 Scalar fugacity(int compIdx) const
163 { return fugacity(0, compIdx); }
164
170 Scalar molarVolume(int phaseIdx) const
171 { return 1.0/molarDensity(phaseIdx); }
172
177 Scalar density(int phaseIdx) const
178 { return density_[phaseIdx]; }
179
188 Scalar molarDensity(int phaseIdx) const
189 { return molarDensity_[phaseIdx]; }
190
194 Scalar temperature(const int phaseIdx) const
195 { return temperature_[phaseIdx]; }
196
201 { return temperatureEquil_ ; }
202
206 Scalar pressure(int phaseIdx) const
207 { return pressure_[phaseIdx]; }
208
213 Scalar partialPressure(int phaseIdx, int compIdx) const
214 {
215 assert(FluidSystem::isGas(phaseIdx));
216 return moleFraction(phaseIdx, compIdx) * pressure(phaseIdx);
217 }
218
222 Scalar enthalpy(int phaseIdx) const
223 { return enthalpy_[phaseIdx]; }
224
232 Scalar internalEnergy(int phaseIdx) const
233 { return enthalpy_[phaseIdx] - pressure(phaseIdx)/density(phaseIdx); }
234
238 Scalar viscosity(int phaseIdx) const
239 { return viscosity_[phaseIdx]; }
240
241 /*****************************************************
242 * Setter methods. Note that these are not part of the
243 * generic FluidState interface but specific for each
244 * implementation...
245 *****************************************************/
249 void setTemperature(int phaseIdx, Scalar value)
250 { temperature_[phaseIdx] = value; }
251
256 { temperatureEquil_ = value; }
257
261 void setPressure(int phaseIdx, Scalar value)
262 { pressure_[phaseIdx] = value; }
263
267 void setSaturation(int phaseIdx, Scalar value)
268 { saturation_[phaseIdx] = value; }
269
273 void setMoleFraction(int phaseIdx, int compIdx, Scalar value)
274 {
275 using std::isfinite;
276 if (isfinite(averageMolarMass_[phaseIdx]))
277 {
278 Scalar delta = value - moleFraction_[phaseIdx][compIdx];
279
280 moleFraction_[phaseIdx][compIdx] = value;
281
282 sumMoleFractions_[phaseIdx] += delta;
283 averageMolarMass_[phaseIdx] += delta*FluidSystem::molarMass(compIdx);
284 }
285 else
286 {
287 moleFraction_[phaseIdx][compIdx] = value;
288
289 // re-calculate the mean molar mass
290 sumMoleFractions_[phaseIdx] = 0.0;
291 averageMolarMass_[phaseIdx] = 0.0;
292 for (int compJIdx = 0; compJIdx < numComponents; ++compJIdx) {
293 sumMoleFractions_[phaseIdx] += moleFraction_[phaseIdx][compJIdx];
294 averageMolarMass_[phaseIdx] += moleFraction_[phaseIdx][compJIdx]*FluidSystem::molarMass(compJIdx);
295 }
296 }
297 }
298
304 void setMassFraction(int phaseIdx, int compIdx, Scalar value)
305 {
306 if (numComponents != 2)
307 DUNE_THROW(Dune::NotImplemented, "This currently only works for 2 components.");
308 else
309 {
310 // calculate average molar mass of the gas phase
311 Scalar M1 = FluidSystem::molarMass(compIdx);
312 Scalar M2 = FluidSystem::molarMass(1-compIdx);
313 Scalar X2 = 1.0-value;
314 Scalar avgMolarMass = M1*M2/(M2 + X2*(M1 - M2));
315
316 moleFraction_[phaseIdx][compIdx] = value * avgMolarMass / M1;
317 moleFraction_[phaseIdx][1-compIdx] = 1.0-moleFraction_[phaseIdx][compIdx];
318
319 // re-calculate the mean molar mass
320 sumMoleFractions_[phaseIdx] = 0.0;
321 averageMolarMass_[phaseIdx] = 0.0;
322 for (int compJIdx = 0; compJIdx < numComponents; ++compJIdx) {
323 sumMoleFractions_[phaseIdx] += moleFraction_[phaseIdx][compJIdx];
324 averageMolarMass_[phaseIdx] += moleFraction_[phaseIdx][compJIdx]*FluidSystem::molarMass(compJIdx);
325 }
326 }
327 }
328
332 void setFugacityCoefficient(int phaseIdx, int compIdx, Scalar value)
333 { fugacityCoefficient_[phaseIdx][compIdx] = value; }
334
338 void setDensity(int phaseIdx, Scalar value)
339 { density_[phaseIdx] = value; }
340
344 void setMolarDensity(int phaseIdx, Scalar value)
345 { molarDensity_[phaseIdx] = value; }
346
350 void setEnthalpy(int phaseIdx, Scalar value)
351 { enthalpy_[phaseIdx] = value; }
352
356 void setViscosity(int phaseIdx, Scalar value)
357 { viscosity_[phaseIdx] = value; }
358
362 void setWettingPhase(int phaseIdx)
363 { wPhaseIdx_ = phaseIdx; }
364
370 template <class FluidState>
371 void assign(const FluidState& fs)
372 {
373 for (int phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx)
374 {
375 averageMolarMass_[phaseIdx] = 0;
376 sumMoleFractions_[phaseIdx] = 0;
377 for (int compIdx = 0; compIdx < numComponents; ++compIdx)
378 {
379 moleFraction_[phaseIdx][compIdx] = fs.moleFraction(phaseIdx, compIdx);
380 fugacityCoefficient_[phaseIdx][compIdx] = fs.fugacityCoefficient(phaseIdx, compIdx);
381 averageMolarMass_[phaseIdx] += moleFraction_[phaseIdx][compIdx]*FluidSystem::molarMass(compIdx);
382 sumMoleFractions_[phaseIdx] += moleFraction_[phaseIdx][compIdx];
383 }
384 pressure_[phaseIdx] = fs.pressure(phaseIdx);
385 saturation_[phaseIdx] = fs.saturation(phaseIdx);
386 density_[phaseIdx] = fs.density(phaseIdx);
387 molarDensity_[phaseIdx] = fs.molarDensity(phaseIdx);
388 enthalpy_[phaseIdx] = fs.enthalpy(phaseIdx);
389 viscosity_[phaseIdx] = fs.viscosity(phaseIdx);
390 temperature_[phaseIdx] = fs.temperature(phaseIdx);
391 }
392 }
393
394protected:
407
408 // For porous medium flow models, here we ... the index of the wetting
409 // phase (needed for vapor pressure evaluation if kelvin equation is used)
411};
412
413} // end namespace Dumux
414
415#endif
make the local view function available whenever we use the grid geometry
Definition: adapt.hh:29
Represents all relevant thermodynamic quantities of a multi-phase, multi-component fluid system witho...
Definition: nonequilibrium.hh:43
Scalar temperatureEquil_
Definition: nonequilibrium.hh:406
ScalarType Scalar
export the scalar type
Definition: nonequilibrium.hh:49
void setPressure(int phaseIdx, Scalar value)
Set the fluid pressure of a phase .
Definition: nonequilibrium.hh:261
Scalar viscosity(int phaseIdx) const
The dynamic viscosity of fluid phase in .
Definition: nonequilibrium.hh:238
Scalar molarDensity_[numPhases]
Definition: nonequilibrium.hh:402
void setMassFraction(int phaseIdx, int compIdx, Scalar value)
Set the mass fraction of a component in a phase and update the average molar mass according to the ...
Definition: nonequilibrium.hh:304
Scalar density(int phaseIdx) const
The mass density of the fluid phase in .
Definition: nonequilibrium.hh:177
Scalar fugacity(int compIdx) const
Definition: nonequilibrium.hh:162
void setDensity(int phaseIdx, Scalar value)
Set the density of a phase .
Definition: nonequilibrium.hh:338
static constexpr int numPhases
Definition: nonequilibrium.hh:45
void setEnthalpy(int phaseIdx, Scalar value)
Set the specific enthalpy of a phase .
Definition: nonequilibrium.hh:350
static constexpr int numComponents
Definition: nonequilibrium.hh:46
int wettingPhase() const
Returns the index of the wetting phase in the fluid-solid configuration (for porous medium systems).
Definition: nonequilibrium.hh:59
void setMolarDensity(int phaseIdx, Scalar value)
Set the molar density of a phase .
Definition: nonequilibrium.hh:344
Scalar pressure(int phaseIdx) const
The pressure of a fluid phase in .
Definition: nonequilibrium.hh:206
Scalar fugacityCoefficient(int phaseIdx, int compIdx) const
The fugacity coefficient of component in fluid phase in .
Definition: nonequilibrium.hh:135
void setTemperature(int phaseIdx, Scalar value)
Set the temperature of a fluid phase.
Definition: nonequilibrium.hh:249
Scalar averageMolarMass(int phaseIdx) const
The average molar mass of phase in .
Definition: nonequilibrium.hh:117
Scalar moleFraction_[numPhases][numComponents]
Definition: nonequilibrium.hh:395
void assign(const FluidState &fs)
Retrieve all parameters from an arbitrary fluid state.
Definition: nonequilibrium.hh:371
void setFugacityCoefficient(int phaseIdx, int compIdx, Scalar value)
Set the fugacity of a component in a phase .
Definition: nonequilibrium.hh:332
void setMoleFraction(int phaseIdx, int compIdx, Scalar value)
Set the mole fraction of a component in a phase .
Definition: nonequilibrium.hh:273
Scalar internalEnergy(int phaseIdx) const
The specific internal energy of a fluid phase in .
Definition: nonequilibrium.hh:232
Scalar temperature(const int phaseIdx) const
The absolute temperature of a fluid phase in .
Definition: nonequilibrium.hh:194
Scalar partialPressure(int phaseIdx, int compIdx) const
The partial pressure of a component in a phase .
Definition: nonequilibrium.hh:213
Scalar pressure_[numPhases]
Definition: nonequilibrium.hh:399
Scalar molarDensity(int phaseIdx) const
The molar density of a fluid phase in .
Definition: nonequilibrium.hh:188
Scalar saturation_[numPhases]
Definition: nonequilibrium.hh:400
Scalar temperature_[numPhases]
Definition: nonequilibrium.hh:405
Scalar molarVolume(int phaseIdx) const
The molar volume of a fluid phase in .
Definition: nonequilibrium.hh:170
Scalar sumMoleFractions_[numPhases]
Definition: nonequilibrium.hh:398
Scalar enthalpy_[numPhases]
Definition: nonequilibrium.hh:403
Scalar density_[numPhases]
Definition: nonequilibrium.hh:401
Scalar fugacity(int phaseIdx, int compIdx) const
The fugacity of component in fluid phase in .
Definition: nonequilibrium.hh:155
void setSaturation(int phaseIdx, Scalar value)
Set the saturation of a phase .
Definition: nonequilibrium.hh:267
void setTemperature(Scalar value)
Set the temperature of all fluid phases.
Definition: nonequilibrium.hh:255
Scalar temperature() const
Get the equilibrium temperature of the fluid phases.
Definition: nonequilibrium.hh:200
Scalar saturation(int phaseIdx) const
Returns the saturation of a fluid phase in .
Definition: nonequilibrium.hh:68
Scalar enthalpy(int phaseIdx) const
The specific enthalpy of a fluid phase in .
Definition: nonequilibrium.hh:222
void setWettingPhase(int phaseIdx)
Set the index of the wetting phase.
Definition: nonequilibrium.hh:362
Scalar averageMolarMass_[numPhases]
Definition: nonequilibrium.hh:397
Scalar molarity(int phaseIdx, int compIdx) const
The molar concentration of component in fluid phase in .
Definition: nonequilibrium.hh:129
int wPhaseIdx_
Definition: nonequilibrium.hh:410
void setViscosity(int phaseIdx, Scalar value)
Set the dynamic viscosity of a phase .
Definition: nonequilibrium.hh:356
Scalar moleFraction(int phaseIdx, int compIdx) const
Returns the molar fraction of the component in fluid phase in .
Definition: nonequilibrium.hh:80
Scalar viscosity_[numPhases]
Definition: nonequilibrium.hh:404
Scalar massFraction(int phaseIdx, int compIdx) const
Returns the mass fraction of component in fluid phase in .
Definition: nonequilibrium.hh:99
Scalar fugacityCoefficient_[numPhases][numComponents]
Definition: nonequilibrium.hh:396