Skip to main content

conspire/constitutive/solid/elastic_hyperviscous/
mod.rs

1//! Elastic-hyperviscous solid 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 canonical;
14
15use super::{
16    super::fluid::hyperviscous::Hyperviscous,
17    viscoelastic::{AppliedLoad, Viscoelastic},
18    *,
19};
20use crate::{
21    math::{
22        ContractWith, Matrix, Quantity, Vector,
23        integrate::ImplicitDaeMinimize,
24        optimize::{EqualityConstraint, Optimization},
25    },
26    units::{Dissipation, Time},
27};
28
29/// Required methods for elastic-hyperviscous solid constitutive models.
30pub trait ElasticHyperviscous
31where
32    Self: Viscoelastic + Hyperviscous,
33{
34    /// Calculates and returns the dissipation potential.
35    ///
36    /// ```math
37    /// \phi(\mathbf{F},\dot{\mathbf{F}}) = \mathbf{P}^e(\mathbf{F}):\dot{\mathbf{F}} + \psi(\mathbf{F},\dot{\mathbf{F}})
38    /// ```
39    fn dissipation_potential(
40        &self,
41        deformation_gradient: &DeformationGradient,
42        deformation_gradient_rate: &DeformationGradientRate,
43    ) -> Result<Quantity<Dissipation>, ConstitutiveError> {
44        Ok(self
45            .first_piola_kirchhoff_stress(deformation_gradient, &DeformationGradientRate::zero())?
46            .contract_with(deformation_gradient_rate)
47            + self.viscous_dissipation(deformation_gradient, deformation_gradient_rate)?)
48    }
49    /// Calculates and returns the internal dissipation.
50    ///
51    /// ```math
52    /// T\dot{s} = \mathbf{P}(\mathbf{F},\dot{\mathbf{F}}):\dot{\mathbf{F}} - \mathbf{P}^e(\mathbf{F}):\dot{\mathbf{F}}
53    /// ```
54    fn internal_dissipation(
55        &self,
56        deformation_gradient: &DeformationGradient,
57        deformation_gradient_rate: &DeformationGradientRate,
58    ) -> Result<Quantity<Dissipation>, ConstitutiveError> {
59        Ok(self
60            .first_piola_kirchhoff_stress(deformation_gradient, deformation_gradient_rate)?
61            .contract_with(deformation_gradient_rate)
62            - self
63                .first_piola_kirchhoff_stress(
64                    deformation_gradient,
65                    &DeformationGradientRate::zero(),
66                )?
67                .contract_with(deformation_gradient_rate))
68    }
69}
70
71/// Optimization methods for elastic-hyperviscous solid constitutive models.
72pub trait Minimize {
73    /// Solve for the unknown components of the deformation gradient and rate under an applied load.
74    ///
75    /// ```math
76    /// \Pi(\mathbf{F},\dot{\mathbf{F}},\boldsymbol{\lambda}) = \mathbf{P}^e(\mathbf{F}):\dot{\mathbf{F}} + \psi(\mathbf{F},\dot{\mathbf{F}}) - \boldsymbol{\lambda}:(\dot{\mathbf{F}} - \dot{\mathbf{F}}_0) - \mathbf{P}_0:\dot{\mathbf{F}}
77    /// ```
78    fn minimize(
79        &self,
80        applied_load: AppliedLoad,
81        integrator: impl ImplicitDaeMinimize<
82            Quantity<Dissipation>,
83            FirstPiolaKirchhoffStress,
84            FirstPiolaKirchhoffRateTangentStiffness,
85            DeformationGradient,
86            DeformationGradients,
87            DeformationGradientRates,
88        >,
89        solver: impl Optimization<
90            Quantity<Dissipation>,
91            FirstPiolaKirchhoffStress,
92            FirstPiolaKirchhoffRateTangentStiffness,
93            DeformationGradientRate,
94        >,
95    ) -> Result<(Times, DeformationGradients, DeformationGradientRates), ConstitutiveError>;
96}
97
98impl<T> Minimize for T
99where
100    T: ElasticHyperviscous,
101{
102    fn minimize(
103        &self,
104        applied_load: AppliedLoad,
105        integrator: impl ImplicitDaeMinimize<
106            Quantity<Dissipation>,
107            FirstPiolaKirchhoffStress,
108            FirstPiolaKirchhoffRateTangentStiffness,
109            DeformationGradient,
110            DeformationGradients,
111            DeformationGradientRates,
112        >,
113        solver: impl Optimization<
114            Quantity<Dissipation>,
115            FirstPiolaKirchhoffStress,
116            FirstPiolaKirchhoffRateTangentStiffness,
117            DeformationGradientRate,
118        >,
119    ) -> Result<(Times, DeformationGradients, DeformationGradientRates), ConstitutiveError> {
120        match applied_load {
121            AppliedLoad::UniaxialStress(deformation_gradient_rate_11, time) => {
122                let mut matrix = Matrix::zero(4, 9);
123                let mut vector = Vector::zero(4);
124                matrix[0][0] = 1.0;
125                matrix[1][1] = 1.0;
126                matrix[2][2] = 1.0;
127                matrix[3][5] = 1.0;
128                integrator.integrate(
129                    |_: Quantity<Time>,
130                     deformation_gradient: &DeformationGradient,
131                     deformation_gradient_rate: &DeformationGradientRate| {
132                        Ok(self.dissipation_potential(
133                            deformation_gradient,
134                            deformation_gradient_rate,
135                        )?)
136                    },
137                    |_: Quantity<Time>,
138                     deformation_gradient: &DeformationGradient,
139                     deformation_gradient_rate: &DeformationGradientRate| {
140                        Ok(self.first_piola_kirchhoff_stress(
141                            deformation_gradient,
142                            deformation_gradient_rate,
143                        )?)
144                    },
145                    |_: Quantity<Time>,
146                     deformation_gradient: &DeformationGradient,
147                     deformation_gradient_rate: &DeformationGradientRate| {
148                        Ok(self.first_piola_kirchhoff_rate_tangent_stiffness(
149                            deformation_gradient,
150                            deformation_gradient_rate,
151                        )?)
152                    },
153                    solver,
154                    time,
155                    DeformationGradient::identity(),
156                    |t: Quantity<Time>| {
157                        vector[0] = deformation_gradient_rate_11(t);
158                        EqualityConstraint::Linear(matrix.clone(), vector.clone())
159                    },
160                    None,
161                )
162            }
163            AppliedLoad::BiaxialStress(
164                deformation_gradient_rate_11,
165                deformation_gradient_rate_22,
166                time,
167            ) => {
168                let mut matrix = Matrix::zero(5, 9);
169                let mut vector = Vector::zero(5);
170                matrix[0][0] = 1.0;
171                matrix[1][1] = 1.0;
172                matrix[2][2] = 1.0;
173                matrix[3][5] = 1.0;
174                matrix[4][4] = 1.0;
175                integrator.integrate(
176                    |_: Quantity<Time>,
177                     deformation_gradient: &DeformationGradient,
178                     deformation_gradient_rate: &DeformationGradientRate| {
179                        Ok(self.dissipation_potential(
180                            deformation_gradient,
181                            deformation_gradient_rate,
182                        )?)
183                    },
184                    |_: Quantity<Time>,
185                     deformation_gradient: &DeformationGradient,
186                     deformation_gradient_rate: &DeformationGradientRate| {
187                        Ok(self.first_piola_kirchhoff_stress(
188                            deformation_gradient,
189                            deformation_gradient_rate,
190                        )?)
191                    },
192                    |_: Quantity<Time>,
193                     deformation_gradient: &DeformationGradient,
194                     deformation_gradient_rate: &DeformationGradientRate| {
195                        Ok(self.first_piola_kirchhoff_rate_tangent_stiffness(
196                            deformation_gradient,
197                            deformation_gradient_rate,
198                        )?)
199                    },
200                    solver,
201                    time,
202                    DeformationGradient::identity(),
203                    |t: Quantity<Time>| {
204                        vector[0] = deformation_gradient_rate_11(t);
205                        vector[4] = deformation_gradient_rate_22(t);
206                        EqualityConstraint::Linear(matrix.clone(), vector.clone())
207                    },
208                    None,
209                )
210            }
211        }
212        .map_err(|error| ConstitutiveError::upstream(error, self))
213    }
214}