12#ifndef DUMUX_SOLIDMECHANICS_HYPERELASTIC_LOCAL_RESIDUAL_HH 
   13#define DUMUX_SOLIDMECHANICS_HYPERELASTIC_LOCAL_RESIDUAL_HH 
   15#include <dune/common/fvector.hh> 
   16#include <dune/common/fmatrix.hh> 
   17#include <dune/common/std/type_traits.hh> 
   18#include <dune/common/exceptions.hh> 
   29namespace Detail::Hyperelastic {
 
   31template <
typename T, 
typename ...Ts>
 
   32using SolidDensityDetector = 
decltype(std::declval<T>().solidDensity(std::declval<Ts>()...));
 
   34template<
class T, 
typename ...Args>
 
   35static constexpr bool hasSolidDensity()
 
   36{ 
return Dune::Std::is_detected<SolidDensityDetector, T, Args...>::value; }
 
   39template <
typename T, 
typename ...Ts>
 
   40using AccelerationDetector = 
decltype(std::declval<T>().acceleration(std::declval<Ts>()...));
 
   42template<
class T, 
typename ...Args>
 
   43static constexpr bool hasAcceleration()
 
   44{ 
return Dune::Std::is_detected<AccelerationDetector, T, Args...>::value; }
 
   52template<
class TypeTag>
 
   64    using FVElementGeometry = 
typename GridGeometry::LocalView;
 
   65    using SubControlVolume = 
typename GridGeometry::SubControlVolume;
 
   66    using SubControlVolumeFace = 
typename GridGeometry::SubControlVolumeFace;
 
   67    using GridView = 
typename GridGeometry::GridView;
 
   68    using Element = 
typename GridView::template Codim<0>::Entity;
 
   71    using Indices = 
typename ModelTraits::Indices;
 
   72    static constexpr int dimWorld = GridView::dimensionworld;
 
   74    using ParentType::ParentType;
 
   77    using ParentType::evalStorage;
 
   83                     const Problem& problem,
 
   84                     const Element& element,
 
   85                     const FVElementGeometry& fvGeometry,
 
   86                     const ElementVolumeVariables& prevElemVolVars,
 
   87                     const ElementVolumeVariables& curElemVolVars,
 
   88                     const SubControlVolume& scv)
 const 
   90        const auto& curVolVars = curElemVolVars[scv];
 
   92        if constexpr (Detail::Hyperelastic::hasSolidDensity<Problem, const Element&, const SubControlVolume&>()
 
   93                      && Detail::Hyperelastic::hasAcceleration<Problem, 
const Element&, 
const SubControlVolume&, Scalar, 
decltype(curVolVars.displacement())>())
 
   95            const auto& curVolVars = curElemVolVars[scv];
 
   96            NumEqVector storage = problem.solidDensity(element, scv)*problem.acceleration(
 
   97                element, scv, this->timeLoop().timeStepSize(), curVolVars.displacement()
 
   99            storage *= curVolVars.extrusionFactor()*Extrusion::volume(fvGeometry, scv);
 
  100            residual[scv.localDofIndex()] += storage;
 
  103            DUNE_THROW(Dune::InvalidStateException,
 
  104                "Problem class must implement solidDensity and acceleration" 
  105                " methods to evaluate storage term" 
 
  114                            const Element& element,
 
  115                            const FVElementGeometry& fvGeometry,
 
  116                            const ElementVolumeVariables& elemVolVars,
 
  117                            const SubControlVolumeFace& scvf,
 
  118                            const ElementFluxVariablesCache& elemFluxVarsCache)
 const 
  123            const auto& fluxVarCache = elemFluxVarsCache[scvf];
 
  124            Dune::FieldMatrix<Scalar, dimWorld, dimWorld> F(0.0);
 
  125            for (
const auto& scv : scvs(fvGeometry))
 
  127                const auto& volVars = elemVolVars[scv];
 
  128                for (
int dir = 0; dir < dimWorld; ++dir)
 
  129                    F[dir].axpy(volVars.displacement(dir), fluxVarCache.gradN(scv.indexInElement()));
 
  132            for (
int dir = 0; dir < dimWorld; ++dir)
 
  139        NumEqVector flux(0.0);
 
  142        const auto P = problem.firstPiolaKirchhoffStressTensor(F);
 
  143        P.mv(scvf.unitOuterNormal(), flux);
 
  144        flux *= -scvf.area();
 
 
  149                              const Element& element,
 
  150                              const FVElementGeometry& fvGeometry,
 
  151                              const ElementVolumeVariables& elemVolVars,
 
  152                              const SubControlVolume& scv)
 const 
  155        return ParentType::computeSource(problem, element, fvGeometry, elemVolVars, scv);
 
 
 
Local residual for the hyperelastic model.
Definition solidmechanics/hyperelastic/localresidual.hh:55
typename ParentType::ElementResidualVector ElementResidualVector
Definition solidmechanics/hyperelastic/localresidual.hh:75
NumEqVector computeSource(const Problem &problem, const Element &element, const FVElementGeometry &fvGeometry, const ElementVolumeVariables &elemVolVars, const SubControlVolume &scv) const
Definition solidmechanics/hyperelastic/localresidual.hh:148
void evalStorage(ElementResidualVector &residual, const Problem &problem, const Element &element, const FVElementGeometry &fvGeometry, const ElementVolumeVariables &prevElemVolVars, const ElementVolumeVariables &curElemVolVars, const SubControlVolume &scv) const
Evaluate the storage contribution to the residual: rho*d^2u/dt^2.
Definition solidmechanics/hyperelastic/localresidual.hh:82
NumEqVector computeFlux(const Problem &problem, const Element &element, const FVElementGeometry &fvGeometry, const ElementVolumeVariables &elemVolVars, const SubControlVolumeFace &scvf, const ElementFluxVariablesCache &elemFluxVarsCache) const
Evaluate the fluxes over a face of a sub control volume flux: -div(P) where P is the 1st Piola-Kircho...
Definition solidmechanics/hyperelastic/localresidual.hh:113
Defines all properties used in Dumux.
The default local operator than can be specialized for each discretization scheme.
Helper classes to compute the integration elements.
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.
typename Detail::DiscretizationDefaultLocalOperator< TypeTag >::type DiscretizationDefaultLocalOperator
Definition defaultlocaloperator.hh:26
typename Extrusion< T >::type Extrusion_t
Convenience alias for obtaining the extrusion type.
Definition extrusion.hh:166
A helper to deduce a vector with the same size as numbers of equations.