version 3.10.0
Loading...
Searching...
No Matches
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// 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_3P3C_LOCAL_RESIDUAL_HH
15#define DUMUX_3P3C_LOCAL_RESIDUAL_HH
16
21
22namespace Dumux
23{
29template<class TypeTag>
32{
37 using FVElementGeometry = typename GridGeometry::LocalView;
38 using SubControlVolume = typename FVElementGeometry::SubControlVolume;
39 using SubControlVolumeFace = typename FVElementGeometry::SubControlVolumeFace;
42 using ElementFluxVariablesCache = typename GetPropType<TypeTag, Properties::GridFluxVariablesCache>::LocalView;
44 using GridView = typename GridGeometry::GridView;
45 using Element = typename GridView::template Codim<0>::Entity;
46 using ElementVolumeVariables = typename GetPropType<TypeTag, Properties::GridVolumeVariables>::LocalView;
50
51 enum {
54
55 contiWEqIdx = Indices::conti0EqIdx + FluidSystem::wPhaseIdx,
56 contiNEqIdx = Indices::conti0EqIdx + FluidSystem::nPhaseIdx,
57 contiGEqIdx = Indices::conti0EqIdx + FluidSystem::gPhaseIdx,
58
59 wPhaseIdx = FluidSystem::wPhaseIdx,
60 nPhaseIdx = FluidSystem::nPhaseIdx,
61 gPhaseIdx = FluidSystem::gPhaseIdx,
62
63 wCompIdx = FluidSystem::wCompIdx,
64 nCompIdx = FluidSystem::nCompIdx,
65 gCompIdx = FluidSystem::gCompIdx
66 };
67
68public:
69
70 using ParentType::ParentType;
71
83 NumEqVector computeStorage(const Problem& problem,
84 const SubControlVolume& scv,
85 const VolumeVariables& volVars) const
86 {
87 NumEqVector storage(0.0);
88
89 // compute storage term of all components within all phases
90 for (int phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx)
91 {
92 for (int compIdx = 0; compIdx < numComponents; ++compIdx)
93 {
94 auto eqIdx = Indices::conti0EqIdx + compIdx;
95 storage[eqIdx] += volVars.porosity()
96 * volVars.saturation(phaseIdx)
97 * volVars.molarDensity(phaseIdx)
98 * volVars.moleFraction(phaseIdx, compIdx);
99 }
100
101 // The energy storage in the fluid phase with index phaseIdx
102 EnergyLocalResidual::fluidPhaseStorage(storage, problem, scv, volVars, phaseIdx);
103 }
104
105 // The energy storage in the solid matrix
106 EnergyLocalResidual::solidPhaseStorage(storage, scv, volVars);
107
108 return storage;
109 }
110
122 NumEqVector computeFlux(const Problem& problem,
123 const Element& element,
124 const FVElementGeometry& fvGeometry,
125 const ElementVolumeVariables& elemVolVars,
126 const SubControlVolumeFace& scvf,
127 const ElementFluxVariablesCache& elemFluxVarsCache) const
128 {
129 FluxVariables fluxVars;
130 fluxVars.init(problem, element, fvGeometry, elemVolVars, scvf, elemFluxVarsCache);
131 static constexpr auto referenceSystemFormulation = FluxVariables::MolecularDiffusionType::referenceSystemFormulation();
132
133 // get upwind weights into local scope
134 NumEqVector flux(0.0);
135
136 // advective fluxes
137 for (int phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx)
138 {
139 for (int compIdx = 0; compIdx < numComponents; ++compIdx)
140 {
141 auto upwindTerm = [phaseIdx, compIdx](const VolumeVariables& volVars)
142 { return volVars.molarDensity(phaseIdx)*volVars.moleFraction(phaseIdx, compIdx)*volVars.mobility(phaseIdx); };
143
144 // get equation index
145 auto eqIdx = Indices::conti0EqIdx + compIdx;
146 flux[eqIdx] += fluxVars.advectiveFlux(phaseIdx, upwindTerm);
147 }
148
149 // Add advective phase energy fluxes. For isothermal model the contribution is zero.
150 EnergyLocalResidual::heatConvectionFlux(flux, fluxVars, phaseIdx);
151 }
152
153 // Add diffusive energy fluxes. For isothermal model the contribution is zero.
154 EnergyLocalResidual::heatConductionFlux(flux, fluxVars);
155
156 // diffusive fluxes
157 const auto diffusionFluxesWPhase = fluxVars.molecularDiffusionFlux(wPhaseIdx);
158 Scalar jGW = diffusionFluxesWPhase[gCompIdx];
159 Scalar jNW = diffusionFluxesWPhase[nCompIdx];
160 Scalar jWW = -(jGW+jNW);
161
162 //check for the reference system and adapt units of the diffusive flux accordingly.
163 if (referenceSystemFormulation == ReferenceSystemFormulation::massAveraged)
164 {
165 jGW /= FluidSystem::molarMass(gCompIdx);
166 jNW /= FluidSystem::molarMass(nCompIdx);
167 jWW /= FluidSystem::molarMass(wCompIdx);
168 }
169
170 const auto diffusionFluxesGPhase = fluxVars.molecularDiffusionFlux(gPhaseIdx);
171 Scalar jWG = diffusionFluxesGPhase[wCompIdx];
172 Scalar jNG = diffusionFluxesGPhase[nCompIdx];
173 Scalar jGG = -(jWG+jNG);
174
175 //check for the reference system and adapt units of the diffusive flux accordingly.
176 if (referenceSystemFormulation == ReferenceSystemFormulation::massAveraged)
177 {
178 jWG /= FluidSystem::molarMass(wCompIdx);
179 jNG /= FluidSystem::molarMass(nCompIdx);
180 jGG /= FluidSystem::molarMass(gCompIdx);
181 }
182
183 // At the moment we do not consider diffusion in the NAPL phase
184 const Scalar jWN = 0.0;
185 const Scalar jGN = 0.0;
186 const Scalar jNN = 0.0;
187
188 flux[contiWEqIdx] += jWW+jWG+jWN;
189 flux[contiNEqIdx] += jNW+jNG+jNN;
190 flux[contiGEqIdx] += jGW+jGG+jGN;
191
192 return flux;
193 }
194};
195
196} // end namespace Dumux
197
198#endif
Element-wise calculation of the Jacobian matrix for problems using the three-phase three-component fu...
Definition porousmediumflow/3p3c/localresidual.hh:32
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:83
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:122
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
@ massAveraged
Definition referencesystemformulation.hh:34
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...