3.1-git
DUNE for Multi-{Phase, Component, Scale, Physics, ...} flow and transport in porous media
porousmediumflow/boxdfm/fluxvariablescache.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 *****************************************************************************/
26#ifndef DUMUX_POROUSMEDIUM_BOXDFM_FLUXVARIABLESCACHE_HH
27#define DUMUX_POROUSMEDIUM_BOXDFM_FLUXVARIABLESCACHE_HH
28
29#include <dune/localfunctions/lagrange/pqkfactory.hh>
30
34
35namespace Dumux {
36
43template<class TypeTag>
45{
50 using FVElementGeometry = typename GetPropType<TypeTag, Properties::GridGeometry>::LocalView;
51 using ElementVolumeVariables = typename GetPropType<TypeTag, Properties::GridVolumeVariables>::LocalView;
52 using SubControlVolumeFace = typename FVElementGeometry::SubControlVolumeFace;
53 using Element = typename GridView::template Codim<0>::Entity;
54 using GridIndexType = typename GridView::IndexSet::IndexType;
55 using Stencil = std::vector<GridIndexType>;
56 using TransmissibilityVector = std::vector<GridIndexType>;
57
58 using CoordScalar = typename GridView::ctype;
59 static const int dim = GridView::dimension;
60 static const int dimWorld = GridView::dimensionworld;
61
62 using FeCache = Dune::PQkLocalFiniteElementCache<CoordScalar, Scalar, dim, 1>;
63 using FeLocalBasis = typename FeCache::FiniteElementType::Traits::LocalBasisType;
64 using ShapeJacobian = typename FeLocalBasis::Traits::JacobianType;
65 using ShapeValue = typename Dune::FieldVector<Scalar, 1>;
66 using JacobianInverseTransposed = typename Element::Geometry::JacobianInverseTransposed;
67 using GlobalPosition = Dune::FieldVector<Scalar, dimWorld>;
68
69public:
70
71 void update(const Problem& problem,
72 const Element& element,
73 const FVElementGeometry& fvGeometry,
74 const ElementVolumeVariables& elemVolVars,
75 const SubControlVolumeFace& scvf)
76 {
77 const auto geometry = element.geometry();
78 const auto& localBasis = fvGeometry.feLocalBasis();
79
80 // evaluate shape functions and gradients at the integration point
81 std::vector<ShapeValue> shapeVals;
82 const auto ipLocal = geometry.local(scvf.ipGlobal());
83 jacInvT_ = geometry.jacobianInverseTransposed(ipLocal);
84 localBasis.evaluateJacobian(ipLocal, shapeJacobian_);
85 localBasis.evaluateFunction(ipLocal, shapeVals);
86
87 // set the shape values
88 shapeValues_.resize(fvGeometry.numScv(), 0.0);
89 if (!scvf.isOnFracture())
90 std::copy(shapeVals.begin(), shapeVals.end(), shapeValues_.begin());
91 else
92 {
93 const auto thisFacetIdx = scvf.facetIndexInElement();
94 for (const auto& scv: scvs(fvGeometry))
95 if (scv.isOnFracture() && scv.facetIndexInElement() == thisFacetIdx)
96 shapeValues_[scv.indexInElement()] = shapeVals[scv.localDofIndex()];
97 }
98
99 // set the shape value gradients
100 gradN_.resize(fvGeometry.numScv(), GlobalPosition(0.0));
101 if (!scvf.isOnFracture())
102 {
103 for (const auto& scv: scvs(fvGeometry))
104 if (!scv.isOnFracture())
105 jacInvT_.mv(shapeJacobian_[scv.localDofIndex()][0], gradN_[scv.indexInElement()]);
106 }
107 else
108 {
109 const auto thisFacetIdx = scvf.facetIndexInElement();
110
111 // first, find all local dofs on this facet
112 std::vector<unsigned int> facetLocalDofs;
113 for (const auto& scv : scvs(fvGeometry))
114 if (scv.isOnFracture() && scv.facetIndexInElement() == thisFacetIdx)
115 facetLocalDofs.push_back(scv.localDofIndex());
116
117 for (const auto& scv: scvs(fvGeometry))
118 {
119 // now, create entries for all fracture scvs on this same facet ...
120 if (scv.isOnFracture() && scv.facetIndexInElement() == thisFacetIdx)
121 jacInvT_.mv(shapeJacobian_[scv.localDofIndex()][0], gradN_[scv.indexInElement()]);
122
123 // ... and those non-fracture scvs that are not on this facet
124 else if (!scv.isOnFracture()
125 && std::find( facetLocalDofs.begin(),
126 facetLocalDofs.end(),
127 scv.localDofIndex() ) == facetLocalDofs.end())
128 {
129 jacInvT_.mv(shapeJacobian_[scv.localDofIndex()][0], gradN_[scv.indexInElement()]);
130 }
131 }
132 }
133 }
134
136 const std::vector<ShapeJacobian>& shapeJacobian() const { return shapeJacobian_; }
138 const std::vector<ShapeValue>& shapeValues() const { return shapeValues_; }
140 const JacobianInverseTransposed& jacInvT() const { return jacInvT_; }
142 const GlobalPosition& gradN(unsigned int scvIdxInElement) const { return gradN_[scvIdxInElement]; }
143
144private:
145 std::vector<GlobalPosition> gradN_;
146 std::vector<ShapeJacobian> shapeJacobian_;
147 std::vector<ShapeValue> shapeValues_;
148 JacobianInverseTransposed jacInvT_;
149};
150
151} // end namespace Dumux
152
153#endif
The available discretization methods in Dumux.
Classes related to flux variables caching.
make the local view function available whenever we use the grid geometry
Definition: adapt.hh:29
typename Properties::Detail::GetPropImpl< TypeTag, Property >::type::type GetPropType
get the type alias defined in the property (equivalent to old macro GET_PROP_TYPE(....
Definition: propertysystem.hh:149
We only store discretization-related quantities for the box method. However, we cannot reuse the cach...
Definition: porousmediumflow/boxdfm/fluxvariablescache.hh:45
const GlobalPosition & gradN(unsigned int scvIdxInElement) const
Returns the shape value gradients corresponding to an scv.
Definition: porousmediumflow/boxdfm/fluxvariablescache.hh:142
const std::vector< ShapeValue > & shapeValues() const
Returns the shape values for all scvs at the integration point.
Definition: porousmediumflow/boxdfm/fluxvariablescache.hh:138
const std::vector< ShapeJacobian > & shapeJacobian() const
Returns the Jacobian of the shape functions at the integration point.
Definition: porousmediumflow/boxdfm/fluxvariablescache.hh:136
void update(const Problem &problem, const Element &element, const FVElementGeometry &fvGeometry, const ElementVolumeVariables &elemVolVars, const SubControlVolumeFace &scvf)
Definition: porousmediumflow/boxdfm/fluxvariablescache.hh:71
const JacobianInverseTransposed & jacInvT() const
Returns the shape value gradients for all scvs at the integration point.
Definition: porousmediumflow/boxdfm/fluxvariablescache.hh:140
Declares all properties used in Dumux.