version 3.10.0
Loading...
Searching...
No Matches
solidmechanics/elastic/localresidual.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// SPDX-FileCopyrightText: Copyright © DuMux Project contributors, see AUTHORS.md in root folder
5// SPDX-License-Identifier: GPL-3.0-or-later
6//
13#ifndef DUMUX_SOLIDMECHANICS_ELASTIC_LOCAL_RESIDUAL_HH
14#define DUMUX_SOLIDMECHANICS_ELASTIC_LOCAL_RESIDUAL_HH
15
21
22namespace Dumux {
23
29template<class TypeTag>
32{
35
36 using GridView = typename GridGeometry::GridView;
37 using Element = typename GridView::template Codim<0>::Entity;
38
41 using FVElementGeometry = typename GridGeometry::LocalView;
42 using SubControlVolume = typename FVElementGeometry::SubControlVolume;
43 using SubControlVolumeFace = typename FVElementGeometry::SubControlVolumeFace;
45 using ElementFluxVariablesCache = typename GetPropType<TypeTag, Properties::GridFluxVariablesCache>::LocalView;
46 using ElementVolumeVariables = typename GetPropType<TypeTag, Properties::GridVolumeVariables>::LocalView;
47 using VolumeVariables = typename ElementVolumeVariables::VolumeVariables;
48
49 // class assembling the stress tensor
51
52public:
53 using ParentType::ParentType;
54
63 NumEqVector computeStorage(const Problem& problem,
64 const SubControlVolume& scv,
65 const VolumeVariables& volVars) const
66 {
67 return NumEqVector(0.0);
68 }
69
81 NumEqVector computeFlux(const Problem& problem,
82 const Element& element,
83 const FVElementGeometry& fvGeometry,
84 const ElementVolumeVariables& elemVolVars,
85 const SubControlVolumeFace& scvf,
86 const ElementFluxVariablesCache& elemFluxVarsCache) const
87 {
88 // obtain force on the face from stress type
89 return StressType::force(problem, element, fvGeometry, elemVolVars, scvf, elemFluxVarsCache);
90 }
91
105 NumEqVector computeSource(const Problem& problem,
106 const Element& element,
107 const FVElementGeometry& fvGeometry,
108 const ElementVolumeVariables& elemVolVars,
109 const SubControlVolume &scv) const
110 {
111 NumEqVector source(0.0);
112
113 // add contributions from volume flux sources
114 source += problem.source(element, fvGeometry, elemVolVars, scv);
115
116 // add contribution from possible point sources
117 source += problem.scvPointSources(element, fvGeometry, elemVolVars, scv);
118
119 // maybe add gravitational acceleration
120 static const bool gravity = getParamFromGroup<bool>(problem.paramGroup(), "Problem.EnableGravity");
121 if (gravity)
122 {
123 const auto& g = problem.spatialParams().gravity(scv.center());
124 for (int dir = 0; dir < GridView::dimensionworld; ++dir)
125 source[Indices::momentum(dir)] += elemVolVars[scv].solidDensity()*g[dir];
126 }
127
128 return source;
129 }
130};
131
132} // end namespace Dumux
133
134#endif
Element-wise calculation of the local residual for problems using the elastic model considering linea...
Definition solidmechanics/elastic/localresidual.hh:32
NumEqVector computeFlux(const Problem &problem, const Element &element, const FVElementGeometry &fvGeometry, const ElementVolumeVariables &elemVolVars, const SubControlVolumeFace &scvf, const ElementFluxVariablesCache &elemFluxVarsCache) const
Evaluates the force in all grid directions acting on a face of a sub-control volume.
Definition solidmechanics/elastic/localresidual.hh:81
NumEqVector computeStorage(const Problem &problem, const SubControlVolume &scv, const VolumeVariables &volVars) const
For the elastic model the storage term is zero since we neglect inertia forces.
Definition solidmechanics/elastic/localresidual.hh:63
NumEqVector computeSource(const Problem &problem, const Element &element, const FVElementGeometry &fvGeometry, const ElementVolumeVariables &elemVolVars, const SubControlVolume &scv) const
Calculate the source term of the equation.
Definition solidmechanics/elastic/localresidual.hh:105
Defines all properties used in Dumux.
The default local operator than can be specialized for each discretization scheme.
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
T getParamFromGroup(Args &&... args)
A free function to get a parameter from the parameter tree singleton with a model group.
Definition parameters.hh:149
typename GetProp< TypeTag, Property >::type GetPropType
get the type alias defined in the property
Definition propertysystem.hh:296
Definition adapt.hh:17
typename Detail::DiscretizationDefaultLocalOperator< TypeTag >::type DiscretizationDefaultLocalOperator
Definition defaultlocaloperator.hh:26
A helper to deduce a vector with the same size as numbers of equations.
The infrastructure to retrieve run-time parameters from Dune::ParameterTrees.