lensspatialparams.hh 6.8 KB
Newer Older
1
2
// -*- mode: C++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*-
// vi: set et ts=4 sw=4 sts=4:
Bernd Flemisch's avatar
Bernd Flemisch committed
3
/*****************************************************************************
4
 *   See the file COPYING for full copying permissions.                      *
Bernd Flemisch's avatar
Bernd Flemisch committed
5
 *                                                                           *
6
 *   This program is free software: you can redistribute it and/or modify    *
Bernd Flemisch's avatar
Bernd Flemisch committed
7
 *   it under the terms of the GNU General Public License as published by    *
8
9
 *   the Free Software Foundation, either version 2 of the License, or       *
 *   (at your option) any later version.                                     *
Bernd Flemisch's avatar
Bernd Flemisch committed
10
 *                                                                           *
11
12
 *   This program is distributed in the hope that it will be useful,         *
 *   but WITHOUT ANY WARRANTY; without even the implied warranty of          *
13
 *   MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the            *
14
15
16
17
 *   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/>.   *
Bernd Flemisch's avatar
Bernd Flemisch committed
18
 *****************************************************************************/
19
20
21
22
/*!
 * \file
 *
 * \brief The spatial parameters for the LensProblem which uses the
23
 *        two-phase fully implicit model
24
 */
25
26
#ifndef DUMUX_LENS_SPATIAL_PARAMS_HH
#define DUMUX_LENS_SPATIAL_PARAMS_HH
Bernd Flemisch's avatar
Bernd Flemisch committed
27

28
#include <dumux/material/spatialparams/implicitspatialparams.hh>
Bernd Flemisch's avatar
Bernd Flemisch committed
29
30
31
32
#include <dumux/material/fluidmatrixinteractions/2p/regularizedvangenuchten.hh>
#include <dumux/material/fluidmatrixinteractions/2p/linearmaterial.hh>
#include <dumux/material/fluidmatrixinteractions/2p/efftoabslaw.hh>

33
#include <dumux/implicit/2p/2pmodel.hh>
Andreas Lauser's avatar
Andreas Lauser committed
34

