subcontrolvolumeface.hh 9.33 KB
Newer Older
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
// -*- 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
Dennis Gläser's avatar
Dennis Gläser committed
21
22
 * \ingroup CCMpfaDiscretization
 * \brief The sub control volume face
23
24
25
26
 */
#ifndef DUMUX_DISCRETIZATION_CC_MPFA_SUBCONTROLVOLUMEFACE_HH
#define DUMUX_DISCRETIZATION_CC_MPFA_SUBCONTROLVOLUMEFACE_HH

Timo Koch's avatar
Timo Koch committed
27
#include <vector>
28
29
30
#include <array>
#include <dune/common/reservedvector.hh>
#include <dune/common/fvector.hh>
Timo Koch's avatar
Timo Koch committed
31
#include <dune/geometry/type.hh>
32
#include <dune/geometry/multilineargeometry.hh>
Dennis Gläser's avatar
Dennis Gläser committed
33

34
namespace Dumux {
35

36
37
38
39
40
41
/*!
 * \ingroup CCMpfaDiscretization
 * \brief Default traits class to be used for the sub-control volume faces
 *        for the cell-centered finite volume scheme using MPFA
 * \tparam GV the type of the grid view
 */
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
template<class GridView>
struct CCMpfaDefaultScvfGeometryTraits
{
    using Grid = typename GridView::Grid;

    static const int dim = Grid::dimension;
    static const int dimWorld = Grid::dimensionworld;

    using Scalar = typename Grid::ctype;
    using GridIndexType = typename Grid::LeafGridView::IndexSet::IndexType;
    using LocalIndexType = unsigned int;
    using OutsideGridIndexStorage = typename std::conditional_t< (dim<dimWorld),
                                                                 std::vector<GridIndexType>,
                                                                 Dune::ReservedVector<GridIndexType, 1> >;

    // we use geometry traits that use static corner vectors to and a fixed geometry type
    template <class ct>
    struct ScvfMLGTraits : public Dune::MultiLinearGeometryTraits<ct>
    {
        // we use static vectors to store the corners as we know
        // the number of corners in advance (2^(dim-1) corners (1<<(dim-1))
        template< int mydim, int cdim >
        struct CornerStorage
        {
            using Type = std::array< Dune::FieldVector< ct, cdim >, (1<<(dim-1)) >;
        };

        // we know all scvfs will have the same geometry type
        template< int dim >
        struct hasSingleGeometryType
        {
            static const bool v = true;
            static const unsigned int topologyId = Dune::Impl::CubeTopology< dim >::type::id;
        };
    };

