version 3.10.0
Loading...
Searching...
No Matches
h2oairxylene.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_XYLENE_FLUID_SYSTEM_HH
13#define DUMUX_H2O_AIR_XYLENE_FLUID_SYSTEM_HH
14
20
24
26
27#include <dumux/io/name.hh>
28
29namespace Dumux::FluidSystems {
30
39template <class Scalar,
40 class H2OType = Components::TabulatedComponent<Components::H2O<Scalar> > >
42 : public Base<Scalar, H2OAirXylene<Scalar, H2OType> >
43{
44 using ThisType = H2OAirXylene<Scalar, H2OType>;
45
46public:
47 using H2O = H2OType;
50
51 static const int numPhases = 3;
52 static const int numComponents = 3;
53
54 static const int wPhaseIdx = 0; // index of the water phase
55 static const int nPhaseIdx = 1; // index of the NAPL phase
56 static const int gPhaseIdx = 2; // index of the gas phase
57
58 static const int H2OIdx = 0;
59 static const int NAPLIdx = 1;
60 static const int AirIdx = 2;
61
62 // export component indices to indicate the main component
63 // of the corresponding phase at atmospheric pressure 1 bar
64 // and room temperature 20°C:
65 static const int wCompIdx = H2OIdx;
66 static const int nCompIdx = NAPLIdx;
67 static const int gCompIdx = AirIdx;
68
75 static void init()
76 {
77 init(/*tempMin=*/273.15,
78 /*tempMax=*/623.15,
79 /*numTemp=*/100,
80 /*pMin=*/0.0,
81 /*pMax=*/20e6,
82 /*numP=*/200);
83 }
84
96 static void init(Scalar tempMin, Scalar tempMax, unsigned nTemp,
97 Scalar pressMin, Scalar pressMax, unsigned nPress)
98 {
99 if (H2O::isTabulated)
100 {
101 H2O::init(tempMin, tempMax, nTemp,
102 pressMin, pressMax, nPress);
103 }
104 }
105
109 static constexpr bool isMiscible()
110 { return true; }
111
117 static constexpr bool isGas(int phaseIdx)
118 {
119 assert(0 <= phaseIdx && phaseIdx < numPhases);
120 return phaseIdx == gPhaseIdx;
121 }
122
129 static bool isIdealGas(int phaseIdx)
130 { return phaseIdx == gPhaseIdx && H2O::gasIsIdeal() && Air::gasIsIdeal() && NAPL::gasIsIdeal(); }
131
146 static bool isIdealMixture(int phaseIdx)
147 {
148 assert(0 <= phaseIdx && phaseIdx < numPhases);
149 // we assume Henry's and Raoult's laws for the water phase and
150 // and no interaction between gas molecules of different
151 // components, so all phases are ideal mixtures!
152 return true;
153 }
154
164 static constexpr bool isCompressible(int phaseIdx)
165 {
166 assert(0 <= phaseIdx && phaseIdx < numPhases);
167 // gases are always compressible
168 if (phaseIdx == gPhaseIdx)
169 return true;
170 else if (phaseIdx == wPhaseIdx)
171 // the water component decides for the water phase...
172 return H2O::liquidIsCompressible();
173
174 // the NAPL component decides for the napl phase...
176 }
177
182 static std::string phaseName(int phaseIdx)
183 {
184 assert(0 <= phaseIdx && phaseIdx < numPhases);
185 switch (phaseIdx)
186 {
187 case wPhaseIdx: return IOName::aqueousPhase();
188 case nPhaseIdx: return IOName::naplPhase();
189 case gPhaseIdx: return IOName::gaseousPhase();
190 }
191 DUNE_THROW(Dune::InvalidStateException, "Invalid phase index " << phaseIdx);
192 }
193
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 // This assumes each gas molecule displaces exactly one
240 // molecule in the liquid.
241 return H2O::liquidMolarDensity(fluidState.temperature(phaseIdx), fluidState.pressure(phaseIdx))
242 * (H2O::molarMass()*fluidState.moleFraction(wPhaseIdx, H2OIdx)
243 + Air::molarMass()*fluidState.moleFraction(wPhaseIdx, AirIdx)
244 + NAPL::molarMass()*fluidState.moleFraction(wPhaseIdx, NAPLIdx));
245 }
246 else if (phaseIdx == nPhaseIdx) {
247 // assume pure NAPL for the NAPL phase
248 Scalar pressure = NAPL::liquidIsCompressible()?fluidState.pressure(phaseIdx):1e100;
249 return NAPL::liquidDensity(fluidState.temperature(phaseIdx), pressure);
250 }
251
252 assert (phaseIdx == gPhaseIdx);
253 Scalar pH2O =
254 fluidState.moleFraction(gPhaseIdx, H2OIdx) *
255 fluidState.pressure(gPhaseIdx);
256 Scalar pAir =
257 fluidState.moleFraction(gPhaseIdx, AirIdx) *
258 fluidState.pressure(gPhaseIdx);
259 Scalar pNAPL =
260 fluidState.moleFraction(gPhaseIdx, NAPLIdx) *
261 fluidState.pressure(gPhaseIdx);
262 return H2O::gasDensity(fluidState.temperature(phaseIdx), pH2O)
263 + Air::gasDensity(fluidState.temperature(phaseIdx), pAir)
264 + NAPL::gasDensity(fluidState.temperature(phaseIdx), pNAPL);
265 }
266
267 using Base<Scalar, ThisType>::molarDensity;
269 template <class FluidState>
270 static Scalar molarDensity(const FluidState &fluidState, int phaseIdx)
271 {
272 Scalar temperature = fluidState.temperature(phaseIdx);
273 Scalar pressure = fluidState.pressure(phaseIdx);
274 if (phaseIdx == nPhaseIdx)
275 {
276 return NAPL::liquidMolarDensity(temperature, pressure);
277 }
278 else if (phaseIdx == wPhaseIdx)
279 { // This assumes each gas molecule displaces exactly one
280 // molecule in the liquid.
281 return H2O::liquidMolarDensity(temperature, pressure);
282 }
283 else
284 {
285 return H2O::gasMolarDensity(temperature, fluidState.partialPressure(gPhaseIdx, H2OIdx))
286 + NAPL::gasMolarDensity(temperature, fluidState.partialPressure(gPhaseIdx, NAPLIdx))
287 + Air::gasMolarDensity(temperature, fluidState.partialPressure(gPhaseIdx, AirIdx));
288 }
289 }
290
291 using Base<Scalar, ThisType>::viscosity;
299 template <class FluidState>
300 static Scalar viscosity(const FluidState &fluidState,
301 int phaseIdx)
302 {
303 if (phaseIdx == wPhaseIdx) {
304 // assume pure water viscosity
305 return H2O::liquidViscosity(fluidState.temperature(phaseIdx),
306 fluidState.pressure(phaseIdx));
307 }
308 else if (phaseIdx == nPhaseIdx) {
309 // assume pure NAPL viscosity
310 return NAPL::liquidViscosity(fluidState.temperature(phaseIdx),
311 fluidState.pressure(phaseIdx));
312 }
313
314 assert (phaseIdx == gPhaseIdx);
315
316 /* Wilke method. See:
317 *
318 * See: R. Reid, et al.: The Properties of Gases and Liquids,
319 * 4th edition, McGraw-Hill, 1987, 407-410
320 * 5th edition, McGraw-Hill, 20001, p. 9.21/22
321 *
322 * in this case, we use a simplified version in order to avoid
323 * computationally costly evaluation of sqrt and pow functions and
324 * divisions
325 * -- compare e.g. with Promo Class p. 32/33
326 */
327 Scalar muResult;
328 const Scalar mu[numComponents] = {
329 h2oGasViscosityInMixture(fluidState.temperature(phaseIdx),fluidState.pressure(phaseIdx)),
330 Air::gasViscosity(fluidState.temperature(phaseIdx), fluidState.pressure(phaseIdx)),
331 NAPL::gasViscosity(fluidState.temperature(phaseIdx), NAPL::vaporPressure(fluidState.temperature(phaseIdx)))
332 };
333 // molar masses
334 const Scalar M[numComponents] = {
335 H2O::molarMass(),
338 };
339
340 Scalar muAW = mu[AirIdx]*fluidState.moleFraction(gPhaseIdx, AirIdx)
341 + mu[H2OIdx]*fluidState.moleFraction(gPhaseIdx, H2OIdx)
342 / (fluidState.moleFraction(gPhaseIdx, AirIdx)
343 + fluidState.moleFraction(gPhaseIdx, H2OIdx));
344 Scalar xAW = fluidState.moleFraction(gPhaseIdx, AirIdx)
345 + fluidState.moleFraction(gPhaseIdx, H2OIdx);
346
347 Scalar MAW = (fluidState.moleFraction(gPhaseIdx, AirIdx)*Air::molarMass()
348 + fluidState.moleFraction(gPhaseIdx, H2OIdx)*H2O::molarMass())
349 / xAW;
350
351 Scalar phiCAW = 0.3; // simplification for this particular system
352 /* actually like this
353 * using std::sqrt;
354 * using std::pow;
355 * Scalar phiCAW = pow(1.+sqrt(mu[NAPLIdx]/muAW)*pow(MAW/M[NAPLIdx],0.25),2)
356 * / sqrt(8.*(1.+M[NAPLIdx]/MAW));
357 */
358 Scalar phiAWC = phiCAW * muAW*M[NAPLIdx]/(mu[NAPLIdx]*MAW);
359
360 muResult = (xAW*muAW)/(xAW+fluidState.moleFraction(gPhaseIdx, NAPLIdx)*phiAWC)
361 + (fluidState.moleFraction(gPhaseIdx, NAPLIdx) * mu[NAPLIdx])
362 / (fluidState.moleFraction(gPhaseIdx, NAPLIdx) + xAW*phiCAW);
363 return muResult;
364 }
365
366
367 using Base<Scalar, ThisType>::diffusionCoefficient;
369 template <class FluidState>
370 static Scalar diffusionCoefficient(const FluidState &fluidState,
371 int phaseIdx,
372 int compIdx)
373 {
374 Scalar temperature = fluidState.temperature(phaseIdx);
375 Scalar pressure = fluidState.pressure(phaseIdx);
376 if (phaseIdx==gPhaseIdx) {
377 Scalar diffAC = BinaryCoeff::Air_Xylene::gasDiffCoeff(temperature, pressure);
378 Scalar diffWC = BinaryCoeff::H2O_Xylene::gasDiffCoeff(temperature, pressure);
379 Scalar diffAW = BinaryCoeff::H2O_Air::gasDiffCoeff(temperature, pressure);
380
381 const Scalar xga = fluidState.moleFraction(gPhaseIdx, AirIdx);
382 const Scalar xgw = fluidState.moleFraction(gPhaseIdx, H2OIdx);
383 const Scalar xgc = fluidState.moleFraction(gPhaseIdx, NAPLIdx);
384
385 if (compIdx==NAPLIdx)
386 return (1.- xgw)/(xga/diffAW + xgc/diffWC);
387 else if (compIdx==H2OIdx)
388 return (1.- xgc)/(xgw/diffWC + xga/diffAC);
389 else if (compIdx==AirIdx)
390 DUNE_THROW(Dune::InvalidStateException,
391 "Diffusivity of Air in the gas phase "
392 "is constraint by sum of diffusive fluxes = 0 !\n");
393 } else if (phaseIdx==wPhaseIdx){
394 Scalar diffACl = BinaryCoeff::Air_Xylene::liquidDiffCoeff(temperature, pressure);
395 Scalar diffWCl = BinaryCoeff::H2O_Xylene::liquidDiffCoeff(temperature, pressure);
396 Scalar diffAWl = BinaryCoeff::H2O_Air::liquidDiffCoeff(temperature, pressure);
397
398 Scalar xwa = fluidState.moleFraction(wPhaseIdx, AirIdx);
399 Scalar xww = fluidState.moleFraction(wPhaseIdx, H2OIdx);
400 Scalar xwc = fluidState.moleFraction(wPhaseIdx, NAPLIdx);
401
402 Scalar diffCont;
403
404 switch (compIdx) {
405 case NAPLIdx:
406 diffCont = (1.- xww)/(xwa/diffAWl + xwc/diffWCl);
407 return diffCont;
408 case AirIdx:
409 diffCont = (1.- xwc)/(xww/diffWCl + xwa/diffACl);
410 return diffCont;
411 case H2OIdx:
412 DUNE_THROW(Dune::InvalidStateException,
413 "Diffusivity of water in the water phase "
414 "is constraint by sum of diffusive fluxes = 0 !\n");
415 }
416 } else if (phaseIdx==nPhaseIdx) {
417 return 0.0;
418 }
419 return 0;
420 }
421
424 template <class FluidState>
425 static Scalar binaryDiffusionCoefficient(const FluidState &fluidState,
426 int phaseIdx,
427 int compIIdx,
428 int compJIdx)
429 {
430 DUNE_THROW(Dune::NotImplemented, "FluidSystems::H2OAirXylene::binaryDiffusionCoefficient()");
431 }
432
433 using Base<Scalar, ThisType>::fugacityCoefficient;
447 template <class FluidState>
448 static Scalar fugacityCoefficient(const FluidState &fluidState,
449 int phaseIdx,
450 int compIdx)
451 {
452 assert(0 <= phaseIdx && phaseIdx < numPhases);
453 assert(0 <= compIdx && compIdx < numComponents);
454
455 Scalar T = fluidState.temperature(phaseIdx);
456 Scalar p = fluidState.pressure(phaseIdx);
457
458 if (phaseIdx == wPhaseIdx) {
459 if (compIdx == H2OIdx)
460 return H2O::vaporPressure(T)/p;
461 else if (compIdx == AirIdx)
462 return BinaryCoeff::H2O_Air::henry(T)/p;
463 else if (compIdx == NAPLIdx)
465 }
466
467 // for the NAPL phase, we assume currently that nothing is
468 // dissolved. this means that the affinity of the NAPL
469 // component to the NAPL phase is much higher than for the
470 // other components, i.e. the fugacity coefficient is much
471 // smaller.
472 if (phaseIdx == nPhaseIdx) {
473 Scalar phiNapl = NAPL::vaporPressure(T)/p;
474 if (compIdx == NAPLIdx)
475 return phiNapl;
476 else if (compIdx == AirIdx)
477 return 1e6*phiNapl;
478 else if (compIdx == H2OIdx)
479 return 1e6*phiNapl;
480 }
481
482 // for the gas phase, assume an ideal gas when it comes to
483 // fugacity (-> fugacity == partial pressure)
484 assert(phaseIdx == gPhaseIdx);
485 return 1.0;
486 }
487
488 template <class FluidState>
489 static Scalar kelvinVaporPressure(const FluidState &fluidState,
490 const int phaseIdx,
491 const int compIdx)
492 {
493 DUNE_THROW(Dune::NotImplemented, "FluidSystems::H2OAirXylene::kelvinVaporPressure()");
494 }
495
496 using Base<Scalar, ThisType>::enthalpy;
506 template <class FluidState>
507 static Scalar enthalpy(const FluidState &fluidState,
508 int phaseIdx)
509 {
510 if (phaseIdx == wPhaseIdx) {
511 return H2O::liquidEnthalpy(fluidState.temperature(phaseIdx), fluidState.pressure(phaseIdx));
512 }
513 else if (phaseIdx == nPhaseIdx) {
514 return NAPL::liquidEnthalpy(fluidState.temperature(phaseIdx), fluidState.pressure(phaseIdx));
515 }
516 else if (phaseIdx == gPhaseIdx) { // gas phase enthalpy depends strongly on composition
517 Scalar hgc = NAPL::gasEnthalpy(fluidState.temperature(phaseIdx),
518 fluidState.pressure(phaseIdx));
519 Scalar hgw = H2O::gasEnthalpy(fluidState.temperature(phaseIdx),
520 fluidState.pressure(phaseIdx));
521 // pressure is only a dummy here (not dependent on pressure, just temperature)
522 Scalar hga = Air::gasEnthalpy(fluidState.temperature(phaseIdx), fluidState.pressure(phaseIdx));
523
524 Scalar result = 0;
525 result += hgw * fluidState.massFraction(gPhaseIdx, H2OIdx);
526 result += hga * fluidState.massFraction(gPhaseIdx, AirIdx);
527 result += hgc * fluidState.massFraction(gPhaseIdx, NAPLIdx);
528
529 return result;
530 }
531 DUNE_THROW(Dune::InvalidStateException, "Invalid phase index " << phaseIdx);
532 }
533
540 template <class FluidState>
541 static Scalar componentEnthalpy(const FluidState& fluidState, int phaseIdx, int componentIdx)
542 {
543 const Scalar T = fluidState.temperature(phaseIdx);
544 const Scalar p = fluidState.pressure(phaseIdx);
545
546 if (phaseIdx == wPhaseIdx)
547 {
548 if (componentIdx == H2OIdx)
549 return H2O::liquidEnthalpy(T, p);
550 else if (componentIdx == NAPLIdx)
551 DUNE_THROW(Dune::NotImplemented, "The component enthalpy for NAPL in water is not implemented.");
552 else if (componentIdx == AirIdx)
553 DUNE_THROW(Dune::NotImplemented, "The component enthalpy for Air in water is not implemented.");
554 DUNE_THROW(Dune::InvalidStateException, "Invalid component index " << componentIdx);
555 }
556 else if (phaseIdx == nPhaseIdx)
557 {
558 if (componentIdx == H2OIdx)
559 DUNE_THROW(Dune::NotImplemented, "The component enthalpy for water in NAPL is not implemented.");
560 else if (componentIdx == NAPLIdx)
561 return NAPL::liquidEnthalpy(T, p);
562 else if (componentIdx == AirIdx)
563 DUNE_THROW(Dune::NotImplemented, "The component enthalpy for air in NAPL is not implemented.");
564 DUNE_THROW(Dune::InvalidStateException, "Invalid component index " << componentIdx);
565 }
566 else if (phaseIdx == gPhaseIdx)
567 {
568 if (componentIdx == H2OIdx)
569 return H2O::gasEnthalpy(T, p);
570 else if (componentIdx == NAPLIdx)
571 return NAPL::gasEnthalpy(T, p);
572 else if (componentIdx == AirIdx)
573 return Air::gasEnthalpy(T,p);
574 DUNE_THROW(Dune::InvalidStateException, "Invalid component index " << componentIdx);
575 }
576 DUNE_THROW(Dune::InvalidStateException, "Invalid phase index " << phaseIdx);
577 }
578
579 using Base<Scalar, ThisType>::heatCapacity;
581 template <class FluidState>
582 static Scalar heatCapacity(const FluidState &fluidState,
583 int phaseIdx)
584 {
585 DUNE_THROW(Dune::NotImplemented, "FluidSystems::H2OAirXylene::heatCapacity()");
586 }
587
588 using Base<Scalar, ThisType>::thermalConductivity;
590 template <class FluidState>
591 static Scalar thermalConductivity(const FluidState &fluidState,
592 int phaseIdx)
593 {
594 const Scalar temperature = fluidState.temperature(phaseIdx) ;
595 const Scalar pressure = fluidState.pressure(phaseIdx);
596 if (phaseIdx == wPhaseIdx)
597 {
598 return H2O::liquidThermalConductivity(temperature, pressure);
599 }
600 else if (phaseIdx == nPhaseIdx)
601 {
602 return NAPL::liquidThermalConductivity(temperature, pressure);
603 }
604 else if (phaseIdx == gPhaseIdx)
605 {
606 return Air::gasThermalConductivity(temperature, pressure);
607 }
608 DUNE_THROW(Dune::InvalidStateException, "Invalid phase index " << phaseIdx);
609 }
610
611private:
612 static Scalar waterPhaseDensity_(Scalar T, Scalar pw, Scalar xww, Scalar xwa, Scalar xwc)
613 {
614 Scalar rholH2O = H2O::liquidDensity(T, pw);
615 Scalar clH2O = rholH2O/H2O::molarMass();
616
617 // this assumes each dissolved molecule displaces exactly one
618 // water molecule in the liquid
619 return
620 clH2O*(xww*H2O::molarMass() + xwa*Air::molarMass() + xwc*NAPL::molarMass());
621 }
622
623 static Scalar gasPhaseDensity_(Scalar T, Scalar pg, Scalar xgw, Scalar xga, Scalar xgc)
624 {
625 return H2O::gasDensity(T, pg*xgw) + Air::gasDensity(T, pg*xga) + NAPL::gasDensity(T, pg*xgc);
626 }
627
628 static Scalar NAPLPhaseDensity_(Scalar T, Scalar pn)
629 {
630 return
631 NAPL::liquidDensity(T, pn);
632 }
633
634};
635
636} // end namespace Dumux::FluidSystems
637
638#endif
A simple class for the air fluid properties.
Binary coefficients for air and xylene.
static Scalar gasDiffCoeff(Scalar temperature, Scalar pressure)
Binary diffusion coefficient for air and xylene. method according to Wilke and Lee see W....
Definition air_xylene.hh:50
static Scalar liquidDiffCoeff(Scalar temperature, Scalar pressure)
Diffusion coefficient for air and xylene in liquid water.
Definition air_xylene.hh:96
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 xylene in liquid water.
Definition h2o_xylene.hh:105
static Scalar henry(Scalar temperature)
Henry coefficient for xylene in liquid water.
Definition h2o_xylene.hh:40
static Scalar gasDiffCoeff(Scalar temperature, Scalar pressure)
Binary diffusion coefficient for molecular water and xylene.
Definition h2o_xylene.hh:57
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
Properties of xylene.
Definition xylene.hh:36
static Scalar gasViscosity(Scalar temp, Scalar pressure)
The dynamic viscosity of xylene vapor.
Definition xylene.hh:315
static std::string name()
A human readable name for the xylene.
Definition xylene.hh:44
static Scalar liquidDensity(Scalar temperature, Scalar pressure)
The density of pure xylene at a given pressure and temperature .
Definition xylene.hh:286
static Scalar liquidEnthalpy(const Scalar temperature, const Scalar pressure)
Specific enthalpy of liquid xylene .
Definition xylene.hh:154
static Scalar gasMolarDensity(Scalar temperature, Scalar pressure)
The molar gas density of xylene gas at a given pressure and temperature.
Definition xylene.hh:249
static constexpr bool gasIsIdeal()
Returns true if the gas phase is assumed to be ideal.
Definition xylene.hh:300
static Scalar vaporPressure(Scalar temperature)
The saturation vapor pressure in of pure xylene at a given temperature according to Antoine after Be...
Definition xylene.hh:93
static Scalar liquidThermalConductivity(Scalar temperature, Scalar pressure)
Thermal conductivity of xylene.
Definition xylene.hh:371
static constexpr bool liquidIsCompressible()
Returns true if the liquid phase is assumed to be compressible.
Definition xylene.hh:306
static Scalar gasEnthalpy(Scalar temperature, Scalar pressure)
Specific enthalpy of xylene vapor .
Definition xylene.hh:225
static constexpr Scalar molarMass()
The molar mass in of xylene.
Definition xylene.hh:50
static Scalar gasDensity(Scalar temperature, Scalar pressure)
The density of xylene gas at a given pressure and temperature.
Definition xylene.hh:236
static Scalar liquidMolarDensity(Scalar temp, Scalar pressure)
The molar liquid density of pure xylene at a given pressure and temperature .
Definition xylene.hh:261
static Scalar liquidViscosity(Scalar temp, Scalar pressure)
The dynamic viscosity of pure xylene.
Definition xylene.hh:343
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 h2oairxylene.hh:43
static Scalar diffusionCoefficient(const FluidState &fluidState, int phaseIdx, int compIdx)
Calculate the binary molecular diffusion coefficient for a component in a fluid phase .
Definition h2oairxylene.hh:370
static Scalar density(const FluidState &fluidState, int phaseIdx)
Given a phase's composition, temperature, pressure, and the partial pressures of all components,...
Definition h2oairxylene.hh:236
static constexpr bool isGas(int phaseIdx)
Return whether a phase is gaseous.
Definition h2oairxylene.hh:117
static const int gCompIdx
Definition h2oairxylene.hh:67
static Scalar viscosity(const FluidState &fluidState, int phaseIdx)
Return the viscosity of a phase .
Definition h2oairxylene.hh:300
static void init()
Initialize the fluid system's static parameters generically.
Definition h2oairxylene.hh:75
static Scalar heatCapacity(const FluidState &fluidState, int phaseIdx)
Specific isobaric heat capacity of a fluid phase .
Definition h2oairxylene.hh:582
static const int nCompIdx
Definition h2oairxylene.hh:66
static Scalar molarMass(int compIdx)
Return the molar mass of a component in .
Definition h2oairxylene.hh:212
static const int AirIdx
Definition h2oairxylene.hh:60
static const int gPhaseIdx
Definition h2oairxylene.hh:56
static std::string componentName(int compIdx)
Return the human readable name of a component (used in indices)
Definition h2oairxylene.hh:198
static std::string phaseName(int phaseIdx)
Return the human readable name of a phase (used in indices)
Definition h2oairxylene.hh:182
static const int numComponents
Definition h2oairxylene.hh:52
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 h2oairxylene.hh:425
static const int H2OIdx
Definition h2oairxylene.hh:58
static const int wCompIdx
Definition h2oairxylene.hh:65
static constexpr bool isCompressible(int phaseIdx)
Returns true if and only if a fluid phase is assumed to be compressible.
Definition h2oairxylene.hh:164
static constexpr bool isMiscible()
Returns whether the fluids are miscible.
Definition h2oairxylene.hh:109
Components::Xylene< Scalar > NAPL
Definition h2oairxylene.hh:48
static bool isIdealGas(int phaseIdx)
Returns true if and only if a fluid phase is assumed to be an ideal gas.
Definition h2oairxylene.hh:129
H2OType H2O
Definition h2oairxylene.hh:47
static const int nPhaseIdx
Definition h2oairxylene.hh:55
static const int numPhases
Definition h2oairxylene.hh:51
static const int NAPLIdx
Definition h2oairxylene.hh:59
static Scalar fugacityCoefficient(const FluidState &fluidState, int phaseIdx, int compIdx)
Returns the fugacity coefficient of a component in a phase.
Definition h2oairxylene.hh:448
static Scalar molarDensity(const FluidState &fluidState, int phaseIdx)
Calculate the molar density of a fluid phase.
Definition h2oairxylene.hh:270
static Scalar componentEnthalpy(const FluidState &fluidState, int phaseIdx, int componentIdx)
Returns the specific enthalpy of a component in a specific phase.
Definition h2oairxylene.hh:541
static Scalar enthalpy(const FluidState &fluidState, int phaseIdx)
Given all mole fractions in a phase, return the specific phase enthalpy .
Definition h2oairxylene.hh:507
static const int wPhaseIdx
Definition h2oairxylene.hh:54
static Scalar thermalConductivity(const FluidState &fluidState, int phaseIdx)
Thermal conductivity of a fluid phase .
Definition h2oairxylene.hh:591
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 h2oairxylene.hh:96
static bool isIdealMixture(int phaseIdx)
Returns true if and only if a fluid phase is assumed to be an ideal mixture.
Definition h2oairxylene.hh:146
static Scalar kelvinVaporPressure(const FluidState &fluidState, const int phaseIdx, const int compIdx)
Definition h2oairxylene.hh:489
Dumux::Components::Air< Scalar > Air
Definition h2oairxylene.hh:49
Fluid system base class.
Material properties of pure water .
Binary coefficients for water and air.
Binary coefficients for water and xylene.
Relations valid for an ideal gas.
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.
Properties of xylene.