3.5-git
DUNE for Multi-{Phase, Component, Scale, Physics, ...} flow and transport in porous media
fluxstencil.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_DISCRETIZATION_FLUXSTENCIL_HH
25#define DUMUX_DISCRETIZATION_FLUXSTENCIL_HH
26
27#include <vector>
28
29#include <dune/common/reservedvector.hh>
32
33namespace Dumux {
34
44template<class FVElementGeometry, class DiscretizationMethod = typename FVElementGeometry::GridGeometry::DiscretizationMethod>
46
47/*
48 * \ingroup Discretization
49 * \brief Flux stencil specialization for the cell-centered tpfa scheme
50 * \tparam FVElementGeometry The local view on the finite volume grid geometry
51 */
52template<class FVElementGeometry>
53class FluxStencil<FVElementGeometry, DiscretizationMethods::CCTpfa>
54{
55 using GridGeometry = typename FVElementGeometry::GridGeometry;
56 using SubControlVolumeFace = typename GridGeometry::SubControlVolumeFace;
57 using GridView = typename GridGeometry::GridView;
58 using Element = typename GridView::template Codim<0>::Entity;
59 using GridIndexType = typename IndexTraits<GridView>::GridIndex;
60
61public:
63 using ScvfStencilIForJ = Dune::ReservedVector<GridIndexType, 1>;
64
66 using Stencil = typename SubControlVolumeFace::Traits::GridIndexStorage;
67
69 static Stencil stencil(const Element& element,
70 const FVElementGeometry& fvGeometry,
71 const SubControlVolumeFace& scvf)
72 {
73 if (scvf.boundary())
74 return Stencil({scvf.insideScvIdx()});
75 else if (scvf.numOutsideScvs() > 1)
76 {
77 Stencil stencil({scvf.insideScvIdx()});
78 for (unsigned int i = 0; i < scvf.numOutsideScvs(); ++i)
79 stencil.push_back(scvf.outsideScvIdx(i));
80 return stencil;
81 }
82 else
83 return Stencil({scvf.insideScvIdx(), scvf.outsideScvIdx()});
84 }
85};
86
87/*
88 * \ingroup Discretization
89 * \brief Flux stencil specialization for the cell-centered mpfa scheme
90 * \tparam FVElementGeometry The local view on the finite volume grid geometry
91 */
92template<class FVElementGeometry>
93class FluxStencil<FVElementGeometry, DiscretizationMethods::CCMpfa>
94{
95 using GridGeometry = typename FVElementGeometry::GridGeometry;
96 using SubControlVolumeFace = typename GridGeometry::SubControlVolumeFace;
97 using GridView = typename GridGeometry::GridView;
98 using Element = typename GridView::template Codim<0>::Entity;
99 using GridIndexType = typename IndexTraits<GridView>::GridIndex;
100
101 // Use the stencil type of the primary interaction volume
102 using NodalIndexSet = typename GridGeometry::GridIVIndexSets::DualGridIndexSet::NodalIndexSet;
103
104public:
106 using ScvfStencilIForJ = std::vector<GridIndexType>;
107
109 using Stencil = typename NodalIndexSet::NodalGridStencilType;
110
112 static const Stencil& stencil(const Element& element,
113 const FVElementGeometry& fvGeometry,
114 const SubControlVolumeFace& scvf)
115 {
116 const auto& gridGeometry = fvGeometry.gridGeometry();
117
118 // return the scv (element) indices in the interaction region
119 if (gridGeometry.vertexUsesSecondaryInteractionVolume(scvf.vertexIndex()))
120 return gridGeometry.gridInteractionVolumeIndexSets().secondaryIndexSet(scvf).gridScvIndices();
121 else
122 return gridGeometry.gridInteractionVolumeIndexSets().primaryIndexSet(scvf).gridScvIndices();
123 }
124};
125
126} // end namespace Dumux
127
128#endif
Defines the index types used for grid and local indices.
The available discretization methods in Dumux.
Definition: adapt.hh:29
typename GridView::IndexSet::IndexType GridIndex
Definition: indextraits.hh:39
The flux stencil specialized for different discretization schemes.
Definition: fluxstencil.hh:45
Dune::ReservedVector< GridIndexType, 1 > ScvfStencilIForJ
Each cell I couples to a cell J always only via one face.
Definition: fluxstencil.hh:63
typename SubControlVolumeFace::Traits::GridIndexStorage Stencil
The flux stencil type.
Definition: fluxstencil.hh:66
static Stencil stencil(const Element &element, const FVElementGeometry &fvGeometry, const SubControlVolumeFace &scvf)
Returns the flux stencil.
Definition: fluxstencil.hh:69
std::vector< GridIndexType > ScvfStencilIForJ
We don't know yet how many faces couple to a neighboring element.
Definition: fluxstencil.hh:106
static const Stencil & stencil(const Element &element, const FVElementGeometry &fvGeometry, const SubControlVolumeFace &scvf)
Returns the indices of the elements required for flux calculation on an scvf.
Definition: fluxstencil.hh:112
typename NodalIndexSet::NodalGridStencilType Stencil
The flux stencil type.
Definition: fluxstencil.hh:109