My Project
Loading...
Searching...
No Matches
Spe5ParameterCache.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 This file is part of the Open Porous Media project (OPM).
5
6 OPM is free software: you can redistribute it and/or modify
7 it under the terms of the GNU General Public License as published by
8 the Free Software Foundation, either version 2 of the License, or
9 (at your option) any later version.
10
11 OPM is distributed in the hope that it will be useful,
12 but WITHOUT ANY WARRANTY; without even the implied warranty of
13 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
14 GNU General Public License for more details.
15
16 You should have received a copy of the GNU General Public License
17 along with OPM. If not, see <http://www.gnu.org/licenses/>.
18
19 Consult the COPYING file in the top-level source directory of this
20 module for the precise wording of the license and the list of
21 copyright holders.
22*/
27#ifndef OPM_SPE5_PARAMETER_CACHE_HPP
28#define OPM_SPE5_PARAMETER_CACHE_HPP
29
32
35
36#include <cassert>
37
38namespace Opm {
39
44template <class Scalar, class FluidSystem>
46 : public ParameterCacheBase<Spe5ParameterCache<Scalar, FluidSystem> >
47{
50
51 typedef ::Opm::PengRobinson<Scalar> PengRobinson;
52
53 enum { numPhases = FluidSystem::numPhases };
54
55 enum { waterPhaseIdx = FluidSystem::waterPhaseIdx };
56 enum { oilPhaseIdx = FluidSystem::oilPhaseIdx };
57 enum { gasPhaseIdx = FluidSystem::gasPhaseIdx };
58
59 static_assert(static_cast<int>(oilPhaseIdx) >= 0, "Oil phase index must be non-negative");
60 static_assert(static_cast<int>(oilPhaseIdx) < static_cast<int>(numPhases),
61 "Oil phase index must be strictly less than FluidSystem's number of phases");
62
63 static_assert(static_cast<int>(gasPhaseIdx) >= 0, "Gas phase index must be non-negative");
64 static_assert(static_cast<int>(gasPhaseIdx) < static_cast<int>(numPhases),
65 "Gas phase index must be strictly less than FluidSystem's number of phases");
66
67 static_assert(static_cast<int>(waterPhaseIdx) >= 0, "Water phase index must be non-negative");
68 static_assert(static_cast<int>(waterPhaseIdx) < static_cast<int>(numPhases),
69 "Water phase index must be strictly less than FluidSystem's number of phases");
70
71public:
73 typedef PengRobinsonParamsMixture<Scalar, FluidSystem, oilPhaseIdx, /*useSpe5=*/true> OilPhaseParams;
75 typedef PengRobinsonParamsMixture<Scalar, FluidSystem, gasPhaseIdx, /*useSpe5=*/true> GasPhaseParams;
76
78 {
79 for (unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
80 VmUpToDate_[phaseIdx] = false;
81 Valgrind::SetUndefined(Vm_[phaseIdx]);
82 }
83 }
84
86 template <class FluidState>
87 void updatePhase(const FluidState& fluidState,
88 unsigned phaseIdx,
89 int exceptQuantities = ParentType::None)
90 {
91 assert ((phaseIdx == static_cast<unsigned int>(oilPhaseIdx)) ||
92 (phaseIdx == static_cast<unsigned int>(gasPhaseIdx)) ||
93 (phaseIdx == static_cast<unsigned int>(waterPhaseIdx)));
94
95 updateEosParams(fluidState, phaseIdx, exceptQuantities);
96
97 // if we don't need to recalculate the molar volume, we exit
98 // here
99 if (VmUpToDate_[phaseIdx])
100 return;
101
102 // update the phase's molar volume
103 updateMolarVolume_(fluidState, phaseIdx);
104 }
105
107 template <class FluidState>
108 void updateSingleMoleFraction(const FluidState& fluidState,
109 unsigned phaseIdx,
110 unsigned compIdx)
111 {
112 assert ((phaseIdx == static_cast<unsigned int>(oilPhaseIdx)) ||
113 (phaseIdx == static_cast<unsigned int>(gasPhaseIdx)) ||
114 (phaseIdx == static_cast<unsigned int>(waterPhaseIdx)));
115
116 if (phaseIdx == oilPhaseIdx)
117 oilPhaseParams_.updateSingleMoleFraction(fluidState, compIdx);
118 else if (phaseIdx == gasPhaseIdx)
119 gasPhaseParams_.updateSingleMoleFraction(fluidState, compIdx);
120
121 // update the phase's molar volume
122 updateMolarVolume_(fluidState, phaseIdx);
123 }
124
130 Scalar a(unsigned phaseIdx) const
131 {
132 switch (phaseIdx)
133 {
134 case oilPhaseIdx: return oilPhaseParams_.a();
135 case gasPhaseIdx: return gasPhaseParams_.a();
136 default:
137 throw std::logic_error("The a() parameter is only defined for "
138 "oil and gas phases");
139 };
140 }
141
147 Scalar b(unsigned phaseIdx) const
148 {
149 switch (phaseIdx)
150 {
151 case oilPhaseIdx: return oilPhaseParams_.b();
152 case gasPhaseIdx: return gasPhaseParams_.b();
153 default:
154 throw std::logic_error("The b() parameter is only defined for "
155 "oil and gas phases");
156 };
157 }
158
167 Scalar aPure(unsigned phaseIdx, unsigned compIdx) const
168 {
169 switch (phaseIdx)
170 {
171 case oilPhaseIdx: return oilPhaseParams_.pureParams(compIdx).a();
172 case gasPhaseIdx: return gasPhaseParams_.pureParams(compIdx).a();
173 default:
174 throw std::logic_error("The a() parameter is only defined for "
175 "oil and gas phases");
176 };
177 }
178
186 Scalar bPure(unsigned phaseIdx, unsigned compIdx) const
187 {
188 switch (phaseIdx)
189 {
190 case oilPhaseIdx: return oilPhaseParams_.pureParams(compIdx).b();
191 case gasPhaseIdx: return gasPhaseParams_.pureParams(compIdx).b();
192 default:
193 throw std::logic_error("The b() parameter is only defined for "
194 "oil and gas phases");
195 };
196 }
197
205 Scalar aCache(unsigned phaseIdx, unsigned compIdx, unsigned compJIdx) const
206 {
207 switch (phaseIdx)
208 {
209 case oilPhaseIdx: return oilPhaseParams_.getaCache(compIdx,compJIdx);
210 case gasPhaseIdx: return gasPhaseParams_.getaCache(compIdx,compJIdx);
211 default:
212 throw std::logic_error("The aCache() parameter is only defined for "
213 "oil and gas phase");
214 };
215 }
216
222 Scalar molarVolume(unsigned phaseIdx) const
223 { assert(VmUpToDate_[phaseIdx]); return Vm_[phaseIdx]; }
224
225
231 { return oilPhaseParams_; }
232
238 { return gasPhaseParams_; }
239
248 template <class FluidState>
249 void updateEosParams(const FluidState& fluidState,
250 unsigned phaseIdx,
251 int exceptQuantities = ParentType::None)
252 {
253 assert ((phaseIdx == static_cast<unsigned int>(oilPhaseIdx)) ||
254 (phaseIdx == static_cast<unsigned int>(gasPhaseIdx)) ||
255 (phaseIdx == static_cast<unsigned int>(waterPhaseIdx)));
256
257 if (!(exceptQuantities & ParentType::Temperature))
258 {
259 updatePure_(fluidState, phaseIdx);
260 updateMix_(fluidState, phaseIdx);
261 VmUpToDate_[phaseIdx] = false;
262 }
263 else if (!(exceptQuantities & ParentType::Composition))
264 {
265 updateMix_(fluidState, phaseIdx);
266 VmUpToDate_[phaseIdx] = false;
267 }
268 else if (!(exceptQuantities & ParentType::Pressure)) {
269 VmUpToDate_[phaseIdx] = false;
270 }
271 }
272
273protected:
280 template <class FluidState>
281 void updatePure_(const FluidState& fluidState, unsigned phaseIdx)
282 {
283 Scalar T = fluidState.temperature(phaseIdx);
284 Scalar p = fluidState.pressure(phaseIdx);
285
286 switch (phaseIdx)
287 {
288 case oilPhaseIdx: oilPhaseParams_.updatePure(T, p); break;
289 case gasPhaseIdx: gasPhaseParams_.updatePure(T, p); break;
290 //case waterPhaseIdx: waterPhaseParams_.updatePure(phaseIdx, temperature, pressure);break;
291 }
292 }
293
301 template <class FluidState>
302 void updateMix_(const FluidState& fluidState, unsigned phaseIdx)
303 {
304 Valgrind::CheckDefined(fluidState.averageMolarMass(phaseIdx));
305 switch (phaseIdx)
306 {
307 case oilPhaseIdx:
308 oilPhaseParams_.updateMix(fluidState);
309 break;
310 case gasPhaseIdx:
311 gasPhaseParams_.updateMix(fluidState);
312 break;
313 case waterPhaseIdx:
314 break;
315 }
316 }
317
318 template <class FluidState>
319 void updateMolarVolume_(const FluidState& fluidState,
320 unsigned phaseIdx)
321 {
322 VmUpToDate_[phaseIdx] = true;
323
324 // calculate molar volume of the phase (we will need this for the
325 // fugacity coefficients and the density anyway)
326 switch (phaseIdx) {
327 case gasPhaseIdx: {
328 // calculate molar volumes for the given composition. although
329 // this isn't a Peng-Robinson parameter strictly speaking, the
330 // molar volume appears in basically every quantity the fluid
331 // system can get queried, so it is okay to calculate it
332 // here...
333 Vm_[gasPhaseIdx] =
335 *this,
336 phaseIdx,
337 /*isGasPhase=*/true);
338 break;
339 }
340 case oilPhaseIdx: {
341 // calculate molar volumes for the given composition. although
342 // this isn't a Peng-Robinson parameter strictly speaking, the
343 // molar volume appears in basically every quantity the fluid
344 // system can get queried, so it is okay to calculate it
345 // here...
346 Vm_[oilPhaseIdx] =
348 *this,
349 phaseIdx,
350 /*isGasPhase=*/false);
351
352 break;
353 }
354 case waterPhaseIdx: {
355 // Density of water in the stock tank (i.e. atmospheric
356 // pressure) is specified as 62.4 lb/ft^3 by the SPE-5
357 // paper. Also 1 lb = 0.4535923 and 1 ft = 0.3048 m.
358 const Scalar stockTankWaterDensity = 62.4 * 0.45359237 / 0.028316847;
359 // Water compressibility is specified as 3.3e-6 per psi
360 // overpressure, where 1 psi = 6894.7573 Pa
361 Scalar overPressure = fluidState.pressure(waterPhaseIdx) - 1.013e5; // [Pa]
362 Scalar waterDensity =
363 stockTankWaterDensity * (1 + 3.3e-6*overPressure/6894.7573);
364
365 // convert water density [kg/m^3] to molar volume [m^3/mol]
366 Vm_[waterPhaseIdx] = fluidState.averageMolarMass(waterPhaseIdx)/waterDensity;
367 break;
368 }
369 }
370 }
371
372 bool VmUpToDate_[numPhases];
373 Scalar Vm_[numPhases];
374
375 OilPhaseParams oilPhaseParams_;
376 GasPhaseParams gasPhaseParams_;
377};
378
379} // namespace Opm
380
381#endif
Material properties of pure water .
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.
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
The mixing rule for the oil and the gas phases of the SPE5 problem.
Definition PengRobinsonParamsMixture.hpp:60
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
Specifies the parameter cache used by the SPE-5 fluid system.
Definition Spe5ParameterCache.hpp:47
void updatePhase(const FluidState &fluidState, unsigned phaseIdx, int exceptQuantities=ParentType::None)
Update all cached parameters of a specific fluid phase.
Definition Spe5ParameterCache.hpp:87
Scalar aPure(unsigned phaseIdx, unsigned compIdx) const
The Peng-Robinson attractive parameter for a pure component given the same temperature and pressure o...
Definition Spe5ParameterCache.hpp:167
void updatePure_(const FluidState &fluidState, unsigned phaseIdx)
Update all parameters of a phase which only depend on temperature and/or pressure.
Definition Spe5ParameterCache.hpp:281
const GasPhaseParams & gasPhaseParams() const
Returns the Peng-Robinson mixture parameters for the gas phase.
Definition Spe5ParameterCache.hpp:237
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 Spe5ParameterCache.hpp:186
const OilPhaseParams & oilPhaseParams() const
Returns the Peng-Robinson mixture parameters for the oil phase.
Definition Spe5ParameterCache.hpp:230
PengRobinsonParamsMixture< Scalar, FluidSystem, gasPhaseIdx, true > GasPhaseParams
The cached parameters for the gas phase.
Definition Spe5ParameterCache.hpp:75
void updateMix_(const FluidState &fluidState, unsigned phaseIdx)
Update all parameters of a phase which depend on the fluid composition.
Definition Spe5ParameterCache.hpp:302
Scalar aCache(unsigned phaseIdx, unsigned compIdx, unsigned compJIdx) const
TODO.
Definition Spe5ParameterCache.hpp:205
Scalar molarVolume(unsigned phaseIdx) const
Returns the molar volume of a phase [m^3/mol].
Definition Spe5ParameterCache.hpp:222
Scalar a(unsigned phaseIdx) const
The Peng-Robinson attractive parameter for a phase.
Definition Spe5ParameterCache.hpp:130
Scalar b(unsigned phaseIdx) const
The Peng-Robinson covolume for a phase.
Definition Spe5ParameterCache.hpp:147
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 Spe5ParameterCache.hpp:249
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 Spe5ParameterCache.hpp:108
PengRobinsonParamsMixture< Scalar, FluidSystem, oilPhaseIdx, true > OilPhaseParams
The cached parameters for the oil phase.
Definition Spe5ParameterCache.hpp:73
This class implements a small container which holds the transmissibility mulitpliers for all the face...
Definition Exceptions.hpp:30