conspire/constitutive/solid/elastic_viscoplastic/
mod.rs

1//! Elastic-viscoplastic solid constitutive models.
2
3use crate::{
4    constitutive::{
5        ConstitutiveError,
6        fluid::viscoplastic::{
7            Viscoplastic, ViscoplasticStateVariables, ViscoplasticStateVariablesHistory,
8        },
9    },
10    math::{
11        Matrix, Rank2, Tensor, TensorArray, Vector,
12        integrate::{ExplicitDaeFirstOrderRoot, ExplicitDaeZerothOrderRoot},
13        optimize::{EqualityConstraint, FirstOrderRootFinding, ZerothOrderRootFinding},
14    },
15    mechanics::{
16        DeformationGradient, DeformationGradients, FirstPiolaKirchhoffStress,
17        FirstPiolaKirchhoffTangentStiffness, Scalar, Times,
18    },
19};
20
21pub use crate::constitutive::solid::elastic_plastic::{AppliedLoad, ElasticPlasticOrViscoplastic};
22
23/// Required methods for elastic-viscoplastic solid constitutive models.
24pub trait ElasticViscoplastic<Y>
25where
26    Self: ElasticPlasticOrViscoplastic + Viscoplastic<Y>,
27    Y: Tensor,
28{
29    /// Calculates and returns the evolution of the state variables.
30    fn state_variables_evolution(
31        &self,
32        deformation_gradient: &DeformationGradient,
33        state_variables: &ViscoplasticStateVariables<Y>,
34    ) -> Result<ViscoplasticStateVariables<Y>, ConstitutiveError> {
35        let deformation_gradient_p = &state_variables.0;
36        let jacobian = self.jacobian(deformation_gradient)?;
37        let deformation_gradient_e = deformation_gradient * deformation_gradient_p.inverse();
38        let cauchy_stress = self.cauchy_stress(deformation_gradient, deformation_gradient_p)?;
39        let mandel_stress = (deformation_gradient_e.transpose()
40            * cauchy_stress
41            * deformation_gradient_e.inverse_transpose())
42            * jacobian;
43        self.plastic_evolution(mandel_stress, state_variables)
44    }
45}
46
47/// Zeroth-order root-finding methods for elastic-viscoplastic solid constitutive models.
48pub trait ZerothOrderRoot<Y>
49where
50    Y: Tensor,
51{
52    /// Solve for the unknown components of the deformation gradients under an applied load.
53    ///
54    /// ```math
55    /// \mathbf{P}(\mathbf{F},\mathbf{F}_\mathrm{p}) - \boldsymbol{\lambda} - \mathbf{P}_0 = \mathbf{0}
56    /// ```
57    fn root(
58        &self,
59        applied_load: AppliedLoad,
60        integrator: impl ExplicitDaeZerothOrderRoot<
61            ViscoplasticStateVariables<Y>,
62            DeformationGradient,
63            ViscoplasticStateVariablesHistory<Y>,
64            DeformationGradients,
65        >,
66        solver: impl ZerothOrderRootFinding<DeformationGradient>,
67    ) -> Result<
68        (
69            Times,
70            DeformationGradients,
71            ViscoplasticStateVariablesHistory<Y>,
72        ),
73        ConstitutiveError,
74    >;
75}
76
77/// First-order root-finding methods for elastic-viscoplastic solid constitutive models.
78pub trait FirstOrderRoot<Y>
79where
80    Y: Tensor,
81{
82    /// Solve for the unknown components of the deformation gradients under an applied load.
83    ///
84    /// ```math
85    /// \mathbf{P}(\mathbf{F},\mathbf{F}_\mathrm{p}) - \boldsymbol{\lambda} - \mathbf{P}_0 = \mathbf{0}
86    /// ```
87    fn root(
88        &self,
89        applied_load: AppliedLoad,
90        integrator: impl ExplicitDaeFirstOrderRoot<
91            FirstPiolaKirchhoffStress,
92            FirstPiolaKirchhoffTangentStiffness,
93            ViscoplasticStateVariables<Y>,
94            DeformationGradient,
95            ViscoplasticStateVariablesHistory<Y>,
96            DeformationGradients,
97        >,
98        solver: impl FirstOrderRootFinding<
99            FirstPiolaKirchhoffStress,
100            FirstPiolaKirchhoffTangentStiffness,
101            DeformationGradient,
102        >,
103    ) -> Result<
104        (
105            Times,
106            DeformationGradients,
107            ViscoplasticStateVariablesHistory<Y>,
108        ),
109        ConstitutiveError,
110    >;
111}
112
113impl<C, Y> ZerothOrderRoot<Y> for C
114where
115    C: ElasticViscoplastic<Y>,
116    Y: Tensor,
117{
118    fn root(
119        &self,
120        applied_load: AppliedLoad,
121        integrator: impl ExplicitDaeZerothOrderRoot<
122            ViscoplasticStateVariables<Y>,
123            DeformationGradient,
124            ViscoplasticStateVariablesHistory<Y>,
125            DeformationGradients,
126        >,
127        solver: impl ZerothOrderRootFinding<DeformationGradient>,
128    ) -> Result<
129        (
130            Times,
131            DeformationGradients,
132            ViscoplasticStateVariablesHistory<Y>,
133        ),
134        ConstitutiveError,
135    > {
136        match match applied_load {
137            AppliedLoad::UniaxialStress(deformation_gradient_11, time) => {
138                let mut matrix = Matrix::zero(4, 9);
139                let mut vector = Vector::zero(4);
140                matrix[0][0] = 1.0;
141                matrix[1][1] = 1.0;
142                matrix[2][2] = 1.0;
143                matrix[3][5] = 1.0;
144                integrator.integrate(
145                    |_: Scalar,
146                     state_variables: &ViscoplasticStateVariables<Y>,
147                     deformation_gradient: &DeformationGradient| {
148                        Ok(self.state_variables_evolution(deformation_gradient, state_variables)?)
149                    },
150                    |_: Scalar,
151                     state_variables: &ViscoplasticStateVariables<Y>,
152                     deformation_gradient: &DeformationGradient| {
153                        let deformation_gradient_p = &state_variables.0;
154                        Ok(self.first_piola_kirchhoff_stress(
155                            deformation_gradient,
156                            deformation_gradient_p,
157                        )?)
158                    },
159                    solver,
160                    time,
161                    (self.initial_state(), DeformationGradient::identity()),
162                    |t: Scalar| {
163                        vector[0] = deformation_gradient_11(t);
164                        EqualityConstraint::Linear(matrix.clone(), vector.clone())
165                    },
166                )
167            }
168        } {
169            Ok((times, state_variables, _, deformation_gradients)) => {
170                Ok((times, deformation_gradients, state_variables))
171            }
172            Err(error) => Err(ConstitutiveError::Upstream(
173                format!("{error}"),
174                format!("{self:?}"),
175            )),
176        }
177    }
178}
179
180impl<C, Y> FirstOrderRoot<Y> for C
181where
182    C: ElasticViscoplastic<Y>,
183    Y: Tensor,
184{
185    fn root(
186        &self,
187        applied_load: AppliedLoad,
188        integrator: impl ExplicitDaeFirstOrderRoot<
189            FirstPiolaKirchhoffStress,
190            FirstPiolaKirchhoffTangentStiffness,
191            ViscoplasticStateVariables<Y>,
192            DeformationGradient,
193            ViscoplasticStateVariablesHistory<Y>,
194            DeformationGradients,
195        >,
196        solver: impl FirstOrderRootFinding<
197            FirstPiolaKirchhoffStress,
198            FirstPiolaKirchhoffTangentStiffness,
199            DeformationGradient,
200        >,
201    ) -> Result<
202        (
203            Times,
204            DeformationGradients,
205            ViscoplasticStateVariablesHistory<Y>,
206        ),
207        ConstitutiveError,
208    > {
209        match match applied_load {
210            AppliedLoad::UniaxialStress(deformation_gradient_11, time) => {
211                let mut matrix = Matrix::zero(4, 9);
212                let mut vector = Vector::zero(4);
213                matrix[0][0] = 1.0;
214                matrix[1][1] = 1.0;
215                matrix[2][2] = 1.0;
216                matrix[3][5] = 1.0;
217                integrator.integrate(
218                    |_: Scalar,
219                     state_variables: &ViscoplasticStateVariables<Y>,
220                     deformation_gradient: &DeformationGradient| {
221                        Ok(self.state_variables_evolution(deformation_gradient, state_variables)?)
222                    },
223                    |_: Scalar,
224                     state_variables: &ViscoplasticStateVariables<Y>,
225                     deformation_gradient: &DeformationGradient| {
226                        let deformation_gradient_p = &state_variables.0;
227                        Ok(self.first_piola_kirchhoff_stress(
228                            deformation_gradient,
229                            deformation_gradient_p,
230                        )?)
231                    },
232                    |_: Scalar,
233                     state_variables: &ViscoplasticStateVariables<Y>,
234                     deformation_gradient: &DeformationGradient| {
235                        let deformation_gradient_p = &state_variables.0;
236                        Ok(self.first_piola_kirchhoff_tangent_stiffness(
237                            deformation_gradient,
238                            deformation_gradient_p,
239                        )?)
240                    },
241                    solver,
242                    time,
243                    (self.initial_state(), DeformationGradient::identity()),
244                    |t: Scalar| {
245                        vector[0] = deformation_gradient_11(t);
246                        EqualityConstraint::Linear(matrix.clone(), vector.clone())
247                    },
248                )
249            }
250        } {
251            Ok((times, state_variables, _, deformation_gradients)) => {
252                Ok((times, deformation_gradients, state_variables))
253            }
254            Err(error) => Err(ConstitutiveError::Upstream(
255                format!("{error}"),
256                format!("{self:?}"),
257            )),
258        }
259    }
260}