3.1-git
DUNE for Multi-{Phase, Component, Scale, Physics, ...} flow and transport in porous media
discretization/box/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 *****************************************************************************/
24#ifndef DUMUX_DISCRETIZATION_BOX_FLUXVARIABLES_CACHE_HH
25#define DUMUX_DISCRETIZATION_BOX_FLUXVARIABLES_CACHE_HH
26
27#include <dune/common/fvector.hh>
28#include <dune/localfunctions/lagrange/pqkfactory.hh>
29
30namespace Dumux {
31
38template< class Scalar, class GridGeometry >
40{
41 using GridView = typename GridGeometry::GridView;
42 using Element = typename GridView::template Codim<0>::Entity;
43 using GlobalPosition = typename Element::Geometry::GlobalCoordinate;
44
45 using FVElementGeometry = typename GridGeometry::LocalView;
46 using SubControlVolumeFace = typename GridGeometry::SubControlVolumeFace;
47
48 static const int dim = GridView::dimension;
49 static const int dimWorld = GridView::dimensionworld;
50
51 using CoordScalar = typename GridView::ctype;
52 using FeCache = Dune::PQkLocalFiniteElementCache<CoordScalar, Scalar, dim, 1>;
53 using FeLocalBasis = typename FeCache::FiniteElementType::Traits::LocalBasisType;
54 using ShapeJacobian = typename FeLocalBasis::Traits::JacobianType;
55 using ShapeValue = typename Dune::FieldVector<Scalar, 1>;
56 using JacobianInverseTransposed = typename Element::Geometry::JacobianInverseTransposed;
57
58public:
59
61 template< class Problem, class ElementVolumeVariables >
62 void update(const Problem& problem,
63 const Element& element,
64 const FVElementGeometry& fvGeometry,
65 const ElementVolumeVariables& elemVolVars,
66 const SubControlVolumeFace& scvf)
67 {
68 update(problem, element, fvGeometry, elemVolVars, scvf.ipGlobal());
69 }
70
72 template< class Problem, class ElementVolumeVariables >
73 void update(const Problem& problem,
74 const Element& element,
75 const FVElementGeometry& fvGeometry,
76 const ElementVolumeVariables& elemVolVars,
77 const GlobalPosition& globalPos)
78 {
79 const auto geometry = element.geometry();
80 const auto& localBasis = fvGeometry.feLocalBasis();
81
82 // evaluate shape functions and gradients at the integration point
83 ipGlobal_ = globalPos;
84 const auto ipLocal = geometry.local(ipGlobal_);
85 jacInvT_ = geometry.jacobianInverseTransposed(ipLocal);
86 localBasis.evaluateJacobian(ipLocal, shapeJacobian_);
87 localBasis.evaluateFunction(ipLocal, shapeValues_); // shape values for rho
88
89 // compute the gradN at for every scv/dof
90 gradN_.resize(fvGeometry.numScv());
91 for (const auto& scv: scvs(fvGeometry))
92 jacInvT_.mv(shapeJacobian_[scv.localDofIndex()][0], gradN_[scv.indexInElement()]);
93 }
94
96 const GlobalPosition& ipGlobal() const { return ipGlobal_; }
98 const std::vector<ShapeJacobian>& shapeJacobian() const { return shapeJacobian_; }
100 const std::vector<ShapeValue>& shapeValues() const { return shapeValues_; }
102 const JacobianInverseTransposed& jacInvT() const { return jacInvT_; }
104 const GlobalPosition& gradN(unsigned int scvIdxInElement) const { return gradN_[scvIdxInElement]; }
105
106private:
107 GlobalPosition ipGlobal_;
108 std::vector<GlobalPosition> gradN_;
109 std::vector<ShapeJacobian> shapeJacobian_;
110 std::vector<ShapeValue> shapeValues_;
111 JacobianInverseTransposed jacInvT_;
112};
113
114} // end namespace Dumux
115
116#endif
make the local view function available whenever we use the grid geometry
Definition: adapt.hh:29
Flux variables cache class for the box scheme. For the box scheme, this class does not contain any ph...
Definition: discretization/box/fluxvariablescache.hh:40
const std::vector< ShapeJacobian > & shapeJacobian() const
returns the shape function gradients in local coordinates at the integration point
Definition: discretization/box/fluxvariablescache.hh:98
const GlobalPosition & gradN(unsigned int scvIdxInElement) const
returns the shape function gradients in global coordinates at the integration point
Definition: discretization/box/fluxvariablescache.hh:104
const std::vector< ShapeValue > & shapeValues() const
returns the shape function values at the integration point
Definition: discretization/box/fluxvariablescache.hh:100
void update(const Problem &problem, const Element &element, const FVElementGeometry &fvGeometry, const ElementVolumeVariables &elemVolVars, const GlobalPosition &globalPos)
update the cache for a given global position
Definition: discretization/box/fluxvariablescache.hh:73
const GlobalPosition & ipGlobal() const
returns the global position for which this cache has been updated
Definition: discretization/box/fluxvariablescache.hh:96
void update(const Problem &problem, const Element &element, const FVElementGeometry &fvGeometry, const ElementVolumeVariables &elemVolVars, const SubControlVolumeFace &scvf)
update the cache for an scvf
Definition: discretization/box/fluxvariablescache.hh:62
const JacobianInverseTransposed & jacInvT() const
returns inverse transposed jacobian at the integration point
Definition: discretization/box/fluxvariablescache.hh:102