28#ifndef OPM_FLASH_INTENSIVE_QUANTITIES_HH
29#define OPM_FLASH_INTENSIVE_QUANTITIES_HH
31#include <dune/common/fmatrix.hh>
32#include <dune/common/fvector.hh>
34#include <opm/material/Constants.hpp>
35#include <opm/material/common/Valgrind.hpp>
36#include <opm/material/fluidstates/CompositionalFluidState.hpp>
54template <
class TypeTag>
55class FlashIntensiveQuantities
56 :
public GetPropType<TypeTag, Properties::DiscIntensiveQuantities>
57 ,
public DiffusionIntensiveQuantities<TypeTag, getPropValue<TypeTag, Properties::EnableDiffusion>() >
58 ,
public EnergyIntensiveQuantities<TypeTag, getPropValue<TypeTag, Properties::EnableEnergy>() >
59 ,
public GetPropType<TypeTag, Properties::FluxModule>::FluxIntensiveQuantities
72 enum { z0Idx = Indices::z0Idx };
77 enum { dimWorld = GridView::dimensionworld };
78 enum { pressure0Idx = Indices::pressure0Idx };
85 using ComponentVector = Dune::FieldVector<Evaluation, numComponents>;
86 using DimMatrix = Dune::FieldMatrix<Scalar, dimWorld, dimWorld>;
88 using FluxIntensiveQuantities =
typename FluxModule::FluxIntensiveQuantities;
111 const auto& problem =
elemCtx.problem();
113 const Scalar
flashTolerance = Parameters::Get<Parameters::FlashTolerance<Scalar>>();
114 const int flashVerbosity = Parameters::Get<Parameters::FlashVerbosity>();
115 const std::string
flashTwoPhaseMethod = Parameters::Get<Parameters::FlashTwoPhaseMethod>();
118 ComponentVector
z(0.);
120 Evaluation
lastZ = 1.0;
125 z[numComponents - 1] =
lastZ;
127 Evaluation
sumz = 0.0;
135 Evaluation
p = priVars.makeEvaluation(pressure0Idx,
timeIdx);
146 const Evaluation&
Ltmp =
hint->fluidState().L();
147 fluidState_.setLvalue(
Ltmp);
156 const Evaluation&
Ltmp =
hint2->fluidState().L();
157 fluidState_.setLvalue(
Ltmp);
161 const Evaluation
Ktmp = fluidState_.wilsonK_(
compIdx);
164 const Evaluation&
Ltmp = -1.0;
165 fluidState_.setLvalue(
Ltmp);
173 std::cout <<
" updating the intensive quantities for Cell " <<
spatialIdx << std::endl;
180 std::cout <<
" \n After flash solve for cell " <<
spatialIdx << std::endl;
181 ComponentVector x,
y;
187 std::cout <<
" x for component: " <<
comp_idx <<
" is:" << std::endl;
188 std::cout << x[
comp_idx] << std::endl;
190 std::cout <<
" y for component: " <<
comp_idx <<
"is:" << std::endl;
193 const Evaluation&
L = fluidState_.L();
194 std::cout <<
" L is:" << std::endl;
195 std::cout <<
L << std::endl;
201 paramCache.updatePhase(fluidState_, FluidSystem::oilPhaseIdx);
203 const Scalar
R = Opm::Constants<Scalar>::R;
204 Evaluation
Z_L = (
paramCache.molarVolume(FluidSystem::oilPhaseIdx) * fluidState_.pressure(FluidSystem::oilPhaseIdx) )/
205 (
R * fluidState_.temperature(FluidSystem::oilPhaseIdx));
206 paramCache.updatePhase(fluidState_, FluidSystem::gasPhaseIdx);
207 Evaluation
Z_V = (
paramCache.molarVolume(FluidSystem::gasPhaseIdx) * fluidState_.pressure(FluidSystem::gasPhaseIdx) )/
208 (
R * fluidState_.temperature(FluidSystem::gasPhaseIdx));
213 Evaluation
L = fluidState_.L();
215 Evaluation Sg =
Opm::max(1 - So, 0.0);
220 fluidState_.setSaturation(0, So);
221 fluidState_.setSaturation(1, Sg);
223 fluidState_.setCompressFactor(0,
Z_L);
224 fluidState_.setCompressFactor(1,
Z_V);
228 std::cout <<
"So = " << So <<std::endl;
229 std::cout <<
"Sg = " << Sg <<std::endl;
234 std::cout <<
"So = " << So <<std::endl;
235 std::cout <<
"Sg = " << Sg <<std::endl;
236 std::cout <<
"Z_L = " <<
Z_L <<std::endl;
237 std::cout <<
"Z_V = " <<
Z_V <<std::endl;
246 MaterialLaw::relativePermeabilities(relativePermeability_,
248 Opm::Valgrind::CheckDefined(relativePermeability_);
259 Opm::Valgrind::CheckDefined(mobility_[
phaseIdx]);
262 fluidState_.setDensity(
phaseIdx, rho);
271 Opm::Valgrind::CheckDefined(porosity_);
290 {
return fluidState_; }
296 {
return intrinsicPerm_; }
302 {
return relativePermeability_[
phaseIdx]; }
316 {
return porosity_; }
319 DimMatrix intrinsicPerm_;
321 Evaluation porosity_;
322 std::array<Evaluation,numPhases> relativePermeability_;
323 std::array<Evaluation,numPhases> mobility_;
Provides the volumetric quantities required for the calculation of molecular diffusive fluxes.
Definition diffusionmodule.hh:141
Provides the volumetric quantities required for the energy equation.
Definition energymodule.hh:532
Contains the intensive quantities of the flash-based compositional multi-phase model.
Definition flashintensivequantities.hh:58
Opm::CompositionalFluidState< Evaluation, FluidSystem, enableEnergy > FluidState
The type of the object returned by the fluidState() method.
Definition flashintensivequantities.hh:91
const Evaluation & mobility(unsigned phaseIdx) const
Returns the effective mobility of a given phase within the control volume.
Definition flashintensivequantities.hh:307
const DimMatrix & intrinsicPermeability() const
Returns the intrinsic permeability tensor a degree of freedom.
Definition flashintensivequantities.hh:295
const FluidState & fluidState() const
Returns the phase state for the control-volume.
Definition flashintensivequantities.hh:289
const Evaluation & porosity() const
Returns the average porosity within the control volume.
Definition flashintensivequantities.hh:315
const Evaluation & relativePermeability(unsigned phaseIdx) const
Returns the relative permeability of a given phase within the control volume.
Definition flashintensivequantities.hh:301
void update(const ElementContext &elemCtx, unsigned dofIdx, unsigned timeIdx)
Definition flashintensivequantities.hh:105
Classes required for molecular diffusion.
Contains the classes required to consider energy as a conservation quantity in a multi-phase module.
Declares the properties required by the compositional multi-phase model based on flash calculations.
This file contains a set of helper functions used by VFPProd / VFPInj.
Definition blackoilboundaryratevector.hh:37
constexpr auto getPropValue()
get the value data member of a property
Definition propertysystem.hh:242
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
Defines the primary variable and equation indices for the compositional multi-phase model based on fl...
Declares the parameters for the compositional multi-phase model based on flash calculations.