volumevariables.hh 9.09 KB
Newer Older
1
2
3
4
5
6
7
// -*- mode: C++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*-
// vi: set et ts=4 sw=4 sts=4:
/*****************************************************************************
 *   See the file COPYING for full copying permissions.                      *
 *                                                                           *
 *   This program is free software: you can redistribute it and/or modify    *
 *   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
10
11
12
13
14
15
16
17
18
19
20
 *   (at your option) any later version.                                     *
 *                                                                           *
 *   This program is distributed in the hope that it will be useful,         *
 *   but WITHOUT ANY WARRANTY; without even the implied warranty of          *
 *   MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the            *
 *   GNU General Public License for more details.                            *
 *                                                                           *
 *   You should have received a copy of the GNU General Public License       *
 *   along with this program.  If not, see <http://www.gnu.org/licenses/>.   *
 *****************************************************************************/
/*!
 * \file
21
 * \ingroup TracerModel
22
 * \brief Quantities required by the tracer model in a control volume.
23
24
25
26
 */
#ifndef DUMUX_TRACER_VOLUME_VARIABLES_HH
#define DUMUX_TRACER_VOLUME_VARIABLES_HH

27
#include <type_traits>
Bernd Flemisch's avatar
Bernd Flemisch committed
28
#include <array>
29
#include <dune/common/std/type_traits.hh>
Bernd Flemisch's avatar
Bernd Flemisch committed
30

31
#include <dumux/porousmediumflow/volumevariables.hh>
Timo Koch's avatar
Timo Koch committed
32
#include <dumux/material/solidstates/updatesolidvolumefractions.hh>
33

