conspire/constitutive/solid/elastic_hyperviscous/
mod.rs

1//! Elastic-hyperviscous constitutive models.
2//!
3//! ---
4//!
5//! Elastic-hyperviscous constitutive models are defined by an elastic stress tensor function and a viscous dissipation function.
6//!
7//! ```math
8//! \mathbf{P}:\dot{\mathbf{F}} - \mathbf{P}^e(\mathbf{F}):\dot{\mathbf{F}} - \phi(\mathbf{F},\dot{\mathbf{F}}) \geq 0
9//! ```
10//! Satisfying the second law of thermodynamics though a minimum viscous dissipation principal yields a relation for the stress.
11//!
12//! ```math
13//! \mathbf{P} = \mathbf{P}^e + \frac{\partial\phi}{\partial\dot{\mathbf{F}}}
14//! ```
15//! Consequently, the rate tangent stiffness associated with the first Piola-Kirchhoff stress is symmetric for these constitutive models.
16//!
17//! ```math
18//! \mathcal{U}_{iJkL} = \mathcal{U}_{kLiJ}
19//! ```
20
21#[cfg(test)]
22pub mod test;
23
24mod almansi_hamel;
25
26pub use almansi_hamel::AlmansiHamel;
27
28use super::{
29    super::fluid::viscous::Viscous,
30    viscoelastic::{AppliedLoad, Viscoelastic},
31    *,
32};
33use crate::math::{
34    Matrix, TensorVec, Vector,
35    integrate::{Explicit, IntegrationError},
36    optimize::{EqualityConstraint, OptimizeError, SecondOrderOptimization},
37};
38use std::fmt::Debug;
39
40/// Required methods for elastic-hyperviscous constitutive models.
41pub trait ElasticHyperviscous
42where
43    Self: Viscoelastic,
44{
45    /// Calculates and returns the dissipation potential.
46    ///
47    /// ```math
48    /// \mathbf{P}^e(\mathbf{F}):\dot{\mathbf{F}} + \phi(\mathbf{F},\dot{\mathbf{F}})
49    /// ```
50    fn dissipation_potential(
51        &self,
52        deformation_gradient: &DeformationGradient,
53        deformation_gradient_rate: &DeformationGradientRate,
54    ) -> Result<Scalar, ConstitutiveError> {
55        Ok(self
56            .first_piola_kirchhoff_stress(deformation_gradient, &ZERO_10)?
57            .full_contraction(deformation_gradient_rate)
58            + self.viscous_dissipation(deformation_gradient, deformation_gradient_rate)?)
59    }
60    /// Calculates and returns the viscous dissipation.
61    ///
62    /// ```math
63    /// \phi = \phi(\mathbf{F},\dot{\mathbf{F}})
64    /// ```
65    fn viscous_dissipation(
66        &self,
67        deformation_gradient: &DeformationGradient,
68        deformation_gradient_rate: &DeformationGradientRate,
69    ) -> Result<Scalar, ConstitutiveError>;
70    /// Solve for the unknown components of the deformation gradient and rate under an applied load.
71    ///
72    /// ```math
73    /// \Pi(\mathbf{F},\dot{\mathbf{F}},\boldsymbol{\lambda}) = \mathbf{P}^e(\mathbf{F}):\dot{\mathbf{F}} + \phi(\mathbf{F},\dot{\mathbf{F}}) - \boldsymbol{\lambda}:(\dot{\mathbf{F}} - \dot{\mathbf{F}}_0) - \mathbf{P}_0:\dot{\mathbf{F}}
74    /// ```
75    fn minimize(
76        &self,
77        applied_load: AppliedLoad,
78        integrator: impl Explicit<DeformationGradientRate, DeformationGradientRates>,
79        solver: impl SecondOrderOptimization<
80            Scalar,
81            FirstPiolaKirchhoffStress,
82            FirstPiolaKirchhoffTangentStiffness,
83            DeformationGradient,
84        >,
85    ) -> Result<(Times, DeformationGradients, DeformationGradientRates), IntegrationError> {
86        let mut solution = DeformationGradientRate::zero();
87        match applied_load {
88            AppliedLoad::UniaxialStress(deformation_gradient_rate_11, time) => integrator
89                .integrate(
90                    |t: Scalar, deformation_gradient: &DeformationGradient| {
91                        solution = self.minimize_uniaxial_inner(
92                            deformation_gradient,
93                            deformation_gradient_rate_11(t),
94                            &solver,
95                            &solution,
96                        )?;
97                        Ok(solution.clone())
98                    },
99                    time,
100                    DeformationGradient::identity(),
101                ),
102            AppliedLoad::BiaxialStress(
103                deformation_gradient_rate_11,
104                deformation_gradient_rate_22,
105                time,
106            ) => integrator.integrate(
107                |t: Scalar, deformation_gradient: &DeformationGradient| {
108                    solution = self.minimize_biaxial_inner(
109                        deformation_gradient,
110                        deformation_gradient_rate_11(t),
111                        deformation_gradient_rate_22(t),
112                        &solver,
113                        &solution,
114                    )?;
115                    Ok(solution.clone())
116                },
117                time,
118                DeformationGradient::identity(),
119            ),
120        }
121    }
122    #[doc(hidden)]
123    fn minimize_uniaxial_inner(
124        &self,
125        deformation_gradient: &DeformationGradient,
126        deformation_gradient_rate_11: Scalar,
127        solver: &impl SecondOrderOptimization<
128            Scalar,
129            FirstPiolaKirchhoffStress,
130            FirstPiolaKirchhoffTangentStiffness,
131            DeformationGradient,
132        >,
133        initial_guess: &DeformationGradientRate,
134    ) -> Result<DeformationGradientRate, OptimizeError> {
135        let mut matrix = Matrix::zero(4, 9);
136        let mut vector = Vector::zero(4);
137        matrix[0][0] = 1.0;
138        matrix[1][1] = 1.0;
139        matrix[2][2] = 1.0;
140        matrix[3][5] = 1.0;
141        vector[0] = deformation_gradient_rate_11;
142        solver.minimize(
143            |deformation_gradient_rate: &DeformationGradientRate| {
144                Ok(self.dissipation_potential(deformation_gradient, deformation_gradient_rate)?)
145            },
146            |deformation_gradient_rate: &DeformationGradientRate| {
147                Ok(self.first_piola_kirchhoff_stress(
148                    deformation_gradient,
149                    deformation_gradient_rate,
150                )?)
151            },
152            |deformation_gradient_rate: &DeformationGradientRate| {
153                Ok(self.first_piola_kirchhoff_rate_tangent_stiffness(
154                    deformation_gradient,
155                    deformation_gradient_rate,
156                )?)
157            },
158            initial_guess.clone(),
159            EqualityConstraint::Linear(matrix, vector),
160            None,
161        )
162    }
163    #[doc(hidden)]
164    fn minimize_biaxial_inner(
165        &self,
166        deformation_gradient: &DeformationGradient,
167        deformation_gradient_rate_11: Scalar,
168        deformation_gradient_rate_22: Scalar,
169        solver: &impl SecondOrderOptimization<
170            Scalar,
171            FirstPiolaKirchhoffStress,
172            FirstPiolaKirchhoffTangentStiffness,
173            DeformationGradient,
174        >,
175        initial_guess: &DeformationGradientRate,
176    ) -> Result<DeformationGradientRate, OptimizeError> {
177        let mut matrix = Matrix::zero(5, 9);
178        let mut vector = Vector::zero(5);
179        matrix[0][0] = 1.0;
180        matrix[1][1] = 1.0;
181        matrix[2][2] = 1.0;
182        matrix[3][5] = 1.0;
183        matrix[4][4] = 1.0;
184        vector[0] = deformation_gradient_rate_11;
185        vector[4] = deformation_gradient_rate_22;
186        solver.minimize(
187            |deformation_gradient_rate: &DeformationGradientRate| {
188                Ok(self.dissipation_potential(deformation_gradient, deformation_gradient_rate)?)
189            },
190            |deformation_gradient_rate: &DeformationGradientRate| {
191                Ok(self.first_piola_kirchhoff_stress(
192                    deformation_gradient,
193                    deformation_gradient_rate,
194                )?)
195            },
196            |deformation_gradient_rate: &DeformationGradientRate| {
197                Ok(self.first_piola_kirchhoff_rate_tangent_stiffness(
198                    deformation_gradient,
199                    deformation_gradient_rate,
200                )?)
201            },
202            initial_guess.clone(),
203            EqualityConstraint::Linear(matrix, vector),
204            None,
205        )
206    }
207}
conspire/constitutive/solid/elastic_hyperviscous/mod.rs

conspire/constitutive/solid/elastic_hyperviscous/
mod.rs
