3.3.0
DUNE for Multi-{Phase, Component, Scale, Physics, ...} flow and transport in porous media
freeflow/compositional/staggered/fluxvariables.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_FREEFLOW_NC_STAGGERED_FLUXVARIABLES_HH
25#define DUMUX_FREEFLOW_NC_STAGGERED_FLUXVARIABLES_HH
26
27#include <numeric>
33
34namespace Dumux {
35
41// forward declaration
42template<class TypeTag, DiscretizationMethod discMethod>
43class FreeflowNCFluxVariablesImpl;
44
45template<class TypeTag>
47: public NavierStokesFluxVariables<TypeTag>
48{
52 using FVElementGeometry = typename GridGeometry::LocalView;
53 using Element = typename GridGeometry::GridView::template Codim<0>::Entity;
54 using SubControlVolumeFace = typename FVElementGeometry::SubControlVolumeFace;
55 using CellCenterPrimaryVariables = GetPropType<TypeTag, Properties::CellCenterPrimaryVariables>;
58
59public:
60 static constexpr auto numComponents = ModelTraits::numFluidComponents();
61 static constexpr bool useMoles = ModelTraits::useMoles();
63
67 template<class ElementVolumeVariables, class ElementFaceVariables, class FluxVariablesCache>
68 CellCenterPrimaryVariables computeMassFlux(const Problem& problem,
69 const Element& element,
70 const FVElementGeometry& fvGeometry,
71 const ElementVolumeVariables& elemVolVars,
72 const ElementFaceVariables& elemFaceVars,
73 const SubControlVolumeFace& scvf,
74 const FluxVariablesCache& fluxVarsCache)
75 {
76 CellCenterPrimaryVariables flux(0.0);
77
78 const auto diffusiveFluxes = MolecularDiffusionType::flux(problem, element, fvGeometry, elemVolVars, scvf);
79
80 static constexpr auto referenceSystemFormulation = MolecularDiffusionType::referenceSystemFormulation();
81
82 for (int compIdx = 0; compIdx < numComponents; ++compIdx)
83 {
84 auto upwindTerm = [compIdx](const auto& volVars)
85 {
86 const auto density = useMoles ? volVars.molarDensity() : volVars.density();
87 const auto fraction = useMoles ? volVars.moleFraction(compIdx) : volVars.massFraction(compIdx);
88 return density * fraction;
89 };
90
91 flux[compIdx] = ParentType::advectiveFluxForCellCenter(problem, elemVolVars, elemFaceVars, scvf, upwindTerm);
92
93 //check for the reference system and adapt units of the diffusive flux accordingly.
94 if (referenceSystemFormulation == ReferenceSystemFormulation::massAveraged)
95 {
96 flux[compIdx] += useMoles ? diffusiveFluxes[compIdx]/FluidSystem::molarMass(compIdx) : diffusiveFluxes[compIdx];
97 }
98 else if (referenceSystemFormulation == ReferenceSystemFormulation::molarAveraged)
99 flux[compIdx] += useMoles ? diffusiveFluxes[compIdx] : diffusiveFluxes[compIdx]*FluidSystem::molarMass(compIdx);
100 else
101 DUNE_THROW(Dune::NotImplemented, "other reference systems than mass and molar averaged are not implemented");
102 }
103
104 // in case one balance is substituted by the total mass balance
105 if (ModelTraits::replaceCompEqIdx() < numComponents)
106 {
107 flux[ModelTraits::replaceCompEqIdx()] = std::accumulate(flux.begin(), flux.end(), 0.0);
108 }
109
110 return flux;
111 }
112};
113
114} // end namespace
115
116#endif
The available discretization methods in Dumux.
Base class for the flux variables living on a sub control volume face.
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
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
std::string density(int phaseIdx) noexcept
I/O name of density for multiphase systems.
Definition: name.hh:65
The flux variables class for the multi-component free-flow model using the staggered grid discretizat...
Definition: freeflow/compositional/fluxvariables.hh:34
GetPropType< TypeTag, Properties::MolecularDiffusionType > MolecularDiffusionType
Definition: freeflow/compositional/staggered/fluxvariables.hh:62
CellCenterPrimaryVariables computeMassFlux(const Problem &problem, const Element &element, const FVElementGeometry &fvGeometry, const ElementVolumeVariables &elemVolVars, const ElementFaceVariables &elemFaceVars, const SubControlVolumeFace &scvf, const FluxVariablesCache &fluxVarsCache)
Computes the flux for the cell center residual.
Definition: freeflow/compositional/staggered/fluxvariables.hh:68
The flux variables class for the Navier-Stokes model using the staggered grid discretization.
Definition: freeflow/navierstokes/fluxvariables.hh:35
Declares all properties used in Dumux.