Bernd Flemisch's avatar
Bernd Flemisch committed
35
36
37
namespace Dumux
{

38
39
//forward declaration
template<class TypeTag>
40
class LensSpatialParams;
41
42
43
44

namespace Properties
{
// The spatial parameters TypeTag
45
NEW_TYPE_TAG(LensSpatialParams);
46
47

// Set the spatial parameters
48
SET_TYPE_PROP(LensSpatialParams, SpatialParams, Dumux::LensSpatialParams<TypeTag>);
49
50

// Set the material Law
51
SET_PROP(LensSpatialParams, MaterialLaw)
52
53
54
55
{
private:
    // define the material law which is parameterized by effective
    // saturations
56
    typedef typename GET_PROP_TYPE(TypeTag, Scalar) Scalar;
57
58
59
60
61
62
    typedef RegularizedVanGenuchten<Scalar> EffectiveLaw;
public:
    // define the material law parameterized by absolute saturations
    typedef EffToAbsLaw<EffectiveLaw> type;
};
}
63
/*!
64
 * \ingroup TwoPModel
65
 * \ingroup ImplicitTestProblems
66
 * \brief The spatial parameters for the LensProblem which uses the
67
 *        two-phase fully implicit model
68
 */
Bernd Flemisch's avatar
Bernd Flemisch committed
69
template<class TypeTag>
70
class LensSpatialParams : public ImplicitSpatialParams<TypeTag>
Bernd Flemisch's avatar
Bernd Flemisch committed
71
{
72
    typedef ImplicitSpatialParams<TypeTag> ParentType;
73
74
75
    typedef typename GET_PROP_TYPE(TypeTag, Grid) Grid;
    typedef typename GET_PROP_TYPE(TypeTag, GridView) GridView;
    typedef typename GET_PROP_TYPE(TypeTag, Scalar) Scalar;
Bernd Flemisch's avatar
Bernd Flemisch committed
76
77
78
79
    typedef typename Grid::ctype CoordScalar;

    enum {
        dim=GridView::dimension,
80
        dimWorld=GridView::dimensionworld
Bernd Flemisch's avatar
Bernd Flemisch committed
81
82
83
84
85
    };

    typedef Dune::FieldVector<CoordScalar,dimWorld> GlobalPosition;

    typedef typename GridView::template Codim<0>::Entity Element;
86
    typedef typename GET_PROP_TYPE(TypeTag, FVElementGeometry) FVElementGeometry;
Bernd Flemisch's avatar
Bernd Flemisch committed
87
88

public:
89
    //get the material law from the property system
90
    typedef typename GET_PROP_TYPE(TypeTag, MaterialLaw) MaterialLaw;
Bernd Flemisch's avatar
Bernd Flemisch committed
91
92
    typedef typename MaterialLaw::Params MaterialLawParams;

93
94
95
96
97
    /*!
     * \brief The constructor
     *
     * \param gridView The grid view
     */
98
    LensSpatialParams(const GridView& gridView)
99
    : ParentType(gridView)
Bernd Flemisch's avatar
Bernd Flemisch committed
100
    {
101
102
103
104
            lensLowerLeft_[0]   = GET_RUNTIME_PARAM(TypeTag, Scalar, SpatialParams.LensLowerLeftX);
            lensLowerLeft_[1]   = GET_RUNTIME_PARAM(TypeTag, Scalar, SpatialParams.LensLowerLeftY);
            lensUpperRight_[0]  = GET_RUNTIME_PARAM(TypeTag, Scalar, SpatialParams.LensUpperRightX);
            lensUpperRight_[1]  = GET_RUNTIME_PARAM(TypeTag, Scalar, SpatialParams.LensUpperRightY);
105

Bernd Flemisch's avatar
Bernd Flemisch committed
106
107
108
109
110
111
        // residual saturations
        lensMaterialParams_.setSwr(0.18);
        lensMaterialParams_.setSnr(0.0);
        outerMaterialParams_.setSwr(0.05);
        outerMaterialParams_.setSnr(0.0);

112
113
114
        // parameters for the Van Genuchten law
        // alpha and n
        lensMaterialParams_.setVgAlpha(0.00045);
115
        lensMaterialParams_.setVgn(7.3);
116
        outerMaterialParams_.setVgAlpha(0.0037);
117
        outerMaterialParams_.setVgn(4.7);
118

119
120
        lensK_ = 9.05e-12;
        outerK_ = 4.6e-10;
Bernd Flemisch's avatar
Bernd Flemisch committed
121
122
123
    }

    /*!
124
     * \brief Returns the scalar intrinsic permeability \f$[m^2]\f$
Bernd Flemisch's avatar
Bernd Flemisch committed
125
     *
126
127
128
     * \param element The finite element
     * \param fvGeometry The finite volume geometry of the element
     * \param scvIdx The local index of the sub-control volume
Bernd Flemisch's avatar
Bernd Flemisch committed
129
     */
130
    Scalar intrinsicPermeability(const Element &element,
131
                                 const FVElementGeometry &fvGeometry,
132
                                 int scvIdx) const
Bernd Flemisch's avatar
Bernd Flemisch committed
133
    {
134
        const GlobalPosition& globalPos = fvGeometry.subContVol[scvIdx].global;
135

Bernd Flemisch's avatar
Bernd Flemisch committed
136
137
138
139
        if (isInLens_(globalPos))
            return lensK_;
        return outerK_;
    }
140

141
    /*!
142
     * \brief Returns the porosity \f$[-]\f$
143
     *
144
145
146
     * \param element The finite element
     * \param fvGeometry The finite volume geometry of the element
     * \param scvIdx The local index of the sub-control volume
147
     */
148
    Scalar porosity(const Element &element,
149
                    const FVElementGeometry &fvGeometry,
150
                    int scvIdx) const
Bernd Flemisch's avatar
Bernd Flemisch committed
151
152
    { return 0.4; }

153
    /*!
154
     * \brief Returns the parameter object for the Brooks-Corey material law
155
     *
156
157
158
     * \param element The finite element
     * \param fvGeometry The finite volume geometry of the element
     * \param scvIdx The local index of the sub-control volume
159
     */
160
    const MaterialLawParams& materialLawParams(const Element &element,
161
                                                const FVElementGeometry &fvGeometry,
162
                                                int scvIdx) const
Bernd Flemisch's avatar
Bernd Flemisch committed
163
    {
164
        const GlobalPosition& globalPos = fvGeometry.subContVol[scvIdx].global;
165

Bernd Flemisch's avatar
Bernd Flemisch committed
166
167
168
169
170
171
172
        if (isInLens_(globalPos))
            return lensMaterialParams_;
        return outerMaterialParams_;
    }


private:
173
    bool isInLens_(const GlobalPosition &globalPos) const
Bernd Flemisch's avatar
Bernd Flemisch committed
174
    {
175
        for (int i = 0; i < dimWorld; ++i) {
176
            if (globalPos[i] < lensLowerLeft_[i] || globalPos[i] > lensUpperRight_[i])
Bernd Flemisch's avatar
Bernd Flemisch committed
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
                return false;
        }
        return true;
    }

    GlobalPosition lensLowerLeft_;
    GlobalPosition lensUpperRight_;

    Scalar lensK_;
    Scalar outerK_;
    MaterialLawParams lensMaterialParams_;
    MaterialLawParams outerMaterialParams_;
};

} // end namespace
#endif