conspire/domain/fem/solid/hyperelastic/
mod.rs1pub mod internal_variables;
2
3use crate::{
4 fem::{
5 Blocks, ElementModel, ElementModelError, Elements, FirstOrderMinimize, Model,
6 NodalCoordinates, SecondOrderMinimize,
7 block::{finalize_node_neighbors, solver_from_neighbors},
8 solid::{NodalForcesSolid, NodalStiffnessesSolidSymmetric, elastic::ElasticElements},
9 },
10 math::{
11 Quantity, Tensor,
12 optimize::{
13 EqualityConstraint, FirstOrderOptimization, OptimizationError, SecondOrderOptimization,
14 },
15 },
16 units::Energy,
17};
18
19pub trait HyperelasticElements<const D: usize>
20where
21 Self: ElasticElements<D>,
22{
23 fn helmholtz_free_energy(
24 &self,
25 nodal_coordinates: &NodalCoordinates<D>,
26 ) -> Result<Quantity<Energy>, ElementModelError>;
27 fn nodal_stiffnesses_symmetric_into(
28 &self,
29 nodal_coordinates: &NodalCoordinates<D>,
30 nodal_stiffnesses: &mut NodalStiffnessesSolidSymmetric<D>,
31 ) -> Result<(), ElementModelError>;
32 fn nodal_stiffnesses_symmetric(
33 &self,
34 nodal_coordinates: &NodalCoordinates<D>,
35 ) -> Result<NodalStiffnessesSolidSymmetric<D>, ElementModelError> {
36 let mut nodal_stiffnesses = NodalStiffnessesSolidSymmetric::zero(nodal_coordinates.len());
37 self.nodal_stiffnesses_symmetric_into(nodal_coordinates, &mut nodal_stiffnesses)?;
38 Ok(nodal_stiffnesses)
39 }
40}
41
42impl<B, const D: usize> HyperelasticElements<D> for Model<B, D>
43where
44 B: HyperelasticElements<D>,
45{
46 fn helmholtz_free_energy(
47 &self,
48 nodal_coordinates: &NodalCoordinates<D>,
49 ) -> Result<Quantity<Energy>, ElementModelError> {
50 self.blocks.helmholtz_free_energy(nodal_coordinates)
51 }
52 fn nodal_stiffnesses_symmetric_into(
53 &self,
54 nodal_coordinates: &NodalCoordinates<D>,
55 nodal_stiffnesses: &mut NodalStiffnessesSolidSymmetric<D>,
56 ) -> Result<(), ElementModelError> {
57 self.blocks
58 .nodal_stiffnesses_symmetric_into(nodal_coordinates, nodal_stiffnesses)
59 }
60}
61
62impl<B1, B2, const D: usize> HyperelasticElements<D> for Blocks<B1, B2>
63where
64 B1: HyperelasticElements<D>,
65 B2: HyperelasticElements<D>,
66{
67 fn helmholtz_free_energy(
68 &self,
69 nodal_coordinates: &NodalCoordinates<D>,
70 ) -> Result<Quantity<Energy>, ElementModelError> {
71 Ok(self.0.helmholtz_free_energy(nodal_coordinates)?
72 + self.1.helmholtz_free_energy(nodal_coordinates)?)
73 }
74 fn nodal_stiffnesses_symmetric_into(
75 &self,
76 nodal_coordinates: &NodalCoordinates<D>,
77 nodal_stiffnesses: &mut NodalStiffnessesSolidSymmetric<D>,
78 ) -> Result<(), ElementModelError> {
79 self.0
80 .nodal_stiffnesses_symmetric_into(nodal_coordinates, nodal_stiffnesses)?;
81 self.1
82 .nodal_stiffnesses_symmetric_into(nodal_coordinates, nodal_stiffnesses)
83 }
84}
85
86impl<B, const D: usize>
87 FirstOrderMinimize<Quantity<Energy>, NodalForcesSolid<D>, NodalCoordinates<D>> for Model<B, D>
88where
89 B: HyperelasticElements<D>,
90{
91 fn minimize(
92 &self,
93 equality_constraint: EqualityConstraint,
94 solver: impl FirstOrderOptimization<Quantity<Energy>, NodalForcesSolid<D>, NodalCoordinates<D>>,
95 ) -> Result<NodalCoordinates<D>, OptimizationError> {
96 solver.minimize(
97 |nodal_coordinates: &NodalCoordinates<D>| {
98 Ok(self.helmholtz_free_energy(nodal_coordinates)?)
99 },
100 |nodal_coordinates: &NodalCoordinates<D>| Ok(self.nodal_forces(nodal_coordinates)?),
101 self.coordinates().clone().into(),
102 equality_constraint,
103 )
104 }
105}
106
107impl<B, const D: usize>
108 SecondOrderMinimize<
109 Quantity<Energy>,
110 NodalForcesSolid<D>,
111 NodalStiffnessesSolidSymmetric<D>,
112 NodalCoordinates<D>,
113 > for Model<B, D>
114where
115 B: HyperelasticElements<D>,
116{
117 fn minimize(
118 &self,
119 equality_constraint: EqualityConstraint,
120 solver: impl SecondOrderOptimization<
121 Quantity<Energy>,
122 NodalForcesSolid<D>,
123 NodalStiffnessesSolidSymmetric<D>,
124 NodalCoordinates<D>,
125 >,
126 ) -> Result<NodalCoordinates<D>, OptimizationError> {
127 let mut neighbors = vec![Vec::new(); self.coordinates().len()];
128 self.node_neighbors(&mut neighbors);
129 finalize_node_neighbors(&mut neighbors);
130 let sparse = solver_from_neighbors(&neighbors, &equality_constraint, D, true);
131 solver.minimize(
132 |nodal_coordinates: &NodalCoordinates<D>| {
133 Ok(self.helmholtz_free_energy(nodal_coordinates)?)
134 },
135 |nodal_coordinates: &NodalCoordinates<D>| Ok(self.nodal_forces(nodal_coordinates)?),
136 |nodal_coordinates: &NodalCoordinates<D>| {
137 Ok(self.nodal_stiffnesses_symmetric(nodal_coordinates)?)
138 },
139 self.coordinates().clone().into(),
140 equality_constraint,
141 Some(sparse),
142 )
143 }
144}