volumevariables.hh 7.35 KB
Newer Older
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
// -*- 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    *
 *   the Free Software Foundation, either version 2 of the License, or       *
 *   (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
 *
 * \brief This file contains the diffusion module for the vertex data
Thomas Fetzer's avatar
Thomas Fetzer committed
23
 *        of the fully coupled MpNc model
24
25
26
27
28
 */
#ifndef DUMUX_MPNC_DIFFUSION_VOLUME_VARIABLES_HH
#define DUMUX_MPNC_DIFFUSION_VOLUME_VARIABLES_HH

#include <dumux/common/valgrind.hh>
29
#include <dumux/implicit/mpnc/mpncproperties.hh>
30
31
32

namespace Dumux {

Thomas Fetzer's avatar
Thomas Fetzer committed
33
34
35
36
/*!
 * \brief Variables for the diffusive fluxes in the MpNc model within 
 *        a finite volume.
 */
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
template<class TypeTag, bool enableDiffusion>
class MPNCVolumeVariablesDiffusion
{
    typedef typename GET_PROP_TYPE(TypeTag, Problem) Problem;
    typedef typename GET_PROP_TYPE(TypeTag, VolumeVariables) VolumeVariables;
    typedef typename GET_PROP_TYPE(TypeTag, FluidSystem) FluidSystem;
    typedef typename GET_PROP_TYPE(TypeTag, Scalar) Scalar;
    typedef typename GET_PROP_TYPE(TypeTag, FluidState) FluidState;
    typedef typename FluidSystem::ParameterCache ParameterCache;

    enum { numComponents = GET_PROP_VALUE(TypeTag, NumComponents) };
    enum { wPhaseIdx = FluidSystem::wPhaseIdx };
    enum { nPhaseIdx = FluidSystem::nPhaseIdx };

public:
Dominik Riesterer's avatar
Dominik Riesterer committed
52
53
54
    /*!
     * \brief The constructor
     */
55
56
    MPNCVolumeVariablesDiffusion()
    {}
Dominik Riesterer's avatar
Dominik Riesterer committed
57
58
59
60
61
62
63
64
    /*!
     * \brief update
     *
     * \param fluidState An arbitrary fluid state
     * \param paramCache Container for cache parameters
     * \param volVars The volume variables
     * \param problem  The problem
     */
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
    void update(FluidState &fluidState,
                ParameterCache &paramCache,
                const VolumeVariables &volVars,
                const Problem &problem)
    {
        Valgrind::SetUndefined(*this);

        // diffusion coefficents in liquid
        diffCoeffL_[0] = 0.0;
        for (int compIdx = 1; compIdx < numComponents; ++compIdx) {
            diffCoeffL_[compIdx] =
                FluidSystem::binaryDiffusionCoefficient(fluidState,
                                                        paramCache,
                                                        wPhaseIdx,
                                                        0,
                                                        compIdx);
        }
        Valgrind::CheckDefined(diffCoeffL_);

        // diffusion coefficents in gas
        for (int compIIdx = 0; compIIdx < numComponents; ++compIIdx) {
            diffCoeffG_[compIIdx][compIIdx] = 0;
            for (int compJIdx = compIIdx + 1; compJIdx < numComponents; ++compJIdx) {
                diffCoeffG_[compIIdx][compJIdx] =
                        FluidSystem::binaryDiffusionCoefficient(fluidState,
                                                                paramCache,
                                                                nPhaseIdx,
                                                                compIIdx,
                                                                compJIdx);

                // fill the symmetric part of the diffusion coefficent
                // matrix
                diffCoeffG_[compJIdx][compIIdx] = diffCoeffG_[compIIdx][compJIdx];
            }
        }
        Valgrind::CheckDefined(diffCoeffG_);
    }

Dominik Riesterer's avatar
Dominik Riesterer committed
103
104
105
    /*!
     * \brief The binary diffusion coefficient for each fluid phase.
     */
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
    Scalar diffCoeff(const unsigned int phaseIdx,
                     const unsigned int compIIdx,
                     const unsigned int compJIdx) const
    {
        if (phaseIdx == nPhaseIdx)
            // TODO: tensorial diffusion coefficients
            return diffCoeffG_[compIIdx][compJIdx];

        const unsigned int i = std::min(compIIdx, compJIdx);
        const unsigned int j = std::max(compIIdx, compJIdx);
        if (i != 0)
            return 0;
        return diffCoeffL_[j];
    }

    /*!
     * \brief If running under valgrind this produces an error message
     *        if some of the object's attributes is undefined.
     */
    void checkDefined() const
    {
        Valgrind::CheckDefined(diffCoeffL_);
        Valgrind::CheckDefined(diffCoeffG_);
    }


protected:
    // the diffusion coefficients for the porous medium for the
    // liquid phase
    Scalar diffCoeffL_[numComponents];

    // the diffusion coefficients for the porous medium for the
    // gas phase
    Scalar diffCoeffG_[numComponents][numComponents];
};

// dummy class for the case where diffusion is disabled
template<class TypeTag>
class MPNCVolumeVariablesDiffusion<TypeTag, false>
{
    typedef typename GET_PROP_TYPE(TypeTag, Problem) Problem;
    typedef typename GET_PROP_TYPE(TypeTag, VolumeVariables) VolumeVariables;
    typedef typename GET_PROP_TYPE(TypeTag, Scalar) Scalar;
    typedef typename GET_PROP_TYPE(TypeTag, FluidSystem) FluidSystem;
    typedef typename GET_PROP_TYPE(TypeTag, FluidState) FluidState;
    typedef typename FluidSystem::ParameterCache ParameterCache;

public:
Dominik Riesterer's avatar
Dominik Riesterer committed
154
155
156
    /*!
     * \brief The constructor
     */
157
158
    MPNCVolumeVariablesDiffusion()
    {}
Dominik Riesterer's avatar
Dominik Riesterer committed
159
160
161
162
163
164
165
166
    /*!
     * \brief update
     *
     * \param fluidState An arbitrary fluid state
     * \param paramCache Container for cache parameters
     * \param volVars The volume variables
     * \param problem  The problem
     */
167
168
169
170
171
    void update(FluidState &fluidState,
                ParameterCache &paramCache,
                const VolumeVariables &volVars,
                const Problem &problem)
    { }
Dominik Riesterer's avatar
Dominik Riesterer committed
172
173
174
175
    /*!
     * \brief The binary diffusion coefficient for each component in the fluid phase.
     * \param compIdx The local index of the components
     */
176
177
    Scalar diffCoeffL(const unsigned int compIdx) const
    { return 0; }
Dominik Riesterer's avatar
Dominik Riesterer committed
178
179
180
181
182
    /*!
     * \brief The binary diffusion coefficient for each component in the gas phase.
     * \param compIIdx The local index of the first component in the phase
     * \param compIJdx The local index of the second component in the phase
     */
183
184
185
186
187
188
189
190
191
192
193
194
195
196
    Scalar diffCoeffG(const unsigned int compIIdx, const unsigned int compJIdx) const
    { return 0; }

    /*!
     * \brief If running under valgrind this produces an error message
     *        if some of the object's attributes is undefined.
     */
    void checkDefined() const
    { }
};

}

#endif // DUMUX_MPNC_DIFFUSION_VOLUME_VARIABLES_HH