version 3.10.0
Loading...
Searching...
No Matches
porousmediumflow/3pwateroil/localresidual.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// SPDX-FileCopyrightText: Copyright © DuMux Project contributors, see AUTHORS.md in root folder
5// SPDX-License-Identifier: GPL-3.0-or-later
6//
13
14#ifndef DUMUX_3P2CNI_LOCAL_RESIDUAL_HH
15#define DUMUX_3P2CNI_LOCAL_RESIDUAL_HH
16
21
22namespace Dumux {
28template<class TypeTag>
31{
32protected:
37 using FVElementGeometry = typename GridGeometry::LocalView;
38 using SubControlVolume = typename FVElementGeometry::SubControlVolume;
39 using SubControlVolumeFace = typename FVElementGeometry::SubControlVolumeFace;
44 using GridView = typename GridGeometry::GridView;
45 using Element = typename GridView::template Codim<0>::Entity;
50
51 enum {
54
55 conti0EqIdx = Indices::conti0EqIdx,
57
58 wPhaseIdx = FluidSystem::wPhaseIdx,
59 nPhaseIdx = FluidSystem::nPhaseIdx,
60 gPhaseIdx = FluidSystem::gPhaseIdx,
61
62 wCompIdx = FluidSystem::wCompIdx,
63 nCompIdx = FluidSystem::nCompIdx,
64 };
65
68public:
69 using ParentType::ParentType;
70
83 const SubControlVolume& scv,
84 const VolumeVariables& volVars) const
85 {
86 NumEqVector storage(0.0);
87
88 const auto massOrMoleDensity = [](const auto& volVars, const int phaseIdx)
89 { return useMoles ? volVars.molarDensity(phaseIdx) : volVars.density(phaseIdx); };
90
91 const auto massOrMoleFraction= [](const auto& volVars, const int phaseIdx, const int compIdx)
92 { return useMoles ? volVars.moleFraction(phaseIdx, compIdx) : volVars.massFraction(phaseIdx, compIdx); };
93
94 // compute storage term of all components within all phases
95 for (int phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx)
96 {
97 for (int compIdx = 0; compIdx < numComponents; ++compIdx)
98 {
99 int eqIdx = (compIdx == wCompIdx) ? conti0EqIdx : conti1EqIdx;
100 storage[eqIdx] += volVars.porosity()
101 * volVars.saturation(phaseIdx)
102 * massOrMoleDensity(volVars, phaseIdx)
103 * massOrMoleFraction(volVars, phaseIdx, compIdx);
104 }
105
107 EnergyLocalResidual::fluidPhaseStorage(storage, problem, scv, volVars, phaseIdx);
108 }
109
111 EnergyLocalResidual::solidPhaseStorage(storage, scv, volVars);
112
113 return storage;
114 }
115
128 const Element& element,
129 const FVElementGeometry& fvGeometry,
130 const ElementVolumeVariables& elemVolVars,
131 const SubControlVolumeFace& scvf,
132 const ElementFluxVariablesCache& elemFluxVarsCache) const
133 {
134 FluxVariables fluxVars;
135 fluxVars.init(problem, element, fvGeometry, elemVolVars, scvf, elemFluxVarsCache);
136
137 // get upwind weights into local scope
138 NumEqVector flux(0.0);
139
140 const auto massOrMoleDensity = [](const auto& volVars, const int phaseIdx)
141 { return useMoles ? volVars.molarDensity(phaseIdx) : volVars.density(phaseIdx); };
142
143 const auto massOrMoleFraction= [](const auto& volVars, const int phaseIdx, const int compIdx)
144 { return useMoles ? volVars.moleFraction(phaseIdx, compIdx) : volVars.massFraction(phaseIdx, compIdx); };
145
146 // advective fluxes
147 for (int phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx)
148 {
149 for (int compIdx = 0; compIdx < numComponents; ++compIdx)
150 {
151 const auto upwindTerm = [&massOrMoleDensity, &massOrMoleFraction, phaseIdx, compIdx] (const auto& volVars)
152 { return massOrMoleDensity(volVars, phaseIdx)*massOrMoleFraction(volVars, phaseIdx, compIdx)*volVars.mobility(phaseIdx); };
153
154 // get equation index
155 auto eqIdx = conti0EqIdx + compIdx;
156 flux[eqIdx] += fluxVars.advectiveFlux(phaseIdx, upwindTerm);
157 }
158
160 EnergyLocalResidual::heatConvectionFlux(flux, fluxVars, phaseIdx);
161 }
162
164 EnergyLocalResidual::heatConductionFlux(flux, fluxVars);
165
166 // diffusive fluxes
167 static constexpr auto referenceSystemFormulation = FluxVariables::MolecularDiffusionType::referenceSystemFormulation();
168 const auto diffusionFluxesWPhase = fluxVars.molecularDiffusionFlux(wPhaseIdx);
169 Scalar jNW = diffusionFluxesWPhase[nCompIdx];
170 Scalar jWW = -jNW;
171 // check for the reference system and adapt units of the diffusive flux accordingly.
172 if (referenceSystemFormulation == ReferenceSystemFormulation::massAveraged)
173 {
174 jNW /= FluidSystem::molarMass(nCompIdx);
175 jWW /= FluidSystem::molarMass(wCompIdx);
176 }
177
178 const auto diffusionFluxesGPhase = fluxVars.molecularDiffusionFlux(gPhaseIdx);
179 Scalar jWG = diffusionFluxesGPhase[wCompIdx];
180 Scalar jNG = -jWG;
181 // check for the reference system and adapt units of the diffusive flux accordingly.
182 if (referenceSystemFormulation == ReferenceSystemFormulation::massAveraged)
183 {
184 jWG /= FluidSystem::molarMass(wCompIdx);
185 jNG /= FluidSystem::molarMass(nCompIdx);
186 }
187
188 const auto diffusionFluxesNPhase = fluxVars.molecularDiffusionFlux(nPhaseIdx);
189 Scalar jWN = diffusionFluxesNPhase[wCompIdx];
190 Scalar jNN = -jWN;
191 // check for the reference system and adapt units of the diffusive flux accordingly.
192 if (referenceSystemFormulation == ReferenceSystemFormulation::massAveraged)
193 {
194 jWN /= FluidSystem::molarMass(wCompIdx);
195 jNN /= FluidSystem::molarMass(nCompIdx);
196 }
197
198 flux[conti0EqIdx] += jWW+jWG+jWN;
199 flux[conti1EqIdx] += jNW+jNG+jNN;
200
201 return flux;
202 }
203};
204
205} // end namespace Dumux
206
207#endif
Element-wise calculation of the local residual for problems using the ThreePWaterOil fully implicit m...
Definition porousmediumflow/3pwateroil/localresidual.hh:31
@ numComponents
Definition porousmediumflow/3pwateroil/localresidual.hh:53
@ conti0EqIdx
Index of the mass conservation equation for the water component.
Definition porousmediumflow/3pwateroil/localresidual.hh:55
@ nPhaseIdx
Definition porousmediumflow/3pwateroil/localresidual.hh:59
@ gPhaseIdx
Definition porousmediumflow/3pwateroil/localresidual.hh:60
@ conti1EqIdx
Index of the mass conservation equation for the contaminant component.
Definition porousmediumflow/3pwateroil/localresidual.hh:56
@ wPhaseIdx
Definition porousmediumflow/3pwateroil/localresidual.hh:58
@ nCompIdx
Definition porousmediumflow/3pwateroil/localresidual.hh:63
@ numPhases
Definition porousmediumflow/3pwateroil/localresidual.hh:52
@ wCompIdx
Definition porousmediumflow/3pwateroil/localresidual.hh:62
NumEqVector computeFlux(const Problem &problem, const Element &element, const FVElementGeometry &fvGeometry, const ElementVolumeVariables &elemVolVars, const SubControlVolumeFace &scvf, const ElementFluxVariablesCache &elemFluxVarsCache) const
Evaluates the total flux of all conservation quantities over a face of a sub-control volume.
Definition porousmediumflow/3pwateroil/localresidual.hh:127
typename GridGeometry::GridView GridView
Definition porousmediumflow/3pwateroil/localresidual.hh:44
typename GetPropType< TypeTag, Properties::GridFluxVariablesCache >::LocalView ElementFluxVariablesCache
Definition porousmediumflow/3pwateroil/localresidual.hh:42
GetPropType< TypeTag, Properties::VolumeVariables > VolumeVariables
Definition porousmediumflow/3pwateroil/localresidual.hh:47
NumEqVector computeStorage(const Problem &problem, const SubControlVolume &scv, const VolumeVariables &volVars) const
Evaluates the amount of all conservation quantities (e.g. phase mass) within a sub-control volume.
Definition porousmediumflow/3pwateroil/localresidual.hh:82
GetPropType< TypeTag, Properties::GridGeometry > GridGeometry
Definition porousmediumflow/3pwateroil/localresidual.hh:33
typename GridView::template Codim< 0 >::Entity Element
Definition porousmediumflow/3pwateroil/localresidual.hh:45
GetPropType< TypeTag, Properties::Problem > Problem
Definition porousmediumflow/3pwateroil/localresidual.hh:35
typename GetPropType< TypeTag, Properties::ModelTraits >::Indices Indices
Definition porousmediumflow/3pwateroil/localresidual.hh:43
typename GridGeometry::LocalView FVElementGeometry
Definition porousmediumflow/3pwateroil/localresidual.hh:37
typename GetPropType< TypeTag, Properties::GridVolumeVariables >::LocalView ElementVolumeVariables
Definition porousmediumflow/3pwateroil/localresidual.hh:46
GetPropType< TypeTag, Properties::FluxVariables > FluxVariables
Definition porousmediumflow/3pwateroil/localresidual.hh:41
Dumux::NumEqVector< GetPropType< TypeTag, Properties::PrimaryVariables > > NumEqVector
Definition porousmediumflow/3pwateroil/localresidual.hh:40
GetPropType< TypeTag, Properties::FluidSystem > FluidSystem
Definition porousmediumflow/3pwateroil/localresidual.hh:49
DiscretizationDefaultLocalOperator< TypeTag > ParentType
Definition porousmediumflow/3pwateroil/localresidual.hh:34
typename FVElementGeometry::SubControlVolumeFace SubControlVolumeFace
Definition porousmediumflow/3pwateroil/localresidual.hh:39
GetPropType< TypeTag, Properties::EnergyLocalResidual > EnergyLocalResidual
Definition porousmediumflow/3pwateroil/localresidual.hh:48
static constexpr bool useMoles
Property that defines whether mole or mass fractions are used.
Definition porousmediumflow/3pwateroil/localresidual.hh:67
typename FVElementGeometry::SubControlVolume SubControlVolume
Definition porousmediumflow/3pwateroil/localresidual.hh:38
GetPropType< TypeTag, Properties::Scalar > Scalar
Definition porousmediumflow/3pwateroil/localresidual.hh:36
Defines all properties used in Dumux.
The default local operator than can be specialized for each discretization scheme.
typename NumEqVectorTraits< PrimaryVariables >::type NumEqVector
A vector with the same size as numbers of equations This is the default implementation and has to be ...
Definition numeqvector.hh:34
VolumeVariables::PrimaryVariables::value_type massOrMoleFraction(const VolumeVariables &volVars, ReferenceSystemFormulation referenceSys, const int phaseIdx, const int compIdx)
returns the mass or mole fraction to be used in Fick's law based on the reference system
Definition referencesystemformulation.hh:54
@ massAveraged
Definition referencesystemformulation.hh:34
constexpr auto getPropValue()
get the value data member of a property
Definition propertysystem.hh:310
typename GetProp< TypeTag, Property >::type GetPropType
get the type alias defined in the property
Definition propertysystem.hh:296
Definition adapt.hh:17
typename Detail::DiscretizationDefaultLocalOperator< TypeTag >::type DiscretizationDefaultLocalOperator
Definition defaultlocaloperator.hh:26
A helper to deduce a vector with the same size as numbers of equations.
The reference frameworks and formulations available for splitting total fluxes into a advective and d...