3.3.0
DUNE for Multi-{Phase, Component, Scale, Physics, ...} flow and transport in porous media
freeflow/nonisothermal/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 *****************************************************************************/
24#ifndef DUMUX_FREE_FLOW_ENERGY_LOCAL_RESIDUAL_HH
25#define DUMUX_FREE_FLOW_ENERGY_LOCAL_RESIDUAL_HH
26
29
30namespace Dumux {
31
32// forward declaration
33template<class GridGeometry, class FluxVariables, DiscretizationMethod discMethod, bool enableEneryBalance, bool isCompositional>
35
41template<class GridGeometry, class FluxVariables, bool enableEneryBalance, bool isCompositional>
44 FluxVariables,
45 GridGeometry::discMethod,
46 enableEneryBalance, isCompositional>;
47
52template<class GridGeometry, class FluxVariables, DiscretizationMethod discMethod, bool isCompositional>
53class FreeFlowEnergyLocalResidualImplementation<GridGeometry, FluxVariables, discMethod, false, isCompositional>
54{
55public:
56
58 template <typename... Args>
59 static void fluidPhaseStorage(Args&&... args)
60 {}
61
63 template <typename... Args>
64 static void heatFlux(Args&&... args)
65 {}
66};
67
72template<class GridGeometry, class FluxVariables>
74 FluxVariables,
76 true, false>
77{
78 using Element = typename GridGeometry::GridView::template Codim<0>::Entity;
79 using FVElementGeometry = typename GridGeometry::LocalView;
80 using SubControlVolumeFace = typename FVElementGeometry::SubControlVolumeFace;
81
82public:
83
85 template<class NumEqVector, class VolumeVariables>
86 static void fluidPhaseStorage(NumEqVector& storage,
87 const VolumeVariables& volVars)
88 {
89 static constexpr auto localEnergyBalanceIdx = NumEqVector::dimension - 1;
90 storage[localEnergyBalanceIdx] += volVars.density() * volVars.internalEnergy();
91 }
92
94 template<class NumEqVector, class Problem, class ElementVolumeVariables, class ElementFaceVariables>
95 static void heatFlux(NumEqVector& flux,
96 const Problem& problem,
97 const Element &element,
98 const FVElementGeometry& fvGeometry,
99 const ElementVolumeVariables& elemVolVars,
100 const ElementFaceVariables& elemFaceVars,
101 const SubControlVolumeFace& scvf)
102 {
103 static constexpr auto localEnergyBalanceIdx = NumEqVector::dimension - 1;
104
105 auto upwindTerm = [](const auto& volVars) { return volVars.density() * volVars.enthalpy(); };
106 flux[localEnergyBalanceIdx] += FluxVariables::advectiveFluxForCellCenter(problem,
107 elemVolVars,
108 elemFaceVars,
109 scvf,
110 upwindTerm);
111
112 flux[localEnergyBalanceIdx] += FluxVariables::HeatConductionType::flux(problem,
113 element,
114 fvGeometry,
115 elemVolVars,
116 scvf);
117 }
118};
119
124template<class GridGeometry, class FluxVariables>
126 FluxVariables,
128 true, true>
129 : public FreeFlowEnergyLocalResidualImplementation<GridGeometry,
130 FluxVariables,
131 DiscretizationMethod::staggered,
132 true, false>
133{
135 FluxVariables,
137 true, false>;
138 using Element = typename GridGeometry::GridView::template Codim<0>::Entity;
139 using FVElementGeometry = typename GridGeometry::LocalView;
140 using SubControlVolumeFace = typename FVElementGeometry::SubControlVolumeFace;
141
142public:
144 template<class NumEqVector, class Problem, class ElementVolumeVariables, class ElementFaceVariables>
145 static void heatFlux(NumEqVector& flux,
146 const Problem& problem,
147 const Element &element,
148 const FVElementGeometry& fvGeometry,
149 const ElementVolumeVariables& elemVolVars,
150 const ElementFaceVariables& elemFaceVars,
151 const SubControlVolumeFace& scvf)
152 {
153 ParentType::heatFlux(flux, problem, element, fvGeometry, elemVolVars, elemFaceVars, scvf);
154
155 static constexpr auto localEnergyBalanceIdx = NumEqVector::dimension - 1;
156 auto diffusiveFlux = FluxVariables::MolecularDiffusionType::flux(problem, element, fvGeometry, elemVolVars, scvf);
157 for (int compIdx = 0; compIdx < FluxVariables::numComponents; ++compIdx)
158 {
159 const bool insideIsUpstream = scvf.directionSign() == sign(diffusiveFlux[compIdx]);
160 const auto& upstreamVolVars = insideIsUpstream ? elemVolVars[scvf.insideScvIdx()] : elemVolVars[scvf.outsideScvIdx()];
161
162 if (FluxVariables::MolecularDiffusionType::referenceSystemFormulation() == ReferenceSystemFormulation::massAveraged)
163 flux[localEnergyBalanceIdx] += diffusiveFlux[compIdx] * upstreamVolVars.componentEnthalpy(compIdx);
164 else
165 flux[localEnergyBalanceIdx] += diffusiveFlux[compIdx] * upstreamVolVars.componentEnthalpy(compIdx)* elemVolVars[scvf.insideScvIdx()].molarMass(compIdx);
166
167 }
168 }
169};
170
171} // end namespace Dumux
172
173#endif
The available discretization methods in Dumux.
The reference frameworks and formulations available for splitting total fluxes into a advective and d...
DiscretizationMethod
The available discretization methods in Dumux.
Definition: method.hh:37
constexpr int sign(const ValueType &value) noexcept
Sign or signum function.
Definition: math.hh:627
Definition: adapt.hh:29
Definition: freeflow/nonisothermal/localresidual.hh:34
static void heatFlux(Args &&... args)
do nothing for the isothermal case
Definition: freeflow/nonisothermal/localresidual.hh:64
static void fluidPhaseStorage(Args &&... args)
do nothing for the isothermal case
Definition: freeflow/nonisothermal/localresidual.hh:59
Specialization for staggered one-phase, non-isothermal models.
Definition: freeflow/nonisothermal/localresidual.hh:77
static void heatFlux(NumEqVector &flux, const Problem &problem, const Element &element, const FVElementGeometry &fvGeometry, const ElementVolumeVariables &elemVolVars, const ElementFaceVariables &elemFaceVars, const SubControlVolumeFace &scvf)
The convective and conductive heat fluxes in the fluid phase.
Definition: freeflow/nonisothermal/localresidual.hh:95
static void fluidPhaseStorage(NumEqVector &storage, const VolumeVariables &volVars)
The energy storage in the fluid phase.
Definition: freeflow/nonisothermal/localresidual.hh:86
static void heatFlux(NumEqVector &flux, const Problem &problem, const Element &element, const FVElementGeometry &fvGeometry, const ElementVolumeVariables &elemVolVars, const ElementFaceVariables &elemFaceVars, const SubControlVolumeFace &scvf)
The convective and conductive heat fluxes in the fluid phase.
Definition: freeflow/nonisothermal/localresidual.hh:145