12#ifndef DUMUX_LINEAR_HELMHOLTZ_OPERATOR_HH
13#define DUMUX_LINEAR_HELMHOLTZ_OPERATOR_HH
26template <
class Traits>
29 using Scalar =
typename Traits::PrimaryVariables::value_type;
33 template<
class ElementSolution,
class Problem,
class Element,
class SubControlVolume>
34 void update(
const ElementSolution& elemSol,
const Problem&,
const Element&,
const SubControlVolume& scv)
35 { priVars_ = elemSol[scv.indexInElement()]; }
37 Scalar
priVar(
const int pvIdx)
const {
return priVars_[pvIdx]; }
42 PrimaryVariables priVars_;
45template<
class TypeTag>
57 using FVElementGeometry =
typename GridGeometry::LocalView;
58 using SubControlVolume =
typename FVElementGeometry::SubControlVolume;
59 using SubControlVolumeFace =
typename FVElementGeometry::SubControlVolumeFace;
60 using GridView =
typename GridGeometry::GridView;
61 using Element =
typename GridView::template Codim<0>::Entity;
62 static constexpr int dimWorld = GridView::dimensionworld;
64 using ParentType::ParentType;
67 const SubControlVolume&,
68 const VolumeVariables&)
const
69 {
return NumEqVector(0.0); }
72 const Element&,
const FVElementGeometry& fvGeometry,
73 const ElementVolumeVariables& elemVolVars,
74 const SubControlVolumeFace& scvf,
75 const ElementFluxVariablesCache& elemFluxVarsCache)
const
77 NumEqVector flux(0.0);
84 const auto& fluxVarCache = elemFluxVarsCache[scvf];
85 Dune::FieldVector<Scalar, dimWorld> gradConcentration(0.0);
86 for (
const auto& scv : scvs(fvGeometry))
88 const auto& volVars = elemVolVars[scv];
89 gradConcentration.axpy(volVars.concentration(), fluxVarCache.gradN(scv.indexInElement()));
93 flux[0] = -1.0*
vtmv(scvf.unitOuterNormal(), problem.a(), gradConcentration)*scvf.area();
99 const auto& insideScv = fvGeometry.scv(scvf.insideScvIdx());
100 const auto& insideV = elemVolVars[insideScv];
101 const auto& outsideV = elemVolVars[scvf.outsideScvIdx()];
103 fvGeometry, scvf, insideScv, problem.a(), insideV.extrusionFactor()
106 const auto tij = [&]{
108 return scvf.area()*ti;
111 const auto& outsideScv = fvGeometry.scv(scvf.outsideScvIdx());
113 fvGeometry, scvf, outsideScv, problem.a(), outsideV.extrusionFactor()
116 return scvf.area()*(ti * tj)/(ti + tj);
120 const auto u = insideV.priVar(0);
121 const auto v = outsideV.priVar(0);
122 flux[0] = tij*(
u - v);
125 DUNE_THROW(Dune::NotImplemented,
"Discretization method " << GridGeometry::discMethod);
131 const Element&,
const FVElementGeometry&,
132 const ElementVolumeVariables& elemVolVars,
133 const SubControlVolume& scv)
const
135 NumEqVector source(0.0);
136 source[0] = -problem.b()*elemVolVars[scv].priVar(0);
141template<
class TypeTag>
154 template<
class GlobalPosition>
156 { BoundaryTypes values; values.setAllNeumann();
return values; }
158 Scalar
a()
const {
return a_; }
159 Scalar
b()
const {
return b_; }
164template<
class G,
class D,
class S>
182template<
class TypeTag,
class P>
184{
using type =
typename P::Scalar; };
187template<
class TypeTag,
class P>
190 using type = Dune::FieldVector<
196template<
class TypeTag,
class P>
203template<
class TypeTag,
class P>
207template<
class TypeTag,
class P>
211template<
class TypeTag,
class P>
218template<
class TypeTag,
class P>
220{
static constexpr bool value =
true; };
222template<
class TypeTag,
class P>
224{
static constexpr bool value =
true; };
226template<
class TypeTag,
class P>
228{
static constexpr bool value =
true; };
230template<
class TypeTag,
class P>
232{
using type =
typename P::Grid; };
251template<
class Discretization,
class Gr
idView,
class Scalar>
258 auto gridGeometry = std::make_shared<GridGeometry>(gridView);
261 auto problem = std::make_shared<Problem>(gridGeometry, a, b);
265 auto gridVariables = std::make_shared<GridVariables>(problem, gridGeometry);
266 SolutionVector x(gridGeometry->numDofs());
267 gridVariables->init(x);
270 auto assembler = std::make_shared<Assembler>(problem, gridGeometry, gridVariables);
271 using A =
typename Assembler::JacobianMatrix;
272 using V =
typename Assembler::ResidualType;
273 auto jacobian = std::make_shared<A>();
274 auto residual = std::make_shared<V>();
275 assembler->setLinearSystem(jacobian, residual);
276 assembler->assembleJacobian(x);
283template<
class Discretization,
class Gr
idView,
class Scalar>
287template<
class LinearOperator,
class Discretization,
class Gr
idView,
class Scalar>
291template<
class LinearOperator,
class Discretization,
class Gr
idView,
class Comm,
class Scalar>
295template<
class LinearOperator,
class Discretization,
class Gr
idView,
class Scalar>
299template<
class LinearOperator,
class Discretization,
class Gr
idView,
class Comm,
class Scalar>
A linear system assembler (residual and Jacobian) for finite volume schemes.
Class to specify the type of a boundary.
Definition common/boundarytypes.hh:26
Definition helmholtzoperator.hh:143
Scalar b() const
Definition helmholtzoperator.hh:159
BoundaryTypes boundaryTypesAtPos(const GlobalPosition &globalPos) const
Definition helmholtzoperator.hh:155
Scalar a() const
Definition helmholtzoperator.hh:158
HelmholtzModelHomogeneousNeumannProblem(std::shared_ptr< const GridGeometry > gridGeometry, Scalar a, Scalar b)
Definition helmholtzoperator.hh:149
Definition helmholtzoperator.hh:48
NumEqVector computeStorage(const Problem &, const SubControlVolume &, const VolumeVariables &) const
Definition helmholtzoperator.hh:66
NumEqVector computeSource(const Problem &problem, const Element &, const FVElementGeometry &, const ElementVolumeVariables &elemVolVars, const SubControlVolume &scv) const
Definition helmholtzoperator.hh:130
NumEqVector computeFlux(const Problem &problem, const Element &, const FVElementGeometry &fvGeometry, const ElementVolumeVariables &elemVolVars, const SubControlVolumeFace &scvf, const ElementFluxVariablesCache &elemFluxVarsCache) const
Definition helmholtzoperator.hh:71
Definition helmholtzoperator.hh:28
typename Traits::PrimaryVariables PrimaryVariables
Definition helmholtzoperator.hh:31
void update(const ElementSolution &elemSol, const Problem &, const Element &, const SubControlVolume &scv)
Definition helmholtzoperator.hh:34
const PrimaryVariables & priVars() const
Definition helmholtzoperator.hh:38
Scalar priVar(const int pvIdx) const
Definition helmholtzoperator.hh:37
Scalar extrusionFactor() const
Definition helmholtzoperator.hh:39
A linear system assembler (residual and Jacobian) for finite volume schemes (box, tpfa,...
Definition assembly/fvassembler.hh:95
FVProblem(std::shared_ptr< const GridGeometry > gridGeometry, const std::string ¶mGroup="")
Constructor.
Definition common/fvproblem.hh:84
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.
The default local operator than can be specialized for each discretization scheme.
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:24
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
Definition gridcapabilities.hh:57
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:296
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:252
auto makeLaplaceMatrix(const GridView &gridView, const Scalar a=1.0)
make a Laplace matrix operator aΔ
Definition helmholtzoperator.hh:284
const Scalar PengRobinsonMixture< Scalar, StaticParameters >::u
Definition pengrobinsonmixture.hh:138
typename Detail::DiscretizationDefaultLocalOperator< TypeTag >::type DiscretizationDefaultLocalOperator
Definition defaultlocaloperator.hh:26
std::shared_ptr< LinearOperator > makeHelmholtzLinearOperator(const GridView &gridView, const Scalar a, const Scalar b)
Definition helmholtzoperator.hh:288
A helper to deduce a vector with the same size as numbers of equations.
Definition helmholtzoperator.hh:166
G Grid
Definition helmholtzoperator.hh:168
D Discretization
Definition helmholtzoperator.hh:169
S Scalar
Definition helmholtzoperator.hh:167
Definition helmholtzoperator.hh:173
std::tuple< typename P::Discretization > InheritsFrom
Definition helmholtzoperator.hh:174
static constexpr bool value
Definition helmholtzoperator.hh:224
static constexpr bool value
Definition helmholtzoperator.hh:228
static constexpr bool value
Definition helmholtzoperator.hh:220
typename P::Grid type
Definition helmholtzoperator.hh:232
Dumux::Detail::HelmholtzOperator::HelmholtzModelLocalResidual< TypeTag > type
Definition helmholtzoperator.hh:205
Traits type
Definition helmholtzoperator.hh:200
Definition helmholtzoperator.hh:199
static constexpr int numEq()
Definition helmholtzoperator.hh:199
Dune::FieldVector< GetPropType< TypeTag, Properties::Scalar >, GetPropType< TypeTag, Properties::ModelTraits >::numEq() > type
Definition helmholtzoperator.hh:190
Dumux::Detail::HelmholtzOperator::HelmholtzModelHomogeneousNeumannProblem< TypeTag > type
Definition helmholtzoperator.hh:209
typename P::Scalar type
Definition helmholtzoperator.hh:184
Dumux::Detail::HelmholtzOperator::HelmholtzModelVolumeVariables< Traits > type
Definition helmholtzoperator.hh:215
Definition helmholtzoperator.hh:214
GetPropType< TypeTag, Properties::PrimaryVariables > PrimaryVariables
Definition helmholtzoperator.hh:214
Free functions to evaluate the transmissibilities associated with flux evaluations across sub-control...