3.1-git
DUNE for Multi-{Phase, Component, Scale, Physics, ...} flow and transport in porous media
porousmediumflow/3p3c/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 * See the file COPYING for full copying permissions. *
5 * *
6 * This program 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 3 of the License, or *
9 * (at your option) any later version. *
10 * *
11 * This program 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 this program. If not, see <http://www.gnu.org/licenses/>. *
18 *****************************************************************************/
26#ifndef DUMUX_3P3C_LOCAL_RESIDUAL_HH
27#define DUMUX_3P3C_LOCAL_RESIDUAL_HH
28
31
32namespace Dumux
33{
39template<class TypeTag>
40class ThreePThreeCLocalResidual : public GetPropType<TypeTag, Properties::BaseLocalResidual>
41{
45 using FVElementGeometry = typename GetPropType<TypeTag, Properties::GridGeometry>::LocalView;
46 using SubControlVolume = typename FVElementGeometry::SubControlVolume;
47 using SubControlVolumeFace = typename FVElementGeometry::SubControlVolumeFace;
50 using ElementFluxVariablesCache = typename GetPropType<TypeTag, Properties::GridFluxVariablesCache>::LocalView;
53 using Element = typename GridView::template Codim<0>::Entity;
54 using ElementVolumeVariables = typename GetPropType<TypeTag, Properties::GridVolumeVariables>::LocalView;
58
59 enum {
62
63 contiWEqIdx = Indices::conti0EqIdx + FluidSystem::wPhaseIdx,
64 contiNEqIdx = Indices::conti0EqIdx + FluidSystem::nPhaseIdx,
65 contiGEqIdx = Indices::conti0EqIdx + FluidSystem::gPhaseIdx,
66
67 wPhaseIdx = FluidSystem::wPhaseIdx,
68 nPhaseIdx = FluidSystem::nPhaseIdx,
69 gPhaseIdx = FluidSystem::gPhaseIdx,
70
71 wCompIdx = FluidSystem::wCompIdx,
72 nCompIdx = FluidSystem::nCompIdx,
73 gCompIdx = FluidSystem::gCompIdx
74 };
75
76public:
77
78 using ParentType::ParentType;
79
91 NumEqVector computeStorage(const Problem& problem,
92 const SubControlVolume& scv,
93 const VolumeVariables& volVars) const
94 {
95 NumEqVector storage(0.0);
96
97 // compute storage term of all components within all phases
98 for (int phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx)
99 {
100 for (int compIdx = 0; compIdx < numComponents; ++compIdx)
101 {
102 auto eqIdx = Indices::conti0EqIdx + compIdx;
103 storage[eqIdx] += volVars.porosity()
104 * volVars.saturation(phaseIdx)
105 * volVars.molarDensity(phaseIdx)
106 * volVars.moleFraction(phaseIdx, compIdx);
107 }
108
109 // The energy storage in the fluid phase with index phaseIdx
110 EnergyLocalResidual::fluidPhaseStorage(storage, scv, volVars, phaseIdx);
111 }
112
113 // The energy storage in the solid matrix
114 EnergyLocalResidual::solidPhaseStorage(storage, scv, volVars);
115
116 return storage;
117 }
118
130 NumEqVector computeFlux(const Problem& problem,
131 const Element& element,
132 const FVElementGeometry& fvGeometry,
133 const ElementVolumeVariables& elemVolVars,
134 const SubControlVolumeFace& scvf,
135 const ElementFluxVariablesCache& elemFluxVarsCache) const
136 {
137 FluxVariables fluxVars;
138 fluxVars.init(problem, element, fvGeometry, elemVolVars, scvf, elemFluxVarsCache);
139 static constexpr auto referenceSystemFormulation = FluxVariables::MolecularDiffusionType::referenceSystemFormulation();
140
141 // get upwind weights into local scope
142 NumEqVector flux(0.0);
143
144 // advective fluxes
145 for (int phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx)
146 {
147 for (int compIdx = 0; compIdx < numComponents; ++compIdx)
148 {
149 auto upwindTerm = [phaseIdx, compIdx](const VolumeVariables& volVars)
150 { return volVars.molarDensity(phaseIdx)*volVars.moleFraction(phaseIdx, compIdx)*volVars.mobility(phaseIdx); };
151
152 // get equation index
153 auto eqIdx = Indices::conti0EqIdx + compIdx;
154 flux[eqIdx] += fluxVars.advectiveFlux(phaseIdx, upwindTerm);
155 }
156
157 // Add advective phase energy fluxes. For isothermal model the contribution is zero.
158 EnergyLocalResidual::heatConvectionFlux(flux, fluxVars, phaseIdx);
159 }
160
161 // Add diffusive energy fluxes. For isothermal model the contribution is zero.
162 EnergyLocalResidual::heatConductionFlux(flux, fluxVars);
163
164 // diffusive fluxes
165 const auto diffusionFluxesWPhase = fluxVars.molecularDiffusionFlux(wPhaseIdx);
166 Scalar jGW = diffusionFluxesWPhase[gCompIdx];
167 Scalar jNW = diffusionFluxesWPhase[nCompIdx];
168
169 //check for the reference system and adapt units of the diffusive flux accordingly.
170 if (referenceSystemFormulation == ReferenceSystemFormulation::massAveraged)
171 {
172 jGW /= FluidSystem::molarMass(gCompIdx);
173 jNW /= FluidSystem::molarMass(nCompIdx);
174 }
175 const Scalar jWW = -(jGW+jNW);
176
177 const auto diffusionFluxesGPhase = fluxVars.molecularDiffusionFlux(gPhaseIdx);
178 Scalar jWG = diffusionFluxesGPhase[wCompIdx];
179 Scalar jNG = diffusionFluxesGPhase[nCompIdx];
180
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 const Scalar jGG = -(jWG+jNG);
188
189 // At the moment we do not consider diffusion in the NAPL phase
190 const Scalar jWN = 0.0;
191 const Scalar jGN = 0.0;
192 const Scalar jNN = 0.0;
193
194 flux[contiWEqIdx] += jWW+jWG+jWN;
195 flux[contiNEqIdx] += jNW+jNG+jNN;
196 flux[contiGEqIdx] += jGW+jGG+jGN;
197
198 return flux;
199 }
200};
201
202} // end namespace Dumux
203
204#endif
The reference frameworks and formulations available for splitting total fluxes into a advective and d...
make the local view function available whenever we use the grid geometry
Definition: adapt.hh:29
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:149
Element-wise calculation of the Jacobian matrix for problems using the three-phase three-component fu...
Definition: porousmediumflow/3p3c/localresidual.hh:41
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/3p3c/localresidual.hh:91
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/3p3c/localresidual.hh:130
Declares all properties used in Dumux.