PolyFEM
Loading...
Searching...
No Matches
HookeLinearElasticity.hpp
Go to the documentation of this file.
1#pragma once
2
7
8// local assembler for HookeLinearElasticity C : (F+F^T)/2, see linear elasticity
9namespace polyfem::assembler
10{
12 {
13 public:
18
20
21 // res is R^{dimĀ²}
22 Eigen::Matrix<double, Eigen::Dynamic, 1, 0, 9, 1>
23 assemble(const LinearAssemblerData &data) const override;
24
25 // compute elastic energy
26 double compute_energy(const NonLinearAssemblerData &data) const override;
27 // neccessary for mixing linear model with non-linear collision response
28 Eigen::MatrixXd assemble_hessian(const NonLinearAssemblerData &data) const override;
29 // compute gradient of elastic energy, as assembler
30 Eigen::VectorXd assemble_gradient(const NonLinearAssemblerData &data) const override;
31
32 VectorNd compute_rhs(const AutodiffHessianPt &pt) const override;
33
34 void set_size(const int size) override;
35
36 // sets the elasticty tensor
37 void add_multimaterial(const int index, const json &params, const Units &units) override;
38
40
41 virtual bool is_linear() const override { return true; }
42 bool allow_inversion() const override { return true; }
43 std::string name() const override { return "HookeLinearElasticity"; }
44 std::map<std::string, ParamFunc> parameters() const override;
45
46 void assign_stress_tensor(const OutputData &data,
47 const int all_size,
48 const ElasticityTensorType &type,
49 Eigen::MatrixXd &all,
50 const std::function<Eigen::MatrixXd(const Eigen::MatrixXd &)> &fun) const override;
51
52 private:
54
55 // aux function that computes energy
56 // double compute_energy is the same with T=double
57 // assemble_gradient is the same with T=DScalar1 and return .getGradient()
58 template <typename T>
59 T compute_energy_aux(const NonLinearAssemblerData &data) const;
60 };
61} // namespace polyfem::assembler
Eigen::Matrix< double, Eigen::Dynamic, 1, 0, 9, 1 > assemble(const LinearAssemblerData &data) const override
local assembly function that defines the bilinear form (LHS) computes and returns a single local stif...
double compute_energy(const NonLinearAssemblerData &data) const override
std::map< std::string, ParamFunc > parameters() const override
const ElasticityTensor & elasticity_tensor() const
void add_multimaterial(const int index, const json &params, const Units &units) override
VectorNd compute_rhs(const AutodiffHessianPt &pt) const override
Eigen::MatrixXd assemble_hessian(const NonLinearAssemblerData &data) const override
T compute_energy_aux(const NonLinearAssemblerData &data) const
Eigen::VectorXd assemble_gradient(const NonLinearAssemblerData &data) const override
void assign_stress_tensor(const OutputData &data, const int all_size, const ElasticityTensorType &type, Eigen::MatrixXd &all, const std::function< Eigen::MatrixXd(const Eigen::MatrixXd &)> &fun) const override
assemble matrix based on the local assembler local assembler is eg Laplace, LinearElasticity etc
void assemble(const bool is_volume, const int n_basis, const std::vector< basis::ElementBases > &bases, const std::vector< basis::ElementBases > &gbases, const AssemblyValsCache &cache, const double t, StiffnessMatrix &stiffness, const bool is_mass=false) const override
assembles the stiffness matrix for the given basis the bilinear form (local assembler) is encoded by ...
double assemble_energy(const bool is_volume, const std::vector< basis::ElementBases > &bases, const std::vector< basis::ElementBases > &gbases, const AssemblyValsCache &cache, const double t, const double dt, const Eigen::MatrixXd &displacement, const Eigen::MatrixXd &displacement_prev) const override
void assemble_gradient(const bool is_volume, const int n_basis, const std::vector< basis::ElementBases > &bases, const std::vector< basis::ElementBases > &gbases, const AssemblyValsCache &cache, const double t, const double dt, const Eigen::MatrixXd &displacement, const Eigen::MatrixXd &displacement_prev, Eigen::MatrixXd &rhs) const override
void assemble_hessian(const bool is_volume, const int n_basis, const bool project_to_psd, const std::vector< basis::ElementBases > &bases, const std::vector< basis::ElementBases > &gbases, const AssemblyValsCache &cache, const double t, const double dt, const Eigen::MatrixXd &displacement, const Eigen::MatrixXd &displacement_prev, utils::MatrixCache &mat_cache, StiffnessMatrix &grad) const override
Used for test only.
Eigen::Matrix< AutodiffScalarHessian, Eigen::Dynamic, 1, 0, 3, 1 > AutodiffHessianPt
Eigen::Matrix< double, Eigen::Dynamic, 1, 0, 3, 1 > VectorNd
Definition Types.hpp:11
nlohmann::json json
Definition Common.hpp:9