    using Geometry = Dune::MultiLinearGeometry<Scalar, dim-1, dimWorld, ScvfMLGTraits<Scalar> >;
    using CornerStorage = typename ScvfMLGTraits<Scalar>::template CornerStorage<dim-1, dimWorld>::Type;
    using GlobalPosition = typename CornerStorage::value_type;
};

83
/*!
Dennis Gläser's avatar
Dennis Gläser committed
84
 * \ingroup CCMpfaDiscretization
85
86
 * \brief Class for a sub control volume face in mpfa methods, i.e a part of the boundary
 *        of a control volume we compute fluxes on.
87
88
 * \tparam GV the type of the grid view
 * \tparam T the scvf geometry traits
89
 */
90
91
template<class GV,
         class T = CCMpfaDefaultScvfGeometryTraits<GV> >
92
class CCMpfaSubControlVolumeFace
93
{
94
95
96
97
98
    using GridIndexType = typename T::GridIndexType;
    using Scalar = typename T::Scalar;
    using CornerStorage = typename T::CornerStorage;
    using OutsideGridIndexStorage = typename T::OutsideGridIndexStorage;
    using Geometry = typename T::Geometry;
99
100

public:
101
    //! export the type used for global coordinates
102
    using GlobalPosition = typename T::GlobalPosition;
103
    //! state the traits public and thus export all types
104
    using Traits = T;
Dennis Gläser's avatar
Dennis Gläser committed
105

106
    /*!
107
     * \brief Constructor
108
     *
Dennis Gläser's avatar
Dennis Gläser committed
109
     * \param helper The helper class for mpfa schemes
110
111
112
     * \param corners The corners of the scv face
     * \param unitOuterNormal The unit outer normal vector of the scvf
     * \param vIdxGlobal The global vertex index the scvf is connected to
113
     * \param vIdxLocal The element-local vertex index the scvf is connected to
114
     * \param scvfIndex The global index of this scv face
115
116
     * \param insideScvIdx The inside scv index connected to this face
     * \param outsideScvIndices The outside scv indices connected to this face
117
118
119
     * \param q The parameterization of the quadrature point on the scvf for flux calculation
     * \param boundary Boolean to specify whether or not the scvf is on a boundary
     */
120
    template<class MpfaHelper>
121
122
123
124
125
126
127
    CCMpfaSubControlVolumeFace(const MpfaHelper& helper,
                               CornerStorage&& corners,
                               GlobalPosition&& unitOuterNormal,
                               GridIndexType vIdxGlobal,
                               unsigned int vIdxLocal,
                               GridIndexType scvfIndex,
                               GridIndexType insideScvIdx,
128
                               const OutsideGridIndexStorage& outsideScvIndices,
129
130
                               Scalar q,
                               bool boundary)
Dennis Gläser's avatar
Dennis Gläser committed
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
    : boundary_(boundary)
    ,   vertexIndex_(vIdxGlobal)
    ,   scvfIndex_(scvfIndex)
    ,   insideScvIdx_(insideScvIdx)
    ,   outsideScvIndices_(outsideScvIndices)
    ,   vIdxInElement_(vIdxLocal)
    ,   corners_(std::move(corners))
    ,   center_(0.0)
    ,   unitOuterNormal_(std::move(unitOuterNormal))
    {
          // compute the center of the scvf
          for (const auto& corner : corners_)
              center_ += corner;
          center_ /= corners_.size();

          // use helper class to obtain area & integration point
          ipGlobal_ = helper.getScvfIntegrationPoint(corners_, q);
          area_ = helper.getScvfArea(corners_);
    }
150
151

    //! The area of the sub control volume face
152
153
    Scalar area() const
    { return area_; }
154
155

    //! returns bolean if the sub control volume face is on the domain boundary
156
157
    bool boundary() const
    { return boundary_; }
158
159

    //! The global index of this sub control volume face
160
161
    GridIndexType index() const
    { return scvfIndex_; }
162
163

    //! Returns the index of the vertex the scvf is connected to
164
165
    GridIndexType vertexIndex() const
    { return vertexIndex_; }
166
167

    //! Returns the element-local vertex index the scvf is connected to
168
169
    unsigned int vertexIndexInElement() const
    { return vIdxInElement_; }
170
171

    //! index of the inside sub control volume
172
173
    GridIndexType insideScvIdx() const
    { return insideScvIdx_; }
174
175

    //! The number of outside scvs connection via this scv face
176
177
    std::size_t numOutsideScvs() const
    { return outsideScvIndices_.size(); }
178
179
180

    //! index of the outside sub control volume or boundary scv index
    //! returns undefined behaviour if index exceeds numOutsideScvs
181
182
    GridIndexType outsideScvIdx(int i = 0) const
    { return outsideScvIndices_[i]; }
183
184

    //! returns the outside scv indices (can be more than one index for dim < dimWorld)
185
186
    const OutsideGridIndexStorage& outsideScvIndices() const
    { return outsideScvIndices_; }
187
188

    //! Returns the number of corners
189
190
    std::size_t corners() const
    { return corners_.size(); }
191
192
193
194
195
196
197
198
199

    //! Returns the corner for a given local index
    const GlobalPosition& corner(unsigned int localIdx) const
    {
        assert(localIdx < corners_.size() && "provided index exceeds the number of corners");
        return corners_[localIdx];
    }

    //! Returns the global position of the vertex the scvf is connected to
200
201
    const GlobalPosition& vertexCorner() const
    { return corners_.back(); }
202
203

    //! Returns the global position of the center of the element facet this scvf is embedded in
204
205
    const GlobalPosition& facetCorner() const
    { return corner(0); }
206
207

    //! The center of the sub control volume face
208
209
    const GlobalPosition& center() const
    { return center_; }
210
211

    //! The integration point for flux evaluations in global coordinates
212
213
    const GlobalPosition& ipGlobal() const
    { return ipGlobal_; }
214
215

    //! returns the unit outer normal vector (assumes non-curved geometries)
216
217
    const GlobalPosition& unitOuterNormal() const
    { return unitOuterNormal_; }
218
219

    //! The geometry of the sub control volume face
220
221
    Geometry geometry() const
    { return Geometry(Dune::GeometryTypes::cube(Geometry::mydimension), corners_); }
222
223
224
225
226
227

private:
    bool boundary_;
    GridIndexType vertexIndex_;
    GridIndexType scvfIndex_;
    GridIndexType insideScvIdx_;
228
    OutsideGridIndexStorage outsideScvIndices_;
229
230
231
232
233
234
235
    unsigned int vIdxInElement_;

    CornerStorage corners_;
    GlobalPosition center_;
    GlobalPosition ipGlobal_;
    GlobalPosition unitOuterNormal_;
    Scalar area_;
236
};
237

238
} // end namespace Dumux
239

240
#endif