conspire/constitutive/solid/hyperelastic/
mod.rs

1//! Hyperelastic constitutive models.
2//!
3//! ---
4//!
5#![doc = include_str!("doc.md")]
6
7#[cfg(feature = "doc")]
8pub mod doc;
9
10#[cfg(test)]
11pub mod test;
12
13mod arruda_boyce;
14mod fung;
15mod gent;
16mod hencky;
17mod mooney_rivlin;
18mod neo_hookean;
19mod saint_venant_kirchhoff;
20mod yeoh;
21
22pub use self::{
23    arruda_boyce::ArrudaBoyce, fung::Fung, gent::Gent, hencky::Hencky, mooney_rivlin::MooneyRivlin,
24    neo_hookean::NeoHookean, saint_venant_kirchhoff::SaintVenantKirchhoff, yeoh::Yeoh,
25};
26use super::{
27    elastic::{AppliedLoad, Elastic},
28    *,
29};
30use crate::math::{
31    Matrix, TensorVec, Vector,
32    optimize::{
33        EqualityConstraint, FirstOrderOptimization, OptimizeError, SecondOrderOptimization,
34    },
35};
36
37/// Required methods for hyperelastic constitutive models.
38pub trait Hyperelastic
39where
40    Self: Elastic,
41{
42    /// Calculates and returns the Helmholtz free energy density.
43    ///
44    /// ```math
45    /// a = a(\mathbf{F})
46    /// ```
47    fn helmholtz_free_energy_density(
48        &self,
49        deformation_gradient: &DeformationGradient,
50    ) -> Result<Scalar, ConstitutiveError>;
51}
52
53/// First-order minimization methods for elastic constitutive models.
54pub trait FirstOrderMinimize {
55    /// Solve for the unknown components of the deformation gradient under an applied load.
56    ///
57    /// ```math
58    /// \Pi(\mathbf{F},\boldsymbol{\lambda}) = a(\mathbf{F}) - \boldsymbol{\lambda}:(\mathbf{F} - \mathbf{F}_0) - \mathbf{P}_0:\mathbf{F}
59    /// ```
60    fn minimize(
61        &self,
62        applied_load: AppliedLoad,
63        solver: impl FirstOrderOptimization<Scalar, DeformationGradient>,
64    ) -> Result<DeformationGradient, OptimizeError>;
65}
66
67/// Second-order minimization methods for elastic constitutive models.
68pub trait SecondOrderMinimize {
69    /// Solve for the unknown components of the deformation gradient under an applied load.
70    ///
71    /// ```math
72    /// \Pi(\mathbf{F},\boldsymbol{\lambda}) = a(\mathbf{F}) - \boldsymbol{\lambda}:(\mathbf{F} - \mathbf{F}_0) - \mathbf{P}_0:\mathbf{F}
73    /// ```
74    fn minimize(
75        &self,
76        applied_load: AppliedLoad,
77        solver: impl SecondOrderOptimization<
78            Scalar,
79            FirstPiolaKirchhoffStress,
80            FirstPiolaKirchhoffTangentStiffness,
81            DeformationGradient,
82        >,
83    ) -> Result<DeformationGradient, OptimizeError>;
84}
85
86impl<T> FirstOrderMinimize for T
87where
88    T: Hyperelastic,
89{
90    fn minimize(
91        &self,
92        applied_load: AppliedLoad,
93        solver: impl FirstOrderOptimization<Scalar, DeformationGradient>,
94    ) -> Result<DeformationGradient, OptimizeError> {
95        match applied_load {
96            AppliedLoad::UniaxialStress(deformation_gradient_11) => {
97                let mut matrix = Matrix::zero(4, 9);
98                let mut vector = Vector::zero(4);
99                matrix[0][0] = 1.0;
100                matrix[1][1] = 1.0;
101                matrix[2][2] = 1.0;
102                matrix[3][5] = 1.0;
103                vector[0] = deformation_gradient_11;
104                solver.minimize(
105                    |deformation_gradient: &DeformationGradient| {
106                        Ok(self.helmholtz_free_energy_density(deformation_gradient)?)
107                    },
108                    |deformation_gradient: &DeformationGradient| {
109                        Ok(self.first_piola_kirchhoff_stress(deformation_gradient)?)
110                    },
111                    DeformationGradient::identity(),
112                    EqualityConstraint::Linear(matrix, vector),
113                )
114            }
115            AppliedLoad::BiaxialStress(deformation_gradient_11, deformation_gradient_22) => {
116                let mut matrix = Matrix::zero(5, 9);
117                let mut vector = Vector::zero(5);
118                matrix[0][0] = 1.0;
119                matrix[1][1] = 1.0;
120                matrix[2][2] = 1.0;
121                matrix[3][5] = 1.0;
122                matrix[4][4] = 1.0;
123                vector[0] = deformation_gradient_11;
124                vector[4] = deformation_gradient_22;
125                solver.minimize(
126                    |deformation_gradient: &DeformationGradient| {
127                        Ok(self.helmholtz_free_energy_density(deformation_gradient)?)
128                    },
129                    |deformation_gradient: &DeformationGradient| {
130                        Ok(self.first_piola_kirchhoff_stress(deformation_gradient)?)
131                    },
132                    DeformationGradient::identity(),
133                    EqualityConstraint::Linear(matrix, vector),
134                )
135            }
136        }
137    }
138}
139
140impl<T> SecondOrderMinimize for T
141where
142    T: Hyperelastic,
143{
144    fn minimize(
145        &self,
146        applied_load: AppliedLoad,
147        solver: impl SecondOrderOptimization<
148            Scalar,
149            FirstPiolaKirchhoffStress,
150            FirstPiolaKirchhoffTangentStiffness,
151            DeformationGradient,
152        >,
153    ) -> Result<DeformationGradient, OptimizeError> {
154        match applied_load {
155            AppliedLoad::UniaxialStress(deformation_gradient_11) => {
156                let mut matrix = Matrix::zero(4, 9);
157                let mut vector = Vector::zero(4);
158                matrix[0][0] = 1.0;
159                matrix[1][1] = 1.0;
160                matrix[2][2] = 1.0;
161                matrix[3][5] = 1.0;
162                vector[0] = deformation_gradient_11;
163                solver.minimize(
164                    |deformation_gradient: &DeformationGradient| {
165                        Ok(self.helmholtz_free_energy_density(deformation_gradient)?)
166                    },
167                    |deformation_gradient: &DeformationGradient| {
168                        Ok(self.first_piola_kirchhoff_stress(deformation_gradient)?)
169                    },
170                    |deformation_gradient: &DeformationGradient| {
171                        Ok(self.first_piola_kirchhoff_tangent_stiffness(deformation_gradient)?)
172                    },
173                    DeformationGradient::identity(),
174                    EqualityConstraint::Linear(matrix, vector),
175                    None,
176                )
177            }
178            AppliedLoad::BiaxialStress(deformation_gradient_11, deformation_gradient_22) => {
179                let mut matrix = Matrix::zero(5, 9);
180                let mut vector = Vector::zero(5);
181                matrix[0][0] = 1.0;
182                matrix[1][1] = 1.0;
183                matrix[2][2] = 1.0;
184                matrix[3][5] = 1.0;
185                matrix[4][4] = 1.0;
186                vector[0] = deformation_gradient_11;
187                vector[4] = deformation_gradient_22;
188                solver.minimize(
189                    |deformation_gradient: &DeformationGradient| {
190                        Ok(self.helmholtz_free_energy_density(deformation_gradient)?)
191                    },
192                    |deformation_gradient: &DeformationGradient| {
193                        Ok(self.first_piola_kirchhoff_stress(deformation_gradient)?)
194                    },
195                    |deformation_gradient: &DeformationGradient| {
196                        Ok(self.first_piola_kirchhoff_tangent_stiffness(deformation_gradient)?)
197                    },
198                    DeformationGradient::identity(),
199                    EqualityConstraint::Linear(matrix, vector),
200                    None,
201                )
202            }
203        }
204    }
205}