version 3.10.0
Loading...
Searching...
No Matches
h2oairmesitylene.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// SPDX-FileCopyrightText: Copyright © DuMux Project contributors, see AUTHORS.md in root folder
5// SPDX-License-Identifier: GPL-3.0-or-later
6//
12#ifndef DUMUX_H2O_AIR_MESITYLENE_FLUID_SYSTEM_HH
13#define DUMUX_H2O_AIR_MESITYLENE_FLUID_SYSTEM_HH
14
21
25
27
28#include <dumux/io/name.hh>
29
30namespace Dumux::FluidSystems {
31
40template <class Scalar,
41 class H2OType = Components::TabulatedComponent<Components::H2O<Scalar> > >
43 : public Base<Scalar, H2OAirMesitylene<Scalar, H2OType> >
44{
45 using ThisType = H2OAirMesitylene<Scalar, H2OType>;
46
47public:
50 using H2O = H2OType;
51
52
53 static const int numPhases = 3;
54 static const int numComponents = 3;
55
56 static const int wPhaseIdx = 0; // index of the water phase
57 static const int nPhaseIdx = 1; // index of the NAPL phase
58 static const int gPhaseIdx = 2; // index of the gas phase
59
60 static const int H2OIdx = 0;
61 static const int NAPLIdx = 1;
62 static const int AirIdx = 2;
63
64 // export component indices to indicate the main component
65 // of the corresponding phase at atmospheric pressure 1 bar
66 // and room temperature 20°C:
67 static const int wCompIdx = H2OIdx;
68 static const int nCompIdx = NAPLIdx;
69 static const int gCompIdx = AirIdx;
70
77 static void init()
78 {
79 init(/*tempMin=*/273.15,
80 /*tempMax=*/623.15,
81 /*numTemp=*/100,
82 /*pMin=*/0.0,
83 /*pMax=*/20e6,
84 /*numP=*/200);
85 }
86
98 static void init(Scalar tempMin, Scalar tempMax, unsigned nTemp,
99 Scalar pressMin, Scalar pressMax, unsigned nPress)
100 {
101 if (H2O::isTabulated)
102 {
103 H2O::init(tempMin, tempMax, nTemp,
104 pressMin, pressMax, nPress);
105 }
106 }
107
111 static constexpr bool isMiscible()
112 { return true; }
113
119 static constexpr bool isGas(int phaseIdx)
120 {
121 assert(0 <= phaseIdx && phaseIdx < numPhases);
122 return phaseIdx == gPhaseIdx;
123 }
124
131 static bool isIdealGas(int phaseIdx)
132 { return phaseIdx == gPhaseIdx && H2O::gasIsIdeal() && Air::gasIsIdeal() && NAPL::gasIsIdeal(); }
133
148 static bool isIdealMixture(int phaseIdx)
149 {
150 assert(0 <= phaseIdx && phaseIdx < numPhases);
151 // we assume Henry's and Raoult's laws for the water phase and
152 // and no interaction between gas molecules of different
153 // components, so all phases are ideal mixtures!
154 return true;
155 }
156
166 static constexpr bool isCompressible(int phaseIdx)
167 {
168 assert(0 <= phaseIdx && phaseIdx < numPhases);
169 // gases are always compressible
170 if (phaseIdx == gPhaseIdx)
171 return true;
172 else if (phaseIdx == wPhaseIdx)
173 // the water component decides for the water phase...
174 return H2O::liquidIsCompressible();
175
176 // the NAPL component decides for the napl phase...
178 }
179
183 static std::string phaseName(int phaseIdx)
184 {
185 assert(0 <= phaseIdx && phaseIdx < numPhases);
186 switch (phaseIdx)
187 {
188 case wPhaseIdx: return IOName::aqueousPhase();
189 case nPhaseIdx: return IOName::naplPhase();
190 case gPhaseIdx: return IOName::gaseousPhase();
191 }
192 DUNE_THROW(Dune::InvalidStateException, "Invalid phase index " << phaseIdx);
193 }
194
198 static std::string componentName(int compIdx)
199 {
200 switch (compIdx) {
201 case H2OIdx: return H2O::name();
202 case AirIdx: return Air::name();
203 case NAPLIdx: return NAPL::name();
204 }
205 DUNE_THROW(Dune::InvalidStateException, "Invalid component index " << compIdx);
206 }
207
212 static Scalar molarMass(int compIdx)
213 {
214 switch (compIdx) {
215 case H2OIdx: return H2O::molarMass();
216 case AirIdx: return Air::molarMass();
217 case NAPLIdx: return NAPL::molarMass();
218 }
219 DUNE_THROW(Dune::InvalidStateException, "Invalid component index " << compIdx);
220 }
221
222 using Base<Scalar, ThisType>::density;
235 template <class FluidState>
236 static Scalar density(const FluidState &fluidState, int phaseIdx)
237 {
238 if (phaseIdx == wPhaseIdx) {
239 // See: Eq. (7) in Class et al. (2002a)
240 // this assumes each dissolved molecule displaces exactly one
241 // water molecule in the liquid
242 return H2O::liquidMolarDensity(fluidState.temperature(phaseIdx), fluidState.pressure(phaseIdx))
243 * (H2O::molarMass()*fluidState.moleFraction(wPhaseIdx, H2OIdx)
244 + Air::molarMass()*fluidState.moleFraction(wPhaseIdx, AirIdx)
245 + NAPL::molarMass()*fluidState.moleFraction(wPhaseIdx, NAPLIdx));
246 }
247 else if (phaseIdx == nPhaseIdx) {
248 // assume pure NAPL for the NAPL phase
249 Scalar pressure = NAPL::liquidIsCompressible()?fluidState.pressure(phaseIdx):1e100;
250 return NAPL::liquidDensity(fluidState.temperature(phaseIdx), pressure);
251 }
252
253 assert (phaseIdx == gPhaseIdx);
254 Scalar pH2O =
255 fluidState.moleFraction(gPhaseIdx, H2OIdx) *
256 fluidState.pressure(gPhaseIdx);
257 Scalar pAir =
258 fluidState.moleFraction(gPhaseIdx, AirIdx) *
259 fluidState.pressure(gPhaseIdx);
260 Scalar pNAPL =
261 fluidState.moleFraction(gPhaseIdx, NAPLIdx) *
262 fluidState.pressure(gPhaseIdx);
263 return H2O::gasDensity(fluidState.temperature(phaseIdx), pH2O)
264 + Air::gasDensity(fluidState.temperature(phaseIdx), pAir)
265 + NAPL::gasDensity(fluidState.temperature(phaseIdx), pNAPL);
266 }
267
268 using Base<Scalar, ThisType>::molarDensity;
270 template <class FluidState>
271 static Scalar molarDensity(const FluidState &fluidState, int phaseIdx)
272 {
273 Scalar temperature = fluidState.temperature(phaseIdx);
274 Scalar pressure = fluidState.pressure(phaseIdx);
275
276 if (phaseIdx == nPhaseIdx)
277 {
278 // assume pure NAPL for the NAPL phase
279 return NAPL::liquidMolarDensity(temperature, pressure);
280 }
281 else if (phaseIdx == wPhaseIdx)
282 {
283 return H2O::liquidMolarDensity(temperature, pressure);
284 }
285 else
286 {
287 return H2O::gasMolarDensity(temperature, fluidState.partialPressure(gPhaseIdx, H2OIdx))
288 + NAPL::gasMolarDensity(temperature, fluidState.partialPressure(gPhaseIdx, NAPLIdx))
289 + Air::gasMolarDensity(temperature, fluidState.partialPressure(gPhaseIdx, AirIdx));
290 }
291 }
292
293 using Base<Scalar, ThisType>::viscosity;
300 template <class FluidState>
301 static Scalar viscosity(const FluidState &fluidState,
302 int phaseIdx)
303 {
304 if (phaseIdx == wPhaseIdx) {
305 // assume pure water viscosity
306 return H2O::liquidViscosity(fluidState.temperature(phaseIdx),
307 fluidState.pressure(phaseIdx));
308 }
309 else if (phaseIdx == nPhaseIdx) {
310 // assume pure NAPL viscosity
311 return NAPL::liquidViscosity(fluidState.temperature(phaseIdx),
312 fluidState.pressure(phaseIdx));
313 }
314
315 assert (phaseIdx == gPhaseIdx);
316
317 /* Wilke method. See:
318 *
319 * See: R. Reid, et al.: The Properties of Gases and Liquids,
320 * 4th edition, McGraw-Hill, 1987, 407-410
321 * 5th edition, McGraw-Hill, 2001, p. 9.21/22
322 *
323 * in this case, we use a simplified version in order to avoid
324 * computationally costly evaluation of sqrt and pow functions and
325 * divisions
326 * -- compare e.g. with Promo Class p. 32/33
327 */
328 Scalar muResult;
329 const Scalar mu[numComponents] = {
330 h2oGasViscosityInMixture(fluidState.temperature(phaseIdx), fluidState.pressure(phaseIdx)),
331 Air::gasViscosity(fluidState.temperature(phaseIdx), fluidState.pressure(phaseIdx)),
332 NAPL::gasViscosity(fluidState.temperature(phaseIdx), NAPL::vaporPressure(fluidState.temperature(phaseIdx)))
333 };
334 // molar masses
335 const Scalar M[numComponents] = {
336 H2O::molarMass(),
339 };
340
341 Scalar muAW = mu[AirIdx]*fluidState.moleFraction(gPhaseIdx, AirIdx)
342 + mu[H2OIdx]*fluidState.moleFraction(gPhaseIdx, H2OIdx)
343 / (fluidState.moleFraction(gPhaseIdx, AirIdx)
344 + fluidState.moleFraction(gPhaseIdx, H2OIdx));
345 Scalar xAW = fluidState.moleFraction(gPhaseIdx, AirIdx)
346 + fluidState.moleFraction(gPhaseIdx, H2OIdx);
347
348 Scalar MAW = (fluidState.moleFraction(gPhaseIdx, AirIdx)*Air::molarMass()
349 + fluidState.moleFraction(gPhaseIdx, H2OIdx)*H2O::molarMass())
350 / xAW;
351
352 Scalar phiCAW = 0.3; // simplification for this particular system
353 /* actually like this
354 * using std::sqrt;
355 * using std::pow;
356 * Scalar phiCAW = pow(1.+sqrt(mu[NAPLIdx]/muAW)*pow(MAW/M[NAPLIdx],0.25),2)
357 * / sqrt(8.*(1.+M[NAPLIdx]/MAW));
358 */
359 Scalar phiAWC = phiCAW * muAW*M[NAPLIdx]/(mu[NAPLIdx]*MAW);
360
361 muResult = (xAW*muAW)/(xAW+fluidState.moleFraction(gPhaseIdx, NAPLIdx)*phiAWC)
362 + (fluidState.moleFraction(gPhaseIdx, NAPLIdx) * mu[NAPLIdx])
363 / (fluidState.moleFraction(gPhaseIdx, NAPLIdx) + xAW*phiCAW);
364 return muResult;
365 }
366
367
368 using Base<Scalar, ThisType>::diffusionCoefficient;
376 template <class FluidState>
377 static Scalar diffusionCoefficient(const FluidState &fluidState,
378 int phaseIdx,
379 int compIdx)
380 {
381 switch (phaseIdx)
382 {
383 case gPhaseIdx:
384 {
385 switch (compIdx)
386 {
387 case NAPLIdx:
388 {
389 Scalar diffWC = BinaryCoeff::H2O_Mesitylene::gasDiffCoeff(fluidState.temperature(phaseIdx), fluidState.pressure(phaseIdx));
390 Scalar diffAW = BinaryCoeff::H2O_Air::gasDiffCoeff(fluidState.temperature(phaseIdx), fluidState.pressure(phaseIdx));
391 const Scalar xga = fluidState.moleFraction(gPhaseIdx, AirIdx);
392 const Scalar xgw = fluidState.moleFraction(gPhaseIdx, H2OIdx);
393 const Scalar xgc = fluidState.moleFraction(gPhaseIdx, NAPLIdx);
394 return (1.- xgw)/(xga/diffAW + xgc/diffWC);
395 }
396 case H2OIdx:
397 {
398 Scalar diffAC = BinaryCoeff::Air_Mesitylene::gasDiffCoeff(fluidState.temperature(phaseIdx), fluidState.pressure(phaseIdx));
399 Scalar diffWC = BinaryCoeff::H2O_Mesitylene::gasDiffCoeff(fluidState.temperature(phaseIdx), fluidState.pressure(phaseIdx));
400 const Scalar xga = fluidState.moleFraction(gPhaseIdx, AirIdx);
401 const Scalar xgw = fluidState.moleFraction(gPhaseIdx, H2OIdx);
402 const Scalar xgc = fluidState.moleFraction(gPhaseIdx, NAPLIdx);
403 return (1.- xgc)/(xgw/diffWC + xga/diffAC);
404 }
405 case AirIdx:
406 DUNE_THROW(Dune::InvalidStateException, "Diffusivity of Air in the gas phase is constraint by sum of diffusive fluxes = 0 !");
407 }
408 }
409 case wPhaseIdx:
410 {
411 Scalar diffACl = BinaryCoeff::Air_Mesitylene::liquidDiffCoeff(fluidState.temperature(phaseIdx), fluidState.pressure(phaseIdx));
412 Scalar diffWCl = BinaryCoeff::H2O_Mesitylene::liquidDiffCoeff(fluidState.temperature(phaseIdx), fluidState.pressure(phaseIdx));
413 Scalar diffAWl = BinaryCoeff::H2O_Air::liquidDiffCoeff(fluidState.temperature(phaseIdx), fluidState.pressure(phaseIdx));
414
415 Scalar xwa = fluidState.moleFraction(wPhaseIdx, AirIdx);
416 Scalar xww = fluidState.moleFraction(wPhaseIdx, H2OIdx);
417 Scalar xwc = fluidState.moleFraction(wPhaseIdx, NAPLIdx);
418
419 switch (compIdx)
420 {
421 case NAPLIdx:
422 return (1.- xww)/(xwa/diffAWl + xwc/diffWCl);
423 case AirIdx:
424 return (1.- xwc)/(xww/diffWCl + xwa/diffACl);
425 case H2OIdx:
426 DUNE_THROW(Dune::InvalidStateException,
427 "Diffusivity of water in the water phase "
428 "is constraint by sum of diffusive fluxes = 0 !\n");
429 }
430 }
431 case nPhaseIdx:
432 {
433 return 0;
434 }
435 }
436 return 0;
437 }
438
441 template <class FluidState>
442 static Scalar binaryDiffusionCoefficient(const FluidState &fluidState,
443 int phaseIdx,
444 int compIIdx,
445 int compJIdx)
446 {
447 DUNE_THROW(Dune::NotImplemented, "FluidSystems::H2OAirMesitylene::binaryDiffusionCoefficient()");
448 }
449
450 using Base<Scalar, ThisType>::fugacityCoefficient;
464 template <class FluidState>
465 static Scalar fugacityCoefficient(const FluidState &fluidState,
466 int phaseIdx,
467 int compIdx)
468 {
469 assert(0 <= phaseIdx && phaseIdx < numPhases);
470 assert(0 <= compIdx && compIdx < numComponents);
471
472 Scalar T = fluidState.temperature(phaseIdx);
473 Scalar p = fluidState.pressure(phaseIdx);
474
475 if (phaseIdx == wPhaseIdx) {
476 if (compIdx == H2OIdx)
477 return H2O::vaporPressure(T)/p;
478 else if (compIdx == AirIdx)
479 return BinaryCoeff::H2O_Air::henry(T)/p;
480 else if (compIdx == NAPLIdx)
482 }
483
484 // for the NAPL phase, we assume currently that nothing is
485 // dissolved. this means that the affinity of the NAPL
486 // component to the NAPL phase is much higher than for the
487 // other components, i.e. the fugacity coefficient is much
488 // smaller.
489 if (phaseIdx == nPhaseIdx) {
490 Scalar phiNapl = NAPL::vaporPressure(T)/p;
491 if (compIdx == NAPLIdx)
492 return phiNapl;
493 else if (compIdx == AirIdx)
494 return 1e6*phiNapl;
495 else if (compIdx == H2OIdx)
496 return 1e6*phiNapl;
497 }
498
499 // for the gas phase, assume an ideal gas when it comes to
500 // fugacity (-> fugacity == partial pressure)
501 assert(phaseIdx == gPhaseIdx);
502 return 1.0;
503 }
504
505 template <class FluidState>
506 static Scalar kelvinVaporPressure(const FluidState &fluidState,
507 const int phaseIdx,
508 const int compIdx)
509 {
510 DUNE_THROW(Dune::NotImplemented, "FluidSystems::H2OAirMesitylene::kelvinVaporPressure()");
511 }
512
513 using Base<Scalar, ThisType>::enthalpy;
523 template <class FluidState>
524 static Scalar enthalpy(const FluidState &fluidState,
525 int phaseIdx)
526 {
527 if (phaseIdx == wPhaseIdx) {
528 return H2O::liquidEnthalpy(fluidState.temperature(phaseIdx), fluidState.pressure(phaseIdx));
529 }
530 else if (phaseIdx == nPhaseIdx) {
531 return NAPL::liquidEnthalpy(fluidState.temperature(phaseIdx), fluidState.pressure(phaseIdx));
532 }
533 else if (phaseIdx == gPhaseIdx) { // gas phase enthalpy depends strongly on composition
534 Scalar hgc = NAPL::gasEnthalpy(fluidState.temperature(phaseIdx),
535 fluidState.pressure(phaseIdx));
536 Scalar hgw = H2O::gasEnthalpy(fluidState.temperature(phaseIdx),
537 fluidState.pressure(phaseIdx));
538 // pressure is only a dummy here (not dependent on pressure, just temperature)
539 Scalar hga = Air::gasEnthalpy(fluidState.temperature(phaseIdx), fluidState.pressure(phaseIdx));
540
541 Scalar result = 0;
542 result += hgw * fluidState.massFraction(gPhaseIdx, H2OIdx);
543 result += hga * fluidState.massFraction(gPhaseIdx, AirIdx);
544 result += hgc * fluidState.massFraction(gPhaseIdx, NAPLIdx);
545
546 return result;
547 }
548 DUNE_THROW(Dune::InvalidStateException, "Invalid phase index " << phaseIdx);
549 }
550
557 template <class FluidState>
558 static Scalar componentEnthalpy(const FluidState& fluidState, int phaseIdx, int componentIdx)
559 {
560 const Scalar T = fluidState.temperature(phaseIdx);
561 const Scalar p = fluidState.pressure(phaseIdx);
562
563 if (phaseIdx == wPhaseIdx)
564 {
565 if (componentIdx == H2OIdx)
566 return H2O::liquidEnthalpy(T, p);
567 else if (componentIdx == NAPLIdx)
568 DUNE_THROW(Dune::NotImplemented, "The component enthalpy for NAPL in water is not implemented.");
569 else if (componentIdx == AirIdx)
570 DUNE_THROW(Dune::NotImplemented, "The component enthalpy for Air in water is not implemented.");
571 DUNE_THROW(Dune::InvalidStateException, "Invalid component index " << componentIdx);
572 }
573 else if (phaseIdx == nPhaseIdx)
574 {
575 if (componentIdx == H2OIdx)
576 DUNE_THROW(Dune::NotImplemented, "The component enthalpy for water in NAPL is not implemented.");
577 else if (componentIdx == NAPLIdx)
578 return NAPL::liquidEnthalpy(T, p);
579 else if (componentIdx == AirIdx)
580 DUNE_THROW(Dune::NotImplemented, "The component enthalpy for air in NAPL is not implemented.");
581 DUNE_THROW(Dune::InvalidStateException, "Invalid component index " << componentIdx);
582 }
583 else if (phaseIdx == gPhaseIdx)
584 {
585 if (componentIdx == H2OIdx)
586 return H2O::gasEnthalpy(T, p);
587 else if (componentIdx == NAPLIdx)
588 return NAPL::gasEnthalpy(T, p);
589 else if (componentIdx == AirIdx)
590 return Air::gasEnthalpy(T,p);
591 DUNE_THROW(Dune::InvalidStateException, "Invalid component index " << componentIdx);
592 }
593 DUNE_THROW(Dune::InvalidStateException, "Invalid phase index " << phaseIdx);
594 }
595
596 using Base<Scalar, ThisType>::heatCapacity;
598 template <class FluidState>
599 static Scalar heatCapacity(const FluidState &fluidState,
600 int phaseIdx)
601 {
602 DUNE_THROW(Dune::NotImplemented, "FluidSystems::H2OAirMesitylene::heatCapacity()");
603 }
604
605 using Base<Scalar, ThisType>::thermalConductivity;
607 template <class FluidState>
608 static Scalar thermalConductivity(const FluidState &fluidState,
609 int phaseIdx)
610 {
611 const Scalar temperature = fluidState.temperature(phaseIdx) ;
612 const Scalar pressure = fluidState.pressure(phaseIdx);
613 if (phaseIdx == wPhaseIdx)
614 {
615 return H2O::liquidThermalConductivity(temperature, pressure);
616 }
617 else if (phaseIdx == nPhaseIdx)
618 {
619 return NAPL::liquidThermalConductivity(temperature, pressure);
620 }
621 else if (phaseIdx == gPhaseIdx)
622 {
623 return Air::gasThermalConductivity(temperature, pressure);
624 }
625 DUNE_THROW(Dune::InvalidStateException, "Invalid phase index " << phaseIdx);
626 }
627
628private:
629 static Scalar waterPhaseDensity_(Scalar T, Scalar pw, Scalar xww, Scalar xwa, Scalar xwc)
630 {
631 Scalar rholH2O = H2O::liquidDensity(T, pw);
632 Scalar clH2O = rholH2O/H2O::molarMass();
633
634 // this assumes each dissolved molecule displaces exactly one
635 // water molecule in the liquid
636 return clH2O*(xww*H2O::molarMass() + xwa*Air::molarMass() + xwc*NAPL::molarMass());
637 }
638
639 static Scalar gasPhaseDensity_(Scalar T, Scalar pg, Scalar xgw, Scalar xga, Scalar xgc)
640 {
641 return H2O::gasDensity(T, pg*xgw) + Air::gasDensity(T, pg*xga) + NAPL::gasDensity(T, pg*xgc);
642 }
643
644 static Scalar NAPLPhaseDensity_(Scalar T, Scalar pn)
645 {
646 return NAPL::liquidDensity(T, pn);
647 }
648
649};
650
651} // end namespace Dumux::FluidSystems
652
653#endif
A simple class for the air fluid properties.
Binary coefficients for air and mesitylene.
static Scalar liquidDiffCoeff(Scalar temperature, Scalar pressure)
Diffusion coefficient for air and mesitylene in liquid water.
Definition air_mesitylene.hh:96
static Scalar gasDiffCoeff(Scalar temperature, Scalar pressure)
Binary diffusion coefficient for air and mesitylene. I used the method according to Wilke and Lee se...
Definition air_mesitylene.hh:46
static Scalar henry(Scalar temperature)
Henry coefficient for air in liquid water.
Definition h2o_air.hh:35
static Scalar gasDiffCoeff(Scalar temperature, Scalar pressure)
Binary diffusion coefficient for molecular water and air.
Definition h2o_air.hh:55
static Scalar liquidDiffCoeff(Scalar temperature, Scalar pressure)
Diffusion coefficient for molecular nitrogen in liquid water.
Definition h2o_air.hh:88
static Scalar liquidDiffCoeff(Scalar temperature, Scalar pressure)
Diffusion coefficient for mesitylene in liquid water.
Definition h2o_mesitylene.hh:104
static Scalar henry(Scalar temperature)
Henry coefficient for mesitylene in liquid water.
Definition h2o_mesitylene.hh:38
static Scalar gasDiffCoeff(Scalar temperature, Scalar pressure)
Binary diffusion coefficient for molecular water and mesitylene.
Definition h2o_mesitylene.hh:54
A class for the air fluid properties.
Definition air.hh:33
static Scalar gasDensity(Scalar temperature, Scalar pressure)
The density of Air at a given pressure and temperature.
Definition air.hh:71
static constexpr Scalar molarMass()
The molar mass in of Air.
Definition air.hh:48
static Scalar gasViscosity(Scalar temperature, Scalar pressure)
The dynamic viscosity of Air at a given pressure and temperature.
Definition air.hh:179
static Scalar gasThermalConductivity(Scalar temperature, Scalar pressure)
Thermal conductivity of air.
Definition air.hh:337
static constexpr bool gasIsIdeal()
Returns true, the gas phase is assumed to be ideal.
Definition air.hh:95
static Scalar gasEnthalpy(Scalar temperature, Scalar pressure)
Specific enthalpy of Air with 273.15 as basis.
Definition air.hh:262
static std::string name()
A human readable name for Air.
Definition air.hh:40
static Scalar gasMolarDensity(Scalar temperature, Scalar pressure)
The molar density of air in , depending on pressure and temperature.
Definition air.hh:83
mesitylene
Definition mesitylene.hh:35
static Scalar gasMolarDensity(Scalar temperature, Scalar pressure)
The molar density of mesitylene in , depending on pressure and temperature.
Definition mesitylene.hh:206
static Scalar gasEnthalpy(Scalar temperature, Scalar pressure)
Specific enthalpy of mesitylene vapor .
Definition mesitylene.hh:182
static std::string name()
A human readable name for the mesitylene.
Definition mesitylene.hh:42
static Scalar liquidMolarDensity(Scalar temperature, Scalar pressure)
The molar density of pure mesitylene at a given pressure and temperature .
Definition mesitylene.hh:229
static constexpr Scalar molarMass()
The molar mass in of mesitylene.
Definition mesitylene.hh:48
static Scalar liquidThermalConductivity(Scalar temperature, Scalar pressure)
Thermal conductivity of mesitylene.
Definition mesitylene.hh:365
static Scalar gasViscosity(Scalar temperature, Scalar pressure)
The dynamic viscosity of mesitylene vapor.
Definition mesitylene.hh:269
static Scalar vaporPressure(Scalar temperature)
The saturation vapor pressure in of pure mesitylene at a given temperature according to Antoine afte...
Definition mesitylene.hh:92
static constexpr bool liquidIsCompressible()
Returns true if the liquid phase is assumed to be compressible.
Definition mesitylene.hh:260
static Scalar gasDensity(Scalar temperature, Scalar pressure)
The density of mesitylene at a given pressure and temperature .
Definition mesitylene.hh:193
static Scalar liquidEnthalpy(const Scalar temperature, const Scalar pressure)
Specific enthalpy of liquid mesitylene .
Definition mesitylene.hh:111
static Scalar liquidDensity(Scalar temperature, Scalar pressure)
The density of pure mesitylene at a given pressure and temperature .
Definition mesitylene.hh:215
static constexpr bool gasIsIdeal()
Returns true if the gas phase is assumed to be ideal.
Definition mesitylene.hh:254
static Scalar liquidViscosity(Scalar temperature, Scalar pressure)
The dynamic viscosity of pure mesitylene.
Definition mesitylene.hh:299
Fluid system base class.
Definition fluidsystems/base.hh:32
ScalarType Scalar
Definition fluidsystems/base.hh:35
A three-phase fluid system featuring gas, NAPL and water as phases and distilled water and air (Pseu...
Definition h2oairmesitylene.hh:44
static const int numComponents
Definition h2oairmesitylene.hh:54
static Scalar thermalConductivity(const FluidState &fluidState, int phaseIdx)
Thermal conductivity of a fluid phase .
Definition h2oairmesitylene.hh:608
static Scalar molarDensity(const FluidState &fluidState, int phaseIdx)
Calculate the molar density of a fluid phase.
Definition h2oairmesitylene.hh:271
static void init()
Initialize the fluid system's static parameters generically.
Definition h2oairmesitylene.hh:77
static const int numPhases
Definition h2oairmesitylene.hh:53
static const int nPhaseIdx
Definition h2oairmesitylene.hh:57
static Scalar componentEnthalpy(const FluidState &fluidState, int phaseIdx, int componentIdx)
Returns the specific enthalpy of a component in a specific phase.
Definition h2oairmesitylene.hh:558
Dumux::Components::Air< Scalar > Air
Definition h2oairmesitylene.hh:49
static const int H2OIdx
Definition h2oairmesitylene.hh:60
static const int AirIdx
Definition h2oairmesitylene.hh:62
static Scalar molarMass(int compIdx)
Return the molar mass of a component in .
Definition h2oairmesitylene.hh:212
static Scalar viscosity(const FluidState &fluidState, int phaseIdx)
Return the viscosity of a phase .
Definition h2oairmesitylene.hh:301
static Scalar enthalpy(const FluidState &fluidState, int phaseIdx)
Given all mole fractions in a phase, return the specific phase enthalpy .
Definition h2oairmesitylene.hh:524
static bool isIdealGas(int phaseIdx)
Returns true if and only if a fluid phase is assumed to be an ideal gas.
Definition h2oairmesitylene.hh:131
static Scalar diffusionCoefficient(const FluidState &fluidState, int phaseIdx, int compIdx)
Given all mole fractions, return the diffusion coefficient in of a component in a phase.
Definition h2oairmesitylene.hh:377
Components::Mesitylene< Scalar > NAPL
Definition h2oairmesitylene.hh:48
static const int wPhaseIdx
Definition h2oairmesitylene.hh:56
static constexpr bool isGas(int phaseIdx)
Return whether a phase is gaseous.
Definition h2oairmesitylene.hh:119
static constexpr bool isCompressible(int phaseIdx)
Returns true if and only if a fluid phase is assumed to be compressible.
Definition h2oairmesitylene.hh:166
static const int wCompIdx
Definition h2oairmesitylene.hh:67
static Scalar fugacityCoefficient(const FluidState &fluidState, int phaseIdx, int compIdx)
Returns the fugacity coefficient of a component in a phase.
Definition h2oairmesitylene.hh:465
static const int gPhaseIdx
Definition h2oairmesitylene.hh:58
static void init(Scalar tempMin, Scalar tempMax, unsigned nTemp, Scalar pressMin, Scalar pressMax, unsigned nPress)
Initialize the fluid system's static parameters using problem specific temperature and pressure range...
Definition h2oairmesitylene.hh:98
static constexpr bool isMiscible()
Returns whether the fluids are miscible.
Definition h2oairmesitylene.hh:111
static Scalar heatCapacity(const FluidState &fluidState, int phaseIdx)
Specific isobaric heat capacity of a fluid phase .
Definition h2oairmesitylene.hh:599
static bool isIdealMixture(int phaseIdx)
Returns true if and only if a fluid phase is assumed to be an ideal mixture.
Definition h2oairmesitylene.hh:148
static std::string componentName(int compIdx)
Return the human readable name of a component (used in indices)
Definition h2oairmesitylene.hh:198
static Scalar density(const FluidState &fluidState, int phaseIdx)
Given a phase's composition, temperature, pressure, and the partial pressures of all components,...
Definition h2oairmesitylene.hh:236
static const int NAPLIdx
Definition h2oairmesitylene.hh:61
H2OType H2O
Definition h2oairmesitylene.hh:50
static const int gCompIdx
Definition h2oairmesitylene.hh:69
static const int nCompIdx
Definition h2oairmesitylene.hh:68
static std::string phaseName(int phaseIdx)
Return the human readable name of a phase (used in indices)
Definition h2oairmesitylene.hh:183
static Scalar kelvinVaporPressure(const FluidState &fluidState, const int phaseIdx, const int compIdx)
Definition h2oairmesitylene.hh:506
static Scalar binaryDiffusionCoefficient(const FluidState &fluidState, int phaseIdx, int compIIdx, int compJIdx)
Given a phase's composition, temperature and pressure, return the binary diffusion coefficient for c...
Definition h2oairmesitylene.hh:442
Fluid system base class.
Material properties of pure water .
Binary coefficients for water and air.
Binary coefficients for water and mesitylene.
Relations valid for an ideal gas.
Properties of mesitylene.
A collection of input/output field names for common physical quantities.
Definition h2o.hh:901
Scalar h2oGasViscosityInMixture(Scalar temperature, Scalar pressure)
The dynamic viscosity of steam in a gas mixture.
Definition h2o.hh:964
std::string gaseousPhase() noexcept
I/O name of gaseous phase.
Definition name.hh:111
std::string naplPhase() noexcept
I/O name of napl phase.
Definition name.hh:119
std::string aqueousPhase() noexcept
I/O name of aqueous phase.
Definition name.hh:115
Tabulates all thermodynamic properties of a given untabulated chemical species.