30#ifndef OPM_LBC_CO2RICH_HPP
31#define OPM_LBC_CO2RICH_HPP
38template <
class Scalar,
class Flu
idSystem>
45 template <
class Flu
idState,
class Params,
class LhsEval =
typename Flu
idState::Scalar>
50 const Scalar MPa_atm = 0.101325;
51 const auto& T = Opm::decay<LhsEval>(fluidState.temperature(phaseIdx));
52 const auto& rho = Opm::decay<LhsEval>(fluidState.density(phaseIdx));
55 LhsEval sumVolume = 0.0;
56 for (
unsigned compIdx = 0; compIdx < FluidSystem::numComponents; ++compIdx) {
57 const Scalar& p_c = FluidSystem::criticalPressure(compIdx)/1e6;
58 const Scalar& T_c = FluidSystem::criticalTemperature(compIdx);
59 const Scalar Mm = FluidSystem::molarMass(compIdx) * 1000;
60 const auto& x = Opm::decay<LhsEval>(fluidState.moleFraction(phaseIdx, compIdx));
61 const Scalar v_c = FluidSystem::criticalVolume(compIdx);
66 LhsEval rho_pc = sumMm/sumVolume;
67 LhsEval rho_r = rho/rho_pc;
71 for (
unsigned i_compIdx = 0; i_compIdx < FluidSystem::numComponents; ++i_compIdx) {
72 const Scalar& T_c_i = FluidSystem::criticalTemperature(i_compIdx);
73 const Scalar& p_c_i = FluidSystem::criticalPressure(i_compIdx)/1e6;
74 const auto& x_i = Opm::decay<LhsEval>(fluidState.moleFraction(phaseIdx, i_compIdx));
75 for (
unsigned j_compIdx = 0; j_compIdx < FluidSystem::numComponents; ++j_compIdx) {
76 const Scalar& T_c_j = FluidSystem::criticalTemperature(j_compIdx);
77 const Scalar& p_c_j = FluidSystem::criticalPressure(j_compIdx)/1e6;
78 const auto& x_j = Opm::decay<LhsEval>(fluidState.moleFraction(phaseIdx, j_compIdx));
80 const Scalar T_c_ij =
std::sqrt(T_c_i*T_c_j);
83 xxT_p += x_i*x_j*T_c_ij/p_c_ij;
84 xxT2_p += x_i*x_j*T_c_ij*T_c_ij/p_c_ij;
88 const LhsEval T_pc = xxT2_p/xxT_p;
89 const LhsEval p_pc = T_pc/xxT_p;
91 LhsEval p_pca = p_pc / MPa_atm;
96 for (
unsigned compIdx = 0; compIdx < FluidSystem::numComponents; ++compIdx) {
97 const Scalar& p_c = FluidSystem::criticalPressure(compIdx)/1e6;
98 const Scalar& T_c = FluidSystem::criticalTemperature(compIdx);
99 const Scalar Mm = FluidSystem::molarMass(compIdx) * 1000;
100 const auto& x = Opm::decay<LhsEval>(fluidState.moleFraction(phaseIdx, compIdx));
101 Scalar p_ca = p_c / MPa_atm;
107 mys = 34.0e-5*
Opm::pow(T_r,0.94)/zeta;
109 mys = 17.78e-5*
Opm::pow(4.58*T_r - 1.67, 0.625)/zeta;
116 std::vector<Scalar>
LBC = {0.10230,
122 LhsEval sumLBC = 0.0;
123 for (
int i = 0; i < 5; ++i) {
127 return (my0 + (
Opm::pow(sumLBC,4.0) - 1e-4)/zeta_tot -1.8366e-8*
Opm::pow(rho_r,13.992))/1e3;
static LhsEval LBCco2rich(const FluidState &fluidState, const Params &, unsigned phaseIdx)
Definition: LBCco2rich.hpp:46
static LhsEval LBC(const FluidState &fluidState, const Params &, unsigned phaseIdx)
Definition: LBC.hpp:47
Definition: Air_Mesitylene.hpp:34
Evaluation sqrt(const Evaluation &value)
Definition: MathToolbox.hpp:399
ReturnEval_< Evaluation1, Evaluation2 >::type pow(const Evaluation1 &base, const Evaluation2 &exp)
Definition: MathToolbox.hpp:416