ncpprimaryvariables.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 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*/
28#ifndef EWOMS_NCP_PRIMARY_VARIABLES_HH
29#define EWOMS_NCP_PRIMARY_VARIABLES_HH
30
31#include "ncpproperties.hh"
32
35
36#include <opm/material/constraintsolvers/NcpFlash.hpp>
37#include <opm/material/fluidstates/CompositionalFluidState.hpp>
38#include <opm/material/densead/Math.hpp>
39
40#include <dune/common/fvector.hh>
41
42namespace Opm {
43
53template <class TypeTag>
55{
56 using ParentType = FvBasePrimaryVariables<TypeTag>;
57
63
65 enum { pressure0Idx = Indices::pressure0Idx };
66 enum { saturation0Idx = Indices::saturation0Idx };
67 enum { fugacity0Idx = Indices::fugacity0Idx };
68
69 enum { numPhases = getPropValue<TypeTag, Properties::NumPhases>() };
70 enum { numComponents = getPropValue<TypeTag, Properties::NumComponents>() };
71 using ComponentVector = Dune::FieldVector<Scalar, numComponents>;
72
73 enum { enableEnergy = getPropValue<TypeTag, Properties::EnableEnergy>() };
75
76 using NcpFlash = Opm::NcpFlash<Scalar, FluidSystem>;
77 using Toolbox = Opm::MathToolbox<Evaluation>;
78
79public:
80 NcpPrimaryVariables() : ParentType()
81 {}
82
86 NcpPrimaryVariables(Scalar value) : ParentType(value)
87 {}
88
95
99 template <class FluidState>
100 void assignMassConservative(const FluidState& fluidState,
101 const MaterialLawParams& matParams,
102 bool isInEquilibrium = false)
103 {
104 using FsToolbox = Opm::MathToolbox<typename FluidState::Scalar>;
105
106#ifndef NDEBUG
107 // make sure the temperature is the same in all fluid phases
108 for (unsigned phaseIdx = 1; phaseIdx < numPhases; ++phaseIdx) {
109 assert(fluidState.temperature(0) == fluidState.temperature(phaseIdx));
110 }
111#endif // NDEBUG
112
113 // for the equilibrium case, we don't need complicated
114 // computations.
115 if (isInEquilibrium) {
116 assignNaive(fluidState);
117 return;
118 }
119
120 // use a flash calculation to calculate a fluid state in
121 // thermodynamic equilibrium
122 typename FluidSystem::template ParameterCache<Scalar> paramCache;
123 Opm::CompositionalFluidState<Scalar, FluidSystem> fsFlash;
124
125 // use the externally given fluid state as initial value for
126 // the flash calculation
127 fsFlash.assign(fluidState);
128
129 // calculate the phase densities
130 paramCache.updateAll(fsFlash);
131 for (unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
132 Scalar rho = FluidSystem::density(fsFlash, paramCache, phaseIdx);
133 fsFlash.setDensity(phaseIdx, rho);
134 }
135
136 // calculate the "global molarities"
137 ComponentVector globalMolarities(0.0);
138 for (unsigned compIdx = 0; compIdx < numComponents; ++compIdx) {
139 for (unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
140 globalMolarities[compIdx] +=
141 FsToolbox::value(fsFlash.saturation(phaseIdx))
142 * FsToolbox::value(fsFlash.molarity(phaseIdx, compIdx));
143 }
144 }
145
146 // run the flash calculation
147 NcpFlash::template solve<MaterialLaw>(fsFlash, matParams, paramCache, globalMolarities);
148
149 // use the result to assign the primary variables
150 assignNaive(fsFlash);
151 }
152
156 template <class FluidState>
157 void assignNaive(const FluidState& fluidState, unsigned refPhaseIdx = 0)
158 {
159 using FsToolbox = Opm::MathToolbox<typename FluidState::Scalar>;
160
161 // assign the phase temperatures. this is out-sourced to
162 // the energy module
163 EnergyModule::setPriVarTemperatures(*this, fluidState);
164
165 // assign fugacities.
166 typename FluidSystem::template ParameterCache<Scalar> paramCache;
167 paramCache.updatePhase(fluidState, refPhaseIdx);
168 Scalar pRef = FsToolbox::value(fluidState.pressure(refPhaseIdx));
169 for (unsigned compIdx = 0; compIdx < numComponents; ++compIdx) {
170 // we always compute the fugacities because they are quite exotic quantities
171 // and this easily forgotten to be specified
172 Scalar fugCoeff =
173 FluidSystem::template fugacityCoefficient<FluidState, Scalar>(fluidState,
174 paramCache,
175 refPhaseIdx,
176 compIdx);
177 (*this)[fugacity0Idx + compIdx] =
178 fugCoeff*fluidState.moleFraction(refPhaseIdx, compIdx)*pRef;
179 }
180
181 // assign pressure of first phase
182 (*this)[pressure0Idx] = FsToolbox::value(fluidState.pressure(/*phaseIdx=*/0));
183
184 // assign first M - 1 saturations
185 for (unsigned phaseIdx = 0; phaseIdx < numPhases - 1; ++phaseIdx)
186 (*this)[saturation0Idx + phaseIdx] = FsToolbox::value(fluidState.saturation(phaseIdx));
187 }
188};
189
190} // namespace Opm
191
192#endif
Provides the auxiliary methods required for consideration of the energy equation.
Definition: energymodule.hh:50
Represents the primary variables used by the a model.
Definition: fvbaseprimaryvariables.hh:52
Represents the primary variables used by the compositional multi-phase NCP model.
Definition: ncpprimaryvariables.hh:55
NcpPrimaryVariables()
Definition: ncpprimaryvariables.hh:80
void assignNaive(const FluidState &fluidState, unsigned refPhaseIdx=0)
Directly retrieve the primary variables from an arbitrary fluid state.
Definition: ncpprimaryvariables.hh:157
NcpPrimaryVariables & operator=(const NcpPrimaryVariables &value)=default
NcpPrimaryVariables(Scalar value)
Constructor with assignment from scalar.
Definition: ncpprimaryvariables.hh:86
void assignMassConservative(const FluidState &fluidState, const MaterialLawParams &matParams, bool isInEquilibrium=false)
Set the primary variables from an arbitrary fluid state in a mass conservative way.
Definition: ncpprimaryvariables.hh:100
NcpPrimaryVariables(const NcpPrimaryVariables &value)=default
Contains the classes required to consider energy as a conservation quantity in a multi-phase module.
Definition: blackoilboundaryratevector.hh:37
typename Properties::Detail::GetPropImpl< TypeTag, Property >::type::type GetPropType
get the type alias defined in the property (equivalent to old macro GET_PROP_TYPE(....
Definition: propertysystem.hh:235
Declares the properties required for the NCP compositional multi-phase model.