My Project
Loading...
Searching...
No Matches
PTFlashParameterCache.hpp
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 Copyright 2022 NORCE.
5 Copyright 2022 SINTEF Digital, Mathematics and Cybernetics.
6
7 This file is part of the Open Porous Media project (OPM).
8
9 OPM is free software: you can redistribute it and/or modify
10 it under the terms of the GNU General Public License as published by
11 the Free Software Foundation, either version 2 of the License, or
12 (at your option) any later version.
13
14 OPM is distributed in the hope that it will be useful,
15 but WITHOUT ANY WARRANTY; without even the implied warranty of
16 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
17 GNU General Public License for more details.
18
19 You should have received a copy of the GNU General Public License
20 along with OPM. If not, see <http://www.gnu.org/licenses/>.
21
22 Consult the COPYING file in the top-level source directory of this
23 module for the precise wording of the license and the list of
24 copyright holders.
25*/
30#ifndef OPM_PTFlash_PARAMETER_CACHE_HPP
31#define OPM_PTFlash_PARAMETER_CACHE_HPP
32
36
37#include <cassert>
38
39namespace Opm {
40
45template <class Scalar, class FluidSystem>
47 : public Opm::ParameterCacheBase<PTFlashParameterCache<Scalar, FluidSystem> >
48{
51 using PengRobinson = Opm::PengRobinson<Scalar, false>; // false refer to discard some old code
52
53 enum { numPhases = FluidSystem::numPhases };
54 enum { oilPhaseIdx = FluidSystem::oilPhaseIdx };
55 enum { gasPhaseIdx = FluidSystem::gasPhaseIdx };
56
57 static_assert(static_cast<int>(oilPhaseIdx) >= 0, "Oil phase index must be non-negative");
58 static_assert(static_cast<int>(oilPhaseIdx) < static_cast<int>(numPhases),
59 "Oil phase index must be strictly less than FluidSystem's number of phases");
60
61 static_assert(static_cast<int>(gasPhaseIdx) >= 0, "Gas phase index must be non-negative");
62 static_assert(static_cast<int>(gasPhaseIdx) < static_cast<int>(numPhases),
63 "Gas phase index must be strictly less than FluidSystem's number of phases");
64
65public:
67 using OilPhaseParams = Opm::PengRobinsonParamsMixture<Scalar, FluidSystem, oilPhaseIdx, /*useChi=*/false>;
69 using GasPhaseParams = Opm::PengRobinsonParamsMixture<Scalar, FluidSystem, gasPhaseIdx, /*useChi=*/false>;
70
72 {
73 VmUpToDate_[oilPhaseIdx] = false;
74 Valgrind::SetUndefined(Vm_[oilPhaseIdx]);
75 VmUpToDate_[gasPhaseIdx] = false;
76 Valgrind::SetUndefined(Vm_[gasPhaseIdx]);
77 }
78
80 template <class FluidState>
81 void updatePhase(const FluidState& fluidState,
82 unsigned phaseIdx,
83 int exceptQuantities = ParentType::None)
84 {
85 assert ((phaseIdx == static_cast<unsigned int>(oilPhaseIdx)) ||
86 (phaseIdx == static_cast<unsigned int>(gasPhaseIdx)));
87
88 updateEosParams(fluidState, phaseIdx, exceptQuantities);
89
90 // update the phase's molar volume
91 updateMolarVolume_(fluidState, phaseIdx);
92 }
93
95 template <class FluidState>
96 void updateSingleMoleFraction(const FluidState& fluidState,
97 unsigned phaseIdx,
98 unsigned compIdx)
99 {
100 if (phaseIdx == oilPhaseIdx)
101 oilPhaseParams_.updateSingleMoleFraction(fluidState, compIdx);
102 if (phaseIdx == gasPhaseIdx)
103 gasPhaseParams_.updateSingleMoleFraction(fluidState, compIdx);
104 else
105 return;
106
107 // update the phase's molar volume
108 updateMolarVolume_(fluidState, phaseIdx);
109 }
110
116 Scalar a(unsigned phaseIdx) const
117 {
118 switch (phaseIdx)
119 {
120 case oilPhaseIdx: return oilPhaseParams_.a();
121 case gasPhaseIdx: return gasPhaseParams_.a();
122 default:
123 throw std::logic_error("The a() parameter is only defined for "
124 "oil and gas phases");
125 };
126 }
127
133 Scalar b(unsigned phaseIdx) const
134 {
135 switch (phaseIdx)
136 {
137 case oilPhaseIdx: return oilPhaseParams_.b();
138 case gasPhaseIdx: return gasPhaseParams_.b();
139 default:
140 throw std::logic_error("The b() parameter is only defined for "
141 "oil and gas phase");
142 };
143 }
144
153 Scalar aPure(unsigned phaseIdx, unsigned compIdx) const
154 {
155 switch (phaseIdx)
156 {
157 case oilPhaseIdx: return oilPhaseParams_.pureParams(compIdx).a();
158 case gasPhaseIdx: return gasPhaseParams_.pureParams(compIdx).a();
159 default:
160 throw std::logic_error("The a() parameter is only defined for "
161 "oil and gas phase");
162 };
163 }
164
172 Scalar bPure(unsigned phaseIdx, unsigned compIdx) const
173 {
174 switch (phaseIdx)
175 {
176 case oilPhaseIdx: return oilPhaseParams_.pureParams(compIdx).b();
177 case gasPhaseIdx: return gasPhaseParams_.pureParams(compIdx).b();
178 default:
179 throw std::logic_error("The b() parameter is only defined for "
180 "oil and gas phase");
181 };
182 }
183
184
185
193 Scalar aCache(unsigned phaseIdx, unsigned compIdx, unsigned compJIdx) const
194 {
195 switch (phaseIdx)
196 {
197 case oilPhaseIdx: return oilPhaseParams_.getaCache(compIdx,compJIdx);
198 case gasPhaseIdx: return gasPhaseParams_.getaCache(compIdx,compJIdx);
199 default:
200 throw std::logic_error("The aCache() parameter is only defined for "
201 "oil and gas phase");
202 };
203 }
204
210 Scalar molarVolume(unsigned phaseIdx) const
211 { assert(VmUpToDate_[phaseIdx]); return Vm_[phaseIdx]; }
212
213
219 { return oilPhaseParams_; }
220
226 // { throw std::invalid_argument("gas phase does not exist");}
227 { return gasPhaseParams_; }
228
237 template <class FluidState>
238 void updateEosParams(const FluidState& fluidState,
239 unsigned phaseIdx,
240 int exceptQuantities = ParentType::None)
241 {
242 assert ((phaseIdx == static_cast<unsigned int>(oilPhaseIdx)) ||
243 (phaseIdx == static_cast<unsigned int>(gasPhaseIdx)));
244
245 if (!(exceptQuantities & ParentType::Temperature))
246 {
247 updatePure_(fluidState, phaseIdx);
248 updateMix_(fluidState, phaseIdx);
249 VmUpToDate_[phaseIdx] = false;
250 }
251 else if (!(exceptQuantities & ParentType::Composition))
252 {
253 updateMix_(fluidState, phaseIdx);
254 VmUpToDate_[phaseIdx] = false;
255 }
256 else if (!(exceptQuantities & ParentType::Pressure)) {
257 VmUpToDate_[phaseIdx] = false;
258 }
259 }
260
261protected:
268 template <class FluidState>
269 void updatePure_(const FluidState& fluidState, unsigned phaseIdx)
270 {
271 Scalar T = decay<Scalar>(fluidState.temperature(phaseIdx));
272 Scalar p = decay<Scalar>(fluidState.pressure(phaseIdx));
273
274 switch (phaseIdx)
275 {
276 case oilPhaseIdx: oilPhaseParams_.updatePure(T, p); break;
277 case gasPhaseIdx: gasPhaseParams_.updatePure(T, p); break;
278 }
279 }
280
288 template <class FluidState>
289 void updateMix_(const FluidState& fluidState, unsigned phaseIdx)
290 {
291 Valgrind::CheckDefined(fluidState.averageMolarMass(phaseIdx));
292 switch (phaseIdx)
293 {
294 case oilPhaseIdx:
295 oilPhaseParams_.updateMix(fluidState);
296 break;
297 case gasPhaseIdx:
298 gasPhaseParams_.updateMix(fluidState);
299 break;
300 }
301 }
302
303 template <class FluidState>
304 void updateMolarVolume_(const FluidState& fluidState,
305 unsigned phaseIdx)
306 {
307 VmUpToDate_[phaseIdx] = true;
308
309 // calculate molar volume of the phase (we will need this for the
310 // fugacity coefficients and the density anyway)
311 switch (phaseIdx) {
312 case gasPhaseIdx: {
313 // calculate molar volumes for the given composition. although
314 // this isn't a Peng-Robinson parameter strictly speaking, the
315 // molar volume appears in basically every quantity the fluid
316 // system can get queried, so it is okay to calculate it
317 // here...
318 Vm_[gasPhaseIdx] = decay<Scalar> (
320 *this,
321 phaseIdx,
322 /*isGasPhase=*/true) );
323 break;
324 }
325 case oilPhaseIdx: {
326 // calculate molar volumes for the given composition. although
327 // this isn't a Peng-Robinson parameter strictly speaking, the
328 // molar volume appears in basically every quantity the fluid
329 // system can get queried, so it is okay to calculate it
330 // here...
331 Vm_[oilPhaseIdx] = decay<Scalar> (
333 *this,
334 phaseIdx,
335 /*isGasPhase=*/false) );
336
337 break;
338 }
339 };
340 }
341
342 bool VmUpToDate_[numPhases];
343 Scalar Vm_[numPhases];
344
345 OilPhaseParams oilPhaseParams_;
346 GasPhaseParams gasPhaseParams_;
347};
348
349} // namespace Opm
350
351#endif
The base class of the parameter caches of fluid systems.
The mixing rule for the oil and the gas phases of the SPE5 problem.
Implements the Peng-Robinson equation of state for liquids and gases.
Specifies the parameter cache used by the SPE-5 fluid system.
Definition PTFlashParameterCache.hpp:48
Opm::PengRobinsonParamsMixture< Scalar, FluidSystem, oilPhaseIdx, false > OilPhaseParams
The cached parameters for the oil phase.
Definition PTFlashParameterCache.hpp:67
const OilPhaseParams & oilPhaseParams() const
Returns the Peng-Robinson mixture parameters for the oil phase.
Definition PTFlashParameterCache.hpp:218
void updatePure_(const FluidState &fluidState, unsigned phaseIdx)
Update all parameters of a phase which only depend on temperature and/or pressure.
Definition PTFlashParameterCache.hpp:269
Scalar molarVolume(unsigned phaseIdx) const
Returns the molar volume of a phase [m^3/mol].
Definition PTFlashParameterCache.hpp:210
void updateSingleMoleFraction(const FluidState &fluidState, unsigned phaseIdx, unsigned compIdx)
Update all cached parameters of a specific fluid phase which depend on the mole fraction of a single ...
Definition PTFlashParameterCache.hpp:96
Scalar bPure(unsigned phaseIdx, unsigned compIdx) const
The Peng-Robinson covolume for a pure component given the same temperature and pressure of the phase.
Definition PTFlashParameterCache.hpp:172
Scalar aCache(unsigned phaseIdx, unsigned compIdx, unsigned compJIdx) const
TODO.
Definition PTFlashParameterCache.hpp:193
Opm::PengRobinsonParamsMixture< Scalar, FluidSystem, gasPhaseIdx, false > GasPhaseParams
The cached parameters for the gas phase.
Definition PTFlashParameterCache.hpp:69
void updateMix_(const FluidState &fluidState, unsigned phaseIdx)
Update all parameters of a phase which depend on the fluid composition.
Definition PTFlashParameterCache.hpp:289
void updatePhase(const FluidState &fluidState, unsigned phaseIdx, int exceptQuantities=ParentType::None)
Update all cached parameters of a specific fluid phase.
Definition PTFlashParameterCache.hpp:81
Scalar a(unsigned phaseIdx) const
The Peng-Robinson attractive parameter for a phase.
Definition PTFlashParameterCache.hpp:116
Scalar aPure(unsigned phaseIdx, unsigned compIdx) const
The Peng-Robinson attractive parameter for a pure component given the same temperature and pressure o...
Definition PTFlashParameterCache.hpp:153
Scalar b(unsigned phaseIdx) const
The Peng-Robinson covolume for a phase.
Definition PTFlashParameterCache.hpp:133
void updateEosParams(const FluidState &fluidState, unsigned phaseIdx, int exceptQuantities=ParentType::None)
Update all parameters required by the equation of state to calculate some quantities for the phase.
Definition PTFlashParameterCache.hpp:238
const GasPhaseParams & gasPhaseParams() const
Returns the Peng-Robinson mixture parameters for the gas phase.
Definition PTFlashParameterCache.hpp:225
The base class of the parameter caches of fluid systems.
Definition ParameterCacheBase.hpp:38
@ Temperature
The temperature has not been modified.
Definition ParameterCacheBase.hpp:48
@ None
All quantities have been (potentially) modified.
Definition ParameterCacheBase.hpp:45
@ Pressure
The pressures have not been modified.
Definition ParameterCacheBase.hpp:51
@ Composition
The compositions have not been modified.
Definition ParameterCacheBase.hpp:54
void updateSingleMoleFraction(const FluidState &fs, unsigned)
Calculates the "a" and "b" Peng-Robinson parameters for the mixture provided that only a single mole ...
Definition PengRobinsonParamsMixture.hpp:192
void updatePure(const FluidState &fluidState)
Update Peng-Robinson parameters for the pure components.
Definition PengRobinsonParamsMixture.hpp:82
Scalar getaCache(unsigned compIIdx, unsigned compJIdx) const
TODO.
Definition PengRobinsonParamsMixture.hpp:73
void updateMix(const FluidState &fs)
Calculates the "a" and "b" Peng-Robinson parameters for the mixture.
Definition PengRobinsonParamsMixture.hpp:143
const PureParams & pureParams(unsigned compIdx) const
Return the Peng-Robinson parameters of a pure substance,.
Definition PengRobinsonParamsMixture.hpp:201
Scalar a() const
Returns the attractive parameter 'a' of the Peng-Robinson fluid.
Definition PengRobinsonParams.hpp:50
Scalar b() const
Returns the repulsive parameter 'b' of the Peng-Robinson fluid.
Definition PengRobinsonParams.hpp:57
Implements the Peng-Robinson equation of state for liquids and gases.
Definition PengRobinson.hpp:60
static FluidState::Scalar computeMolarVolume(const FluidState &fs, Params &params, unsigned phaseIdx, bool isGasPhase)
Computes molar volumes where the Peng-Robinson EOS is true.
Definition PengRobinson.hpp:148
This class implements a small container which holds the transmissibility mulitpliers for all the face...
Definition Exceptions.hpp:30