13#ifndef DUMUX_SOLIDMECHANICS_ELASTIC_LOCAL_RESIDUAL_HH 
   14#define DUMUX_SOLIDMECHANICS_ELASTIC_LOCAL_RESIDUAL_HH 
   29template<
class TypeTag>
 
   36    using GridView = 
typename GridGeometry::GridView;
 
   37    using Element = 
typename GridView::template Codim<0>::Entity;
 
   41    using FVElementGeometry = 
typename GridGeometry::LocalView;
 
   42    using SubControlVolume = 
typename FVElementGeometry::SubControlVolume;
 
   43    using SubControlVolumeFace = 
typename FVElementGeometry::SubControlVolumeFace;
 
   47    using VolumeVariables = 
typename ElementVolumeVariables::VolumeVariables;
 
   53    using ParentType::ParentType;
 
   64                               const SubControlVolume& scv,
 
   65                               const VolumeVariables& volVars)
 const 
   67        return NumEqVector(0.0);
 
 
   82                            const Element& element,
 
   83                            const FVElementGeometry& fvGeometry,
 
   84                            const ElementVolumeVariables& elemVolVars,
 
   85                            const SubControlVolumeFace& scvf,
 
   86                            const ElementFluxVariablesCache& elemFluxVarsCache)
 const 
   89        return StressType::force(problem, element, fvGeometry, elemVolVars, scvf, elemFluxVarsCache);
 
 
  106                              const Element& element,
 
  107                              const FVElementGeometry& fvGeometry,
 
  108                              const ElementVolumeVariables& elemVolVars,
 
  109                              const SubControlVolume &scv)
 const 
  111        NumEqVector source(0.0);
 
  114        source += problem.source(element, fvGeometry, elemVolVars, scv);
 
  117        source += problem.scvPointSources(element, fvGeometry, elemVolVars, scv);
 
  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];
 
 
 
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
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.