My Project
Loading...
Searching...
No Matches
GenericOilGasWaterFluidSystem.hpp
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 Copyright 2023 SINTEF Digital, Mathematics and Cybernetics.
5
6 This file is part of the Open Porous Media project (OPM).
7
8 OPM is free software: you can redistribute it and/or modify
9 it under the terms of the GNU General Public License as published by
10 the Free Software Foundation, either version 2 of the License, or
11 (at your option) any later version.
12
13 OPM is distributed in the hope that it will be useful,
14 but WITHOUT ANY WARRANTY; without even the implied warranty of
15 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
16 GNU General Public License for more details.
17
18 You should have received a copy of the GNU General Public License
19 along with OPM. If not, see <http://www.gnu.org/licenses/>.
20
21 Consult the COPYING file in the top-level source directory of this
22 module for the precise wording of the license and the list of
23 copyright holders.
24*/
25
26#ifndef OPM_GENERIC_OIL_GAS_WATER_FLUIDSYSTEM_HPP
27#define OPM_GENERIC_OIL_GAS_WATER_FLUIDSYSTEM_HPP
28
29#include <opm/common/OpmLog/OpmLog.hpp>
30
31#if HAVE_ECL_INPUT
32#include <opm/input/eclipse/EclipseState/EclipseState.hpp>
33#include <opm/input/eclipse/Schedule/Schedule.hpp>
34#endif
35
36#include <opm/material/eos/CubicEOS.hpp>
39#include <opm/material/fluidsystems/PTFlashParameterCache.hpp> // TODO: this is something else need to check
41
42#include <cassert>
43#include <cstddef>
44#include <string>
45#include <string_view>
46
47#include <fmt/format.h>
48
49namespace Opm {
50
59 template<class Scalar, int NumComp, bool enableWater>
61 : public BaseFluidSystem<Scalar, GenericOilGasWaterFluidSystem<Scalar, NumComp, enableWater> > {
62 public:
63 // TODO: I do not think these should be constant in fluidsystem, will try to make it non-constant later
64 static constexpr bool waterEnabled = enableWater;
65 static constexpr int numPhases = enableWater ? 3 : 2;
66 static constexpr int numComponents = NumComp;
67 static constexpr int numMisciblePhases = 2;
68 // \Note: not totally sure when we should distinguish numMiscibleComponents and numComponents.
69 // Possibly when with a dummy phase like water?
70 static constexpr int numMiscibleComponents = NumComp;
71 // TODO: phase location should be more general
72 static constexpr int oilPhaseIdx = 0;
73 static constexpr int gasPhaseIdx = 1;
74 static constexpr int waterPhaseIdx = 2;
75
76 static constexpr int waterCompIdx = -1;
77 static constexpr int oilCompIdx = 0;
78 static constexpr int gasCompIdx = 1;
79 static constexpr int compositionSwitchIdx = -1; // equil initializer
80
81 template <class ValueType>
86
88 std::string name;
89 Scalar molar_mass; // unit: g/mol
90 Scalar critic_temp; // unit: K
91 Scalar critic_pres; // unit: parscal
92 Scalar critic_vol; // unit: m^3/kmol
93 Scalar acentric_factor; // unit: dimension less
94
95 ComponentParam(const std::string_view name_, const Scalar molar_mass_, const Scalar critic_temp_,
96 const Scalar critic_pres_, const Scalar critic_vol_, const Scalar acentric_factor_)
97 : name(name_),
98 molar_mass(molar_mass_),
99 critic_temp(critic_temp_),
100 critic_pres(critic_pres_),
101 critic_vol(critic_vol_),
102 acentric_factor(acentric_factor_)
103 {}
104 };
105
106 static bool phaseIsActive(unsigned phaseIdx)
107 {
108 if constexpr (enableWater) {
109 assert(phaseIdx < numPhases);
110 return true;
111 }
112 else {
113 if (phaseIdx == waterPhaseIdx)
114 return false;
115 assert(phaseIdx < numPhases);
116 return true;
117 }
118 }
119
120 template<typename Param>
121 static void addComponent(const Param& param)
122 {
123 // Check if the current size is less than the maximum allowed components.
124 if (component_param_.size() < numComponents) {
125 component_param_.push_back(param);
126 } else {
127 // Adding another component would exceed the limit.
128 const std::string msg = fmt::format("The fluid system has reached its maximum capacity of {} components,"
129 "the component '{}' will not be added.", NumComp, param.name);
130 OpmLog::note(msg);
131 // Optionally, throw an exception?
132 }
133 }
134
135#if HAVE_ECL_INPUT
139 static void initFromState(const EclipseState& eclState, const Schedule& schedule)
140 {
141 // TODO: we are not considering the EOS region for now
142 const auto& comp_config = eclState.compositionalConfig();
143 // how should we utilize the numComps from the CompositionalConfig?
144 using FluidSystem = GenericOilGasWaterFluidSystem<Scalar, NumComp, enableWater>;
145 const std::size_t num_comps = comp_config.numComps();
146 // const std::size_t num_eos_region = comp_config.
147 assert(num_comps == NumComp);
148 const auto& names = comp_config.compName();
149 // const auto& eos_type = comp_config.eosType(0);
150 const auto& molar_weight = comp_config.molecularWeights(0);
151 const auto& acentric_factor = comp_config.acentricFactors(0);
152 const auto& critic_pressure = comp_config.criticalPressure(0);
153 const auto& critic_temp = comp_config.criticalTemperature(0);
154 const auto& critic_volume = comp_config.criticalVolume(0);
155 FluidSystem::init();
156 using CompParm = typename FluidSystem::ComponentParam;
157 for (std::size_t c = 0; c < num_comps; ++c) {
158 // we use m^3/kmol for the critic volume in the flash calculation, so we multiply 1.e3 for the critic volume
159 FluidSystem::addComponent(CompParm{names[c], molar_weight[c], critic_temp[c], critic_pressure[c],
160 critic_volume[c] * 1.e3, acentric_factor[c]});
161 }
162 FluidSystem::printComponentParams();
163 interaction_coefficients_ = comp_config.binaryInteractionCoefficient(0);
164
165 // Init. water pvt from deck
166 waterPvt_->initFromState(eclState, schedule);
167
168 }
169#endif // HAVE_ECL_INPUT
170
171 static void init()
172 {
173 waterPvt_ = std::make_shared<WaterPvt>();
174 component_param_.reserve(numComponents);
175 }
176
180 static void setWaterPvt(std::shared_ptr<WaterPvt> pvtObj)
181 { waterPvt_ = std::move(pvtObj); }
182
188 static Scalar acentricFactor(unsigned compIdx)
189 {
190 assert(isConsistent());
191 assert(compIdx < numComponents);
192
193 return component_param_[compIdx].acentric_factor;
194 }
200 static Scalar criticalTemperature(unsigned compIdx)
201 {
202 assert(isConsistent());
203 assert(compIdx < numComponents);
204
205 return component_param_[compIdx].critic_temp;
206 }
212 static Scalar criticalPressure(unsigned compIdx)
213 {
214 assert(isConsistent());
215 assert(compIdx < numComponents);
216
217 return component_param_[compIdx].critic_pres;
218 }
224 static Scalar criticalVolume(unsigned compIdx)
225 {
226 assert(isConsistent());
227 assert(compIdx < numComponents);
228
229 return component_param_[compIdx].critic_vol;
230 }
231
233 static Scalar molarMass(unsigned compIdx)
234 {
235 assert(isConsistent());
236 assert(compIdx < numComponents);
237
238 return component_param_[compIdx].molar_mass;
239 }
240
245 static Scalar interactionCoefficient(unsigned comp1Idx, unsigned comp2Idx)
246 {
247 assert(isConsistent());
248 assert(comp1Idx < numComponents);
249 assert(comp2Idx < numComponents);
250 if (interaction_coefficients_.empty() || comp2Idx == comp1Idx) {
251 return 0.0;
252 }
253 // make sure row is the bigger value compared to column number
254 const auto [column, row] = std::minmax(comp1Idx, comp2Idx);
255 const unsigned index = (row * (row - 1) / 2 + column); // it is the current understanding
256 return interaction_coefficients_[index];
257 }
258
260 static std::string_view phaseName(unsigned phaseIdx)
261 {
262 static const std::string_view name[] = {"o", // oleic phase
263 "g", // gas phase
264 "w"}; // aqueous phase
265
266 assert(phaseIdx < numPhases);
267 return name[phaseIdx];
268 }
269
271 static std::string_view componentName(unsigned compIdx)
272 {
273 assert(isConsistent());
274 assert(compIdx < numComponents);
275
276 return component_param_[compIdx].name;
277 }
278
282 template <class FluidState, class LhsEval = typename FluidState::Scalar, class ParamCacheEval = LhsEval>
283 static LhsEval density(const FluidState& fluidState,
284 const ParameterCache<ParamCacheEval>& paramCache,
285 unsigned phaseIdx)
286 {
287 assert(isConsistent());
288 assert(phaseIdx < numPhases);
289
290 if (phaseIdx == oilPhaseIdx || phaseIdx == gasPhaseIdx) {
291 return decay<LhsEval>(fluidState.averageMolarMass(phaseIdx) / paramCache.molarVolume(phaseIdx));
292 }
293 else {
294 const LhsEval& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
295 const LhsEval& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
296 const Scalar& rho_w_ref = waterPvt_->waterReferenceDensity(0);
297 const LhsEval& bw = waterPvt_->inverseFormationVolumeFactor(0, T, p, LhsEval(0.0), LhsEval(0.0));
298 return rho_w_ref * bw;
299 }
300
301 return {};
302 }
303
305 template <class FluidState, class LhsEval = typename FluidState::Scalar, class ParamCacheEval = LhsEval>
306 static LhsEval viscosity(const FluidState& fluidState,
307 const ParameterCache<ParamCacheEval>& paramCache,
308 unsigned phaseIdx)
309 {
310 assert(isConsistent());
311 assert(phaseIdx < numPhases);
312
313 if (phaseIdx == oilPhaseIdx || phaseIdx == gasPhaseIdx) {
314 // Use LBC method to calculate viscosity
315 return decay<LhsEval>(ViscosityModel::LBC(fluidState, paramCache, phaseIdx));
316 }
317 else {
318 const LhsEval& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
319 const LhsEval& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
320 return waterPvt_->viscosity(0, T, p, LhsEval(0.0), LhsEval(0.0));
321 }
322 }
323
325 template <class FluidState, class LhsEval = typename FluidState::Scalar, class ParamCacheEval = LhsEval>
326 static LhsEval fugacityCoefficient(const FluidState& fluidState,
327 const ParameterCache<ParamCacheEval>& paramCache,
328 unsigned phaseIdx,
329 unsigned compIdx)
330 {
331 if (waterEnabled && phaseIdx == static_cast<unsigned int>(waterPhaseIdx))
332 return LhsEval(0.0);
333
334 assert(isConsistent());
335 assert(phaseIdx < numPhases);
336 assert(compIdx < numComponents);
337
338 return decay<LhsEval>(CubicEOS::computeFugacityCoefficient(fluidState, paramCache, phaseIdx, compIdx));
339 }
340
341 // TODO: the following interfaces are needed by function checkFluidSystem()
343 static bool isCompressible([[maybe_unused]] unsigned phaseIdx)
344 {
345 assert(phaseIdx < numPhases);
346
347 return true;
348 }
349
351 static bool isIdealMixture([[maybe_unused]] unsigned phaseIdx)
352 {
353 assert(phaseIdx < numPhases);
354
355 return false;
356 }
357
359 static bool isLiquid(unsigned phaseIdx)
360 {
361 assert(phaseIdx < numPhases);
362
363 return (phaseIdx == oilPhaseIdx);
364 }
365
367 static bool isIdealGas(unsigned phaseIdx)
368 {
369 assert(phaseIdx < numPhases);
370
371 return (phaseIdx == gasPhaseIdx);
372 }
373 // the following funcitons are needed to compile the GenericOutputBlackoilModule
374 // not implemented for this FluidSystem yet
375 template <class LhsEval>
376 static LhsEval convertXwGToxwG(const LhsEval&, unsigned)
377 {
378 assert(false && "convertXwGToxwG not implemented for GenericOilGasWaterFluidSystem!");
379 return 0.;
380 }
381 template <class LhsEval>
382 static LhsEval convertXoGToxoG(const LhsEval&, unsigned)
383 {
384 assert(false && "convertXoGToxoG not implemented for GenericOilGasWaterFluidSystem!");
385 return 0.;
386 }
387
388 template <class LhsEval>
389 static LhsEval convertxoGToXoG(const LhsEval&, unsigned)
390 {
391 assert(false && "convertxoGToXoG not implemented for GenericOilGasWaterFluidSystem!");
392 return 0.;
393 }
394
395 template <class LhsEval>
396 static LhsEval convertXgOToxgO(const LhsEval&, unsigned)
397 {
398 assert(false && "convertXgOToxgO not implemented for GenericOilGasWaterFluidSystem!");
399 return 0.;
400 }
401
402 template <class LhsEval>
403 static LhsEval convertRswToXwG(const LhsEval&, unsigned)
404 {
405 assert(false && "convertRswToXwG not implemented for GenericOilGasWaterFluidSystem!");
406 return 0.;
407 }
408
409 template <class LhsEval>
410 static LhsEval convertRvwToXgW(const LhsEval&, unsigned)
411 {
412 assert(false && "convertRvwToXgW not implemented for GenericOilGasWaterFluidSystem!");
413 return 0.;
414 }
415
416 template <class LhsEval>
417 static LhsEval convertXgWToxgW(const LhsEval&, unsigned)
418 {
419 assert(false && "convertXgWToxgW not implemented for GenericOilGasWaterFluidSystem!");
420 return 0.;
421 }
422
423 static bool enableDissolvedGas()
424 {
425 return false;
426 }
427
428 static bool enableDissolvedGasInWater()
429 {
430 return false;
431 }
432
433 static bool enableVaporizedWater()
434 {
435 return false;
436 }
437
438 static bool enableVaporizedOil()
439 {
440 return false;
441 }
442
443 private:
444 static bool isConsistent() {
445 return component_param_.size() == NumComp;
446 }
447
448 static std::vector<ComponentParam> component_param_;
449 static std::vector<Scalar> interaction_coefficients_;
450 static std::shared_ptr<WaterPvt> waterPvt_;
451
452 public:
453 static std::string printComponentParams() {
454 std::string result = "Components Information:\n";
455 for (const auto& param : component_param_) {
456 result += fmt::format("Name: {}\n", param.name);
457 result += fmt::format("Molar Mass: {} g/mol\n", param.molar_mass);
458 result += fmt::format("Critical Temperature: {} K\n", param.critic_temp);
459 result += fmt::format("Critical Pressure: {} Pascal\n", param.critic_pres);
460 result += fmt::format("Critical Volume: {} m^3/kmol\n", param.critic_vol);
461 result += fmt::format("Acentric Factor: {}\n", param.acentric_factor);
462 result += "---------------------------------\n";
463 }
464 return result;
465 }
466 };
467
468 template <class Scalar, int NumComp, bool enableWater>
469 std::vector<typename GenericOilGasWaterFluidSystem<Scalar, NumComp, enableWater>::ComponentParam>
470 GenericOilGasWaterFluidSystem<Scalar, NumComp, enableWater>::component_param_;
471
472 template <class Scalar, int NumComp, bool enableWater>
473 std::vector<Scalar>
474 GenericOilGasWaterFluidSystem<Scalar, NumComp, enableWater>::interaction_coefficients_;
475
476 template <class Scalar, int NumComp, bool enableWater>
477 std::shared_ptr<WaterPvtMultiplexer<Scalar> >
478 GenericOilGasWaterFluidSystem<Scalar, NumComp, enableWater>::waterPvt_;
479
480}
481#endif // OPM_GENERIC_OIL_GAS_WATER_FLUIDSYSTEM_HPP
The base class for all fluid systems.
Specifies the parameter cache used by the SPE-5 fluid system.
This class represents the Pressure-Volume-Temperature relations of the water phase in the black-oil m...
The base class for all fluid systems.
Definition BaseFluidSystem.hpp:43
Scalar Scalar
The type used for scalar quantities.
Definition BaseFluidSystem.hpp:48
Definition CubicEOS.hpp:34
A two phase system that can contain NumComp components.
Definition GenericOilGasWaterFluidSystem.hpp:61
static bool isCompressible(unsigned phaseIdx)
Returns true if and only if a fluid phase is assumed to be compressible.
Definition GenericOilGasWaterFluidSystem.hpp:343
static LhsEval fugacityCoefficient(const FluidState &fluidState, const ParameterCache< ParamCacheEval > &paramCache, unsigned phaseIdx, unsigned compIdx)
Calculate the fugacity coefficient [Pa] of an individual component in a fluid phase.
Definition GenericOilGasWaterFluidSystem.hpp:326
static Scalar acentricFactor(unsigned compIdx)
The acentric factor of a component [].
Definition GenericOilGasWaterFluidSystem.hpp:188
static std::string_view phaseName(unsigned phaseIdx)
Return the human readable name of a fluid phase.
Definition GenericOilGasWaterFluidSystem.hpp:260
static bool isIdealGas(unsigned phaseIdx)
Returns true if and only if a fluid phase is assumed to be an ideal gas.
Definition GenericOilGasWaterFluidSystem.hpp:367
static Scalar criticalTemperature(unsigned compIdx)
Critical temperature of a component [K].
Definition GenericOilGasWaterFluidSystem.hpp:200
static LhsEval density(const FluidState &fluidState, const ParameterCache< ParamCacheEval > &paramCache, unsigned phaseIdx)
Calculate the density [kg/m^3] of a fluid phase.
Definition GenericOilGasWaterFluidSystem.hpp:283
static Scalar interactionCoefficient(unsigned comp1Idx, unsigned comp2Idx)
Returns the interaction coefficient for two components.
Definition GenericOilGasWaterFluidSystem.hpp:245
static Scalar criticalPressure(unsigned compIdx)
Critical pressure of a component [Pa].
Definition GenericOilGasWaterFluidSystem.hpp:212
static bool isIdealMixture(unsigned phaseIdx)
Returns true if and only if a fluid phase is assumed to be an ideal mixture.
Definition GenericOilGasWaterFluidSystem.hpp:351
static Scalar molarMass(unsigned compIdx)
Return the molar mass of a component in [kg/mol].
Definition GenericOilGasWaterFluidSystem.hpp:233
static std::string_view componentName(unsigned compIdx)
Return the human readable name of a component.
Definition GenericOilGasWaterFluidSystem.hpp:271
static bool isLiquid(unsigned phaseIdx)
Return whether a phase is liquid.
Definition GenericOilGasWaterFluidSystem.hpp:359
static Scalar criticalVolume(unsigned compIdx)
Critical volume of a component [m3].
Definition GenericOilGasWaterFluidSystem.hpp:224
static void setWaterPvt(std::shared_ptr< WaterPvt > pvtObj)
Set the pressure-volume-saturation (PVT) relations for the water phase.
Definition GenericOilGasWaterFluidSystem.hpp:180
static LhsEval viscosity(const FluidState &fluidState, const ParameterCache< ParamCacheEval > &paramCache, unsigned phaseIdx)
Calculate the dynamic viscosity of a fluid phase [Pa*s].
Definition GenericOilGasWaterFluidSystem.hpp:306
Specifies the parameter cache used by the SPE-5 fluid system.
Definition PTFlashParameterCache.hpp:51
Scalar molarVolume(unsigned phaseIdx) const
Returns the molar volume of a phase [m^3/mol].
Definition PTFlashParameterCache.hpp:283
Definition LBC.hpp:40
This class represents the Pressure-Volume-Temperature relations of the water phase in the black-oil m...
Definition WaterPvtMultiplexer.hpp:90
This class implements a small container which holds the transmissibility mulitpliers for all the face...
Definition Exceptions.hpp:30
Definition GenericOilGasWaterFluidSystem.hpp:87