6#include <dune/localfefunctions/derivativetransformators.hh>
7#include <dune/localfefunctions/meta.hh>
15template <
typename PreFE,
typename FE>
28 template <
typename PreFE,
typename FE>
44template <
typename PreFE,
typename FE>
59 : volumeLoad_{pre.volumeLoad} {}
62 template <
template <
typename,
int,
int>
class RT>
63 requires Dune::AlwaysFalse<RT<double, 1, 1>>::value
65 Dune::PriorityTag<0>)
const {}
67 template <
typename ST>
70 const std::optional<std::reference_wrapper<
const Eigen::VectorX<ST>>>& dx = std::nullopt)
const -> ST {
74 const auto uFunction = underlying().displacementFunction(par, dx);
75 const auto geo = underlying().geometry();
76 const auto& lambda = par.parameter();
78 for (
const auto& [gpIndex, gp] : uFunction.viewOverIntegrationPoints()) {
79 const auto uVal = uFunction.evaluate(gpIndex);
80 Eigen::Vector<double, worldDim> fext = volumeLoad_(geo.global(gp.position()), lambda);
81 energy -= uVal.dot(fext) * geo.integrationElement(gp.position()) * gp.weight();
86 template <
typename ST>
89 const std::optional<std::reference_wrapper<
const Eigen::VectorX<ST>>>& dx = std::nullopt)
const {
92 using namespace Dune::DerivativeDirections;
94 const auto uFunction = underlying().displacementFunction(par, dx);
95 const auto geo = underlying().geometry();
98 for (
const auto& [gpIndex, gp] : uFunction.viewOverIntegrationPoints()) {
99 const Eigen::Vector<double, worldDim> fext = volumeLoad_(geo.global(gp.position()), lambda);
100 const double intElement = geo.integrationElement(gp.position()) * gp.weight();
101 for (
size_t i = 0; i < underlying().numberOfNodes(); ++i) {
102 const auto udCi = uFunction.evaluateDerivative(gpIndex, wrt(coeff(i)));
103 force.template segment<worldDim>(
worldDim * i) -= udCi * fext * intElement;
108 template <
typename ST>
111 const std::optional<std::reference_wrapper<
const Eigen::VectorX<ST>>>& dx = std::nullopt)
const {}
116 const auto& underlying()
const {
return static_cast<const FE&
>(*this); }
117 auto& underlying() {
return static_cast<FE&
>(*this); }
126template <
int worldDim>
128 const double&)>& f) {
Contains stl-like type traits.
Definition of the LinearElastic class for finite element mechanics computations.
Definition: assemblermanipulatorbuildingblocks.hh:22
auto volumeLoad(const std::function< Eigen::Vector< double, worldDim >(const Dune::FieldVector< double, worldDim > &, const double &)> &f)
A helper function to create a volume load skill.
Definition: volume.hh:127
MatrixAffordance
A strongly typed enum class representing the matrix affordance.
Definition: ferequirements.hh:64
VolumeLoadPre(F f) -> VolumeLoadPre< traits::FunctionTraits< F >::return_type::RowsAtCompileTime >
VectorAffordance
A strongly typed enum class representing the vector affordance.
Definition: ferequirements.hh:49
ScalarAffordance
A strongly typed enum class representing the scalar affordance.
Definition: ferequirements.hh:38
Definition: utils/dirichletvalues.hh:30
FE class is a base class for all finite elements.
Definition: febase.hh:79
FETraits< BH, useEigenRef, useFlat > Traits
Definition: febase.hh:38
Class representing the requirements for finite element calculations.
Definition: ferequirements.hh:223
const PM & parameter() const
Get the parameter value.
Definition: ferequirements.hh:336
Traits for handling finite elements.
Definition: fetraits.hh:25
static constexpr int worlddim
Dimension of the world space.
Definition: fetraits.hh:60
VolumeLoad class represents distributed volume load that can be applied.
Definition: volume.hh:46
VolumeLoad(const Pre &pre={})
Constructor for the Loads class.
Definition: volume.hh:58
static constexpr int worldDim
Definition: volume.hh:50
auto calculateScalarImpl(const Requirement &par, ScalarAffordance affordance, const std::optional< std::reference_wrapper< const Eigen::VectorX< ST > > > &dx=std::nullopt) const -> ST
Definition: volume.hh:68
::value auto calculateAtImpl(const Requirement &req, const Dune::FieldVector< double, Traits::mydim > &local, Dune::PriorityTag< 0 >) const
Definition: volume.hh:64
void calculateVectorImpl(const Requirement &par, VectorAffordance affordance, typename Traits::template VectorType< ST > force, const std::optional< std::reference_wrapper< const Eigen::VectorX< ST > > > &dx=std::nullopt) const
Definition: volume.hh:87
void calculateMatrixImpl(const Requirement &par, MatrixAffordance, typename Traits::template MatrixType<> K, const std::optional< std::reference_wrapper< const Eigen::VectorX< ST > > > &dx=std::nullopt) const
Definition: volume.hh:109
A PreFE struct for volume load skill.
Definition: volume.hh:24
std::function< Eigen::Vector< double, worldDim >(const Dune::FieldVector< double, worldDim > &, const double &)> volumeLoad
Definition: volume.hh:26
static constexpr int worldDim
Definition: volume.hh:25
Definition: utils/dirichletvalues.hh:32