12#ifndef DUMUX_LINEAR_HELMHOLTZ_OPERATOR_HH
13#define DUMUX_LINEAR_HELMHOLTZ_OPERATOR_HH
25template <
class Traits>
28 using Scalar =
typename Traits::PrimaryVariables::value_type;
32 template<
class ElementSolution,
class Problem,
class Element,
class SubControlVolume>
33 void update(
const ElementSolution& elemSol,
const Problem&,
const Element&,
const SubControlVolume& scv)
34 { priVars_ = elemSol[scv.indexInElement()]; }
36 Scalar
priVar(
const int pvIdx)
const {
return priVars_[pvIdx]; }
44template<
class TypeTag>
46:
public GetPropType<TypeTag, Properties::BaseLocalResidual>
56 using FVElementGeometry =
typename GridGeometry::LocalView;
57 using SubControlVolume =
typename FVElementGeometry::SubControlVolume;
58 using SubControlVolumeFace =
typename FVElementGeometry::SubControlVolumeFace;
59 using GridView =
typename GridGeometry::GridView;
60 using Element =
typename GridView::template Codim<0>::Entity;
61 static constexpr int dimWorld = GridView::dimensionworld;
63 using ParentType::ParentType;
66 const SubControlVolume&,
67 const VolumeVariables&)
const
68 {
return NumEqVector(0.0); }
71 const Element&,
const FVElementGeometry& fvGeometry,
72 const ElementVolumeVariables& elemVolVars,
73 const SubControlVolumeFace& scvf,
74 const ElementFluxVariablesCache& elemFluxVarsCache)
const
76 NumEqVector flux(0.0);
83 const auto& fluxVarCache = elemFluxVarsCache[scvf];
84 Dune::FieldVector<Scalar, dimWorld> gradConcentration(0.0);
85 for (
const auto& scv : scvs(fvGeometry))
87 const auto& volVars = elemVolVars[scv];
88 gradConcentration.axpy(volVars.concentration(), fluxVarCache.gradN(scv.indexInElement()));
92 flux[0] = -1.0*
vtmv(scvf.unitOuterNormal(), problem.a(), gradConcentration)*scvf.area();
98 const auto& insideScv = fvGeometry.scv(scvf.insideScvIdx());
99 const auto& insideV = elemVolVars[insideScv];
100 const auto& outsideV = elemVolVars[scvf.outsideScvIdx()];
102 fvGeometry, scvf, insideScv, problem.a(), insideV.extrusionFactor()
105 const auto tij = [&]{
107 return scvf.area()*ti;
110 const auto& outsideScv = fvGeometry.scv(scvf.outsideScvIdx());
112 fvGeometry, scvf, outsideScv, problem.a(), outsideV.extrusionFactor()
115 return scvf.area()*(ti * tj)/(ti + tj);
119 const auto u = insideV.priVar(0);
120 const auto v = outsideV.priVar(0);
121 flux[0] = tij*(u - v);
124 DUNE_THROW(Dune::NotImplemented,
"Discretization method " << GridGeometry::discMethod);
130 const Element&,
const FVElementGeometry&,
131 const ElementVolumeVariables& elemVolVars,
132 const SubControlVolume& scv)
const
134 NumEqVector source(0.0);
135 source[0] = -problem.b()*elemVolVars[scv].priVar(0);
140template<
class TypeTag>
153 template<
class GlobalPosition>
157 Scalar
a()
const {
return a_; }
158 Scalar
b()
const {
return b_; }
163template<
class G,
class D,
class S>
181template<
class TypeTag,
class P>
183{
using type =
typename P::Scalar; };
186template<
class TypeTag,
class P>
189 using type = Dune::FieldVector<
195template<
class TypeTag,
class P>
198 struct Traits {
static constexpr int numEq() {
return 1; } };
202template<
class TypeTag,
class P>
206template<
class TypeTag,
class P>
210template<
class TypeTag,
class P>
217template<
class TypeTag,
class P>
219{
static constexpr bool value =
true; };
221template<
class TypeTag,
class P>
223{
static constexpr bool value =
true; };
225template<
class TypeTag,
class P>
227{
static constexpr bool value =
true; };
229template<
class TypeTag,
class P>
231{
using type =
typename P::Grid; };
250template<
class Discretization,
class Gr
idView,
class Scalar>
257 auto gridGeometry = std::make_shared<GridGeometry>(gridView);
260 auto problem = std::make_shared<Problem>(gridGeometry, a, b);
264 auto gridVariables = std::make_shared<GridVariables>(problem, gridGeometry);
265 SolutionVector x(gridGeometry->numDofs());
266 gridVariables->init(x);
269 auto assembler = std::make_shared<Assembler>(problem, gridGeometry, gridVariables);
270 using A =
typename Assembler::JacobianMatrix;
271 using V =
typename Assembler::ResidualType;
272 auto jacobian = std::make_shared<A>();
273 auto residual = std::make_shared<V>();
274 assembler->setLinearSystem(jacobian, residual);
275 assembler->assembleJacobian(x);
282template<
class Discretization,
class Gr
idView,
class Scalar>
284{
return makeHelmholtzMatrix<Discretization>(gridView, a, 0.0); };
286template<
class LinearOperator,
class Discretization,
class Gr
idView,
class Scalar>
288{
return std::make_shared<LinearOperator>(makeHelmholtzMatrix<Discretization>(gridView, a, b)); }
290template<
class LinearOperator,
class Discretization,
class Gr
idView,
class Comm,
class Scalar>
292{
return std::make_shared<LinearOperator>(makeHelmholtzMatrix<Discretization>(gridView, a, b), comm); }
294template<
class LinearOperator,
class Discretization,
class Gr
idView,
class Scalar>
296{
return std::make_shared<LinearOperator>(makeLaplaceMatrix<Discretization>(gridView, a)); }
298template<
class LinearOperator,
class Discretization,
class Gr
idView,
class Comm,
class Scalar>
300{
return std::make_shared<LinearOperator>(makeLaplaceMatrix<Discretization>(gridView, a), comm); }
A linear system assembler (residual and Jacobian) for finite volume schemes.
Class to specify the type of a boundary.
Definition: common/boundarytypes.hh:26
void setAllNeumann()
Set all boundary conditions to Neumann.
Definition: common/boundarytypes.hh:90
Definition: helmholtzoperator.hh:142
Scalar b() const
Definition: helmholtzoperator.hh:158
BoundaryTypes boundaryTypesAtPos(const GlobalPosition &globalPos) const
Definition: helmholtzoperator.hh:154
Scalar a() const
Definition: helmholtzoperator.hh:157
HelmholtzModelHomogeneousNeumannProblem(std::shared_ptr< const GridGeometry > gridGeometry, Scalar a, Scalar b)
Definition: helmholtzoperator.hh:148
Definition: helmholtzoperator.hh:47
NumEqVector computeStorage(const Problem &, const SubControlVolume &, const VolumeVariables &) const
Definition: helmholtzoperator.hh:65
NumEqVector computeSource(const Problem &problem, const Element &, const FVElementGeometry &, const ElementVolumeVariables &elemVolVars, const SubControlVolume &scv) const
Definition: helmholtzoperator.hh:129
NumEqVector computeFlux(const Problem &problem, const Element &, const FVElementGeometry &fvGeometry, const ElementVolumeVariables &elemVolVars, const SubControlVolumeFace &scvf, const ElementFluxVariablesCache &elemFluxVarsCache) const
Definition: helmholtzoperator.hh:70
Definition: helmholtzoperator.hh:27
typename Traits::PrimaryVariables PrimaryVariables
Definition: helmholtzoperator.hh:30
void update(const ElementSolution &elemSol, const Problem &, const Element &, const SubControlVolume &scv)
Definition: helmholtzoperator.hh:33
const PrimaryVariables & priVars() const
Definition: helmholtzoperator.hh:37
Scalar priVar(const int pvIdx) const
Definition: helmholtzoperator.hh:36
Scalar extrusionFactor() const
Definition: helmholtzoperator.hh:38
A linear system assembler (residual and Jacobian) for finite volume schemes (box, tpfa,...
Definition: assembly/fvassembler.hh:94
Base class for all finite-volume problems.
Definition: common/fvproblem.hh:43
const GridGeometry & gridGeometry() const
The finite volume grid geometry.
Definition: common/fvproblem.hh:520
Class to specify the type of a boundary.
Base class for all finite volume problems.
Defines all properties used in Dumux.
Tensor::field_type computeTpfaTransmissibility(const FVElementGeometry &fvGeometry, const typename FVElementGeometry::SubControlVolumeFace &scvf, const typename FVElementGeometry::SubControlVolume &scv, const Tensor &T, typename FVElementGeometry::SubControlVolume::Traits::Scalar extrusionFactor)
Free function to evaluate the Tpfa transmissibility associated with the flux (in the form of flux = T...
Definition: tpfa/computetransmissibility.hh:36
Dune::DenseMatrix< MAT >::value_type vtmv(const Dune::DenseVector< V1 > &v1, const Dune::DenseMatrix< MAT > &M, const Dune::DenseVector< V2 > &v2)
Evaluates the scalar product of a vector v2, projected by a matrix M, with a vector v1.
Definition: math.hh:880
typename NumEqVectorTraits< PrimaryVariables >::type NumEqVector
A vector with the same size as numbers of equations This is the default implementation and has to be ...
Definition: numeqvector.hh:34
typename GetProp< TypeTag, Property >::type GetPropType
get the type alias defined in the property
Definition: propertysystem.hh:296
Define some often used mathematical functions.
Definition: helmholtzoperator.hh:23
constexpr FCDiamond fcdiamond
Definition: method.hh:152
constexpr CCTpfa cctpfa
Definition: method.hh:145
constexpr Box box
Definition: method.hh:147
constexpr PQ1Bubble pq1bubble
Definition: method.hh:148
The energy balance equation for a porous solid.
Definition: common/properties.hh:26
std::shared_ptr< LinearOperator > makeLaplaceLinearOperator(const GridView &gridView, const Scalar a)
Definition: helmholtzoperator.hh:295
auto makeHelmholtzMatrix(const GridView &gridView, const Scalar a=1.0, const Scalar b=1.0)
make a Helmholtz matrix operator (aΔ + bI)
Definition: helmholtzoperator.hh:251
auto makeLaplaceMatrix(const GridView &gridView, const Scalar a=1.0)
make a Laplace matrix operator aΔ
Definition: helmholtzoperator.hh:283
std::shared_ptr< LinearOperator > makeHelmholtzLinearOperator(const GridView &gridView, const Scalar a, const Scalar b)
Definition: helmholtzoperator.hh:287
A helper to deduce a vector with the same size as numbers of equations.
Definition: helmholtzoperator.hh:165
G Grid
Definition: helmholtzoperator.hh:167
D Discretization
Definition: helmholtzoperator.hh:168
S Scalar
Definition: helmholtzoperator.hh:166
Definition: helmholtzoperator.hh:172
std::tuple< typename P::Discretization > InheritsFrom
Definition: helmholtzoperator.hh:173
typename P::Grid type
Definition: helmholtzoperator.hh:231
Traits type
Definition: helmholtzoperator.hh:199
static constexpr int numEq()
Definition: helmholtzoperator.hh:198
Dune::FieldVector< GetPropType< TypeTag, Properties::Scalar >, GetPropType< TypeTag, Properties::ModelTraits >::numEq() > type
Definition: helmholtzoperator.hh:192
typename P::Scalar type
Definition: helmholtzoperator.hh:183
GetPropType< TypeTag, Properties::PrimaryVariables > PrimaryVariables
Definition: helmholtzoperator.hh:213
Free functions to evaluate the transmissibilities associated with flux evaluations across sub-control...