PolyFEM
Loading...
Searching...
No Matches
GenericElastic.hpp
Go to the documentation of this file.
1#pragma once
2
5
6// non linear NeoHookean material model
7namespace polyfem::assembler
8{
9 template <typename T>
10 using DefGradMatrix = Eigen::Matrix<T, Eigen::Dynamic, Eigen::Dynamic, 0, 3, 3>;
11
12 template <typename Derived>
14 {
15 public:
19
21 virtual ~GenericElastic() = default;
22
23 // energy, gradient, and hessian used in newton method
24 double compute_energy(const NonLinearAssemblerData &data) const override;
25 Eigen::MatrixXd assemble_hessian(const NonLinearAssemblerData &data) const override;
26 Eigen::VectorXd assemble_gradient(const NonLinearAssemblerData &data) const override;
27
28 void assign_stress_tensor(const OutputData &data,
29 const int all_size,
30 const ElasticityTensorType &type,
31 Eigen::MatrixXd &all,
32 const std::function<Eigen::MatrixXd(const Eigen::MatrixXd &)> &fun) const override;
33
35 const Eigen::MatrixXd &mat,
36 Eigen::MatrixXd &stress,
37 Eigen::MatrixXd &result) const override;
38
40 Eigen::MatrixXd &stress,
41 Eigen::MatrixXd &result) const override;
42
44 const Eigen::MatrixXd &vect,
45 Eigen::MatrixXd &stress,
46 Eigen::MatrixXd &result) const override;
47
49 Derived &derived() { return static_cast<Derived &>(*this); }
51 const Derived &derived() const { return static_cast<const Derived &>(*this); }
52
53 // sets material params
54 virtual void add_multimaterial(const int index, const json &params, const Units &units) override = 0;
55
56 bool allow_inversion() const override { return true; }
57
58 private:
59 // utility function that computes energy, the template is used for double, DScalar1, and DScalar2 in energy, gradient and hessian
60 template <typename T>
62 {
63 typedef Eigen::Matrix<T, Eigen::Dynamic, 1> AutoDiffVect;
64 typedef Eigen::Matrix<T, Eigen::Dynamic, Eigen::Dynamic, 0, 3, 3> AutoDiffGradMat;
65
66 AutoDiffVect local_disp;
67 get_local_disp(data, size(), local_disp);
68
69 AutoDiffGradMat def_grad(size(), size());
70
71 T energy = T(0.0);
72
73 const int n_pts = data.da.size();
74 for (long p = 0; p < n_pts; ++p)
75 {
76 compute_disp_grad_at_quad(data, local_disp, p, size(), def_grad);
77
78 // Id + grad d
79 for (int d = 0; d < size(); ++d)
80 def_grad(d, d) += T(1);
81
82 const T val = derived().elastic_energy(data.vals.val.row(p), data.t, data.vals.element_id, def_grad);
83
84 energy += val * data.da(p);
85 }
86 return energy;
87 }
88 };
89} // namespace polyfem::assembler
double val
Definition Assembler.cpp:86
void compute_stress_grad_multiply_stress(const OptAssemblerData &data, Eigen::MatrixXd &stress, Eigen::MatrixXd &result) 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
Derived & derived()
Returns this as a reference to derived class.
const Derived & derived() const
Returns this as a const reference to derived class.
Eigen::MatrixXd assemble_hessian(const NonLinearAssemblerData &data) const override
virtual void add_multimaterial(const int index, const json &params, const Units &units) override=0
void compute_stress_grad_multiply_vect(const OptAssemblerData &data, const Eigen::MatrixXd &vect, Eigen::MatrixXd &stress, Eigen::MatrixXd &result) const override
Eigen::VectorXd assemble_gradient(const NonLinearAssemblerData &data) const override
bool allow_inversion() const override
double compute_energy(const NonLinearAssemblerData &data) const override
T compute_energy_aux(const NonLinearAssemblerData &data) const
void compute_stress_grad_multiply_mat(const OptAssemblerData &data, const Eigen::MatrixXd &mat, Eigen::MatrixXd &stress, Eigen::MatrixXd &result) const override
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< T, Eigen::Dynamic, Eigen::Dynamic, 0, 3, 3 > DefGradMatrix
Eigen::Matrix< T, Eigen::Dynamic, Eigen::Dynamic, 0, 3, 3 > AutoDiffGradMat
Eigen::Matrix< T, Eigen::Dynamic, 1 > AutoDiffVect
void compute_disp_grad_at_quad(const assembler::NonLinearAssemblerData &data, const AutoDiffVect &local_disp, const int p, const int size, AutoDiffGradMat &def_grad)
nlohmann::json json
Definition Common.hpp:9
void get_local_disp(const assembler::NonLinearAssemblerData &data, const int size, AutoDiffVect &local_disp)