34
namespace Dumux {
35

36
37
38
39
40
41
42
43
44
45
46
47
namespace Detail {
// helper structs and functions detecting if the user-defined spatial params class
// has user-specified functions saturation() for multi-phase tracer.
template <typename T, typename ...Ts>
using SaturationDetector = decltype(std::declval<T>().spatialParams().saturation(std::declval<Ts>()...));

template<class T, typename ...Args>
static constexpr bool hasSaturation()
{ return Dune::Std::is_detected<SaturationDetector, T, Args...>::value; }

} // end namespace Detail

48
49
50
51
52
/*!
 * \ingroup TracerModel
 * \brief Contains the quantities which are constant within a
 *        finite volume for the tracer model.
 */
53
54
template <class Traits>
class TracerVolumeVariables
55
: public PorousMediumFlowVolumeVariables<Traits>
56
{
57
    using ParentType = PorousMediumFlowVolumeVariables<Traits>;
58
59
    using Scalar = typename Traits::PrimaryVariables::value_type;
    static constexpr bool useMoles = Traits::ModelTraits::useMoles();
60
61

public:
62
    //! Export fluid system type
63
    using FluidSystem = typename Traits::FluidSystem;
64
    using SolidState = typename Traits::SolidState;
65
66

    /*!
67
     * \brief Updates all quantities for a given control volume.
68
69
70
71
72
73
     *
     * \param elemSol A vector containing all primary variables connected to the element
     * \param problem The object specifying the problem which ought to
     *                be simulated
     * \param element An element which contains part of the control volume
     * \param scv The sub-control volume
74
     */
75
76
    template<class ElemSol, class Problem, class Element, class Scv>
    void update(const ElemSol &elemSol,
77
78
                const Problem &problem,
                const Element &element,
79
                const Scv &scv)
80
81
82
83
    {
        // update parent type sets primary variables
        ParentType::update(elemSol, problem, element, scv);

84
        updateSolidVolumeFractions(elemSol, problem, element, scv, solidState_, ParentType::numFluidComponents());
85
        // dispersivity_ = problem.spatialParams().dispersivity(element, scv, elemSol);
86
87
88
89

        // the spatial params special to the tracer model
        fluidDensity_ = problem.spatialParams().fluidDensity(element, scv);
        fluidMolarMass_ = problem.spatialParams().fluidMolarMass(element, scv);
90
        fluidSaturation_ = saturation_(problem, element, scv);
91

92
        for (int compIdx = 0; compIdx < ParentType::numFluidComponents(); ++compIdx)
93
        {
94
            moleOrMassFraction_[compIdx] = this->priVars()[compIdx];
95
            diffCoeff_[compIdx] = FluidSystem::binaryDiffusionCoefficient(compIdx, problem, element, scv);
96
97
98
99
        }
    }

    /*!
100
     * \brief Returns the density \f$\mathrm{[kg/m^3]}\f$ the of the fluid phase.
101
102
     *
     * We always forward to the fluid state with the phaseIdx property (see class description).
103
     *
104
     * \param phaseIdx The phase index
105
     */
106
    Scalar density(int phaseIdx = 0) const
107
108
    { return fluidDensity_; }

109
110
111
112
113
114
    /*!
     * \brief Returns the phase state for the control volume.
     */
    const SolidState &solidState() const
    { return solidState_; }

115
    /*!
116
     * \brief Returns the saturation.
117
118
     *
     * This method is here for compatibility reasons with other models. The saturation
119
120
121
122
     * is always 1.0 in a one-phasic context, if two-phases or richards are considered,
     * the spatialParams serve as way to pass the saturation from the main-file to the
     * volVars and then to the localresidual for the tracer model.

123
     * \param phaseIdx The phase index
124
     */
125
    Scalar saturation(int phaseIdx = 0) const
126
    { return fluidSaturation_ ; }
127
128

    /*!
129
     * \brief Returns the mobility.
130
131
132
     *
     * This method is here for compatibility reasons with other models. The mobility is always 1
     * for one-phasic models where the velocity field is given
133
     *
134
     * \param phaseIdx The phase index
135
     */
136
    Scalar mobility(int phaseIdx = 0) const
137
138
139
    { return 1.0; }

    /*!
140
     * \brief Returns the molar density \f$\mathrm{[mol/m^3]}\f$ the of the fluid phase.
141
     *
142
     * \param phaseIdx The phase index
143
     */
144
    Scalar molarDensity(int phaseIdx = 0) const
145
146
147
    { return fluidDensity_/fluidMolarMass_; }

    /*!
148
     * \brief Returns the mole fraction \f$\mathrm{[mol/mol]}\f$ of a component in the phase.
149
     *
150
     * \param phaseIdx The phase index
151
152
     * \param compIdx The index of the component
     */
153
    Scalar moleFraction(int phaseIdx, int compIdx) const
154
    { return useMoles ? moleOrMassFraction_[compIdx] : moleOrMassFraction_[compIdx]/FluidSystem::molarMass(compIdx)*fluidMolarMass_; }
155
156

    /*!
157
     * \brief Returns the mass fraction \f$\mathrm{[kg/kg]}\f$ of a component in the phase.
158
     *
159
     * \param phaseIdx The phase index
160
161
     * \param compIdx The index of the component
     */
162
    Scalar massFraction(int phaseIdx, int compIdx) const
163
    { return useMoles ? moleOrMassFraction_[compIdx]*FluidSystem::molarMass(compIdx)/fluidMolarMass_ : moleOrMassFraction_[compIdx]; }
164
165

    /*!
166
     * \brief Returns the concentration \f$\mathrm{[mol/m^3]}\f$  of a component in the phase.
167
     *
168
     * \param phaseIdx The phase index
169
170
     * \param compIdx The index of the component
     */
171
172
    Scalar molarity(int phaseIdx, int compIdx) const
    { return moleFraction(phaseIdx, compIdx)*molarDensity(); }
173
174

    /*!
175
     * \brief Returns the binary diffusion coefficient \f$\mathrm{[m^2/s]}\f$ in the fluid.
176
     *
177
     * \param phaseIdx The phase index
178
     * \param compIdx The index of the component
179
     */
180
    Scalar diffusionCoefficient(int phaseIdx, int compIdx) const
181
182
    { return diffCoeff_[compIdx]; }

183
184
185
186
187
188
    // /*!
    //  * \brief Returns the dispersivity of the fluid's streamlines.
    //  * \todo implement me
    //  */
    // const DispersivityType &dispersivity() const
    // { return dispersivity_; }
189
190
191
192
193

    /*!
     * \brief Return the average porosity \f$\mathrm{[-]}\f$ within the control volume.
     */
    Scalar porosity() const
194
    { return solidState_.porosity(); }
195
196

protected:
197
    SolidState solidState_;
198
    Scalar fluidDensity_, fluidMolarMass_;
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
    Scalar fluidSaturation_ = 1.0;
    /*!
     * \brief Gets the saturation in an scv.
     *
     * \param problem the problem to solve
     * \param element the element (codim-0-entity) the scv belongs to
     * \param scv the sub control volume
     * \note this gets selected if the user uses the multiphase tracer
     */
     template<class Problem, class Element, class Scv,
              std::enable_if_t<Detail::hasSaturation<Problem, Element, Scv>(), int> = 0>
     Scalar saturation_(const Problem& problem,
                        const Element& element,
                        const Scv& scv)
     { return problem.spatialParams().saturation(element, scv); }

    /*!
     * \brief Gets the saturation in an scv.
     *
     * \param problem the problem to solve
     * \param element the element (codim-0-entity) the scv belongs to
     * \param scv the sub control volume
     * \note this gets selected if the user a single phase tracer
     */
    template<class Problem, class Element, class Scv,
             std::enable_if_t<!Detail::hasSaturation<Problem, Element, Scv>(), int> = 0>
    Scalar saturation_(const Problem& problem,
                       const Element &element,
                       const Scv &scv)
    { return 1.0; }

230
    // DispersivityType dispersivity_;
231
232
    std::array<Scalar, ParentType::numFluidComponents()> diffCoeff_;
    std::array<Scalar, ParentType::numFluidComponents()> moleOrMassFraction_;
233
234
235
236
237
};

} // end namespace Dumux

#endif