conspire/constitutive/hybrid/elastic/multiplicative/elastic/
mod.rs

1#[cfg(test)]
2mod test;
3
4use crate::{
5    constitutive::{
6        ConstitutiveError,
7        hybrid::ElasticMultiplicative,
8        solid::{
9            Solid,
10            elastic::{Elastic, internal_variables::ElasticIV},
11        },
12    },
13    math::{ContractThirdIndexWithFirstIndexOf, Rank2, TensorArray, TensorRank4},
14    mechanics::{
15        CauchyStress, CauchyTangentStiffness, CauchyTangentStiffness1, DeformationGradient,
16        DeformationGradient2, FirstPiolaKirchhoffStress, FirstPiolaKirchhoffStress1,
17        FirstPiolaKirchhoffStress2, FirstPiolaKirchhoffTangentStiffness2, Scalar,
18        SecondPiolaKirchhoffStress,
19    },
20};
21
22impl<C1, C2> Solid for ElasticMultiplicative<C1, C2>
23where
24    C1: Elastic,
25    C2: Elastic,
26{
27    fn bulk_modulus(&self) -> Scalar {
28        todo!()
29    }
30    fn shear_modulus(&self) -> Scalar {
31        todo!()
32    }
33}
34
35impl<C1, C2>
36    ElasticIV<
37        DeformationGradient2,
38        TensorRank4<3, 2, 0, 1, 0>,
39        TensorRank4<3, 1, 0, 2, 0>,
40        FirstPiolaKirchhoffTangentStiffness2,
41    > for ElasticMultiplicative<C1, C2>
42where
43    C1: Elastic,
44    C2: Elastic,
45{
46    /// Calculates and returns the Cauchy stress.
47    ///
48    /// ```math
49    /// \boldsymbol{\sigma} = \frac{1}{J_2}\,\boldsymbol{\sigma}_1
50    /// ```
51    fn cauchy_stress(
52        &self,
53        deformation_gradient: &DeformationGradient,
54        deformation_gradient_2: &DeformationGradient2,
55    ) -> Result<CauchyStress, ConstitutiveError> {
56        let (deformation_gradient_2_inverse, jacobian_2) =
57            deformation_gradient_2.inverse_and_determinant();
58        let deformation_gradient_1 = deformation_gradient * &deformation_gradient_2_inverse;
59        Ok(self.0.cauchy_stress(&deformation_gradient_1.into())? / jacobian_2)
60    }
61    /// Calculates and returns the tangent stiffness associated with the Cauchy stress.
62    ///
63    /// ```math
64    /// \boldsymbol{\mathcal{T}} = \frac{1}{J_2}\,\boldsymbol{\mathcal{T}}_1\cdot\mathbf{F}_2^{-T}
65    /// ```
66    fn cauchy_tangent_stiffness(
67        &self,
68        deformation_gradient: &DeformationGradient,
69        deformation_gradient_2: &DeformationGradient2,
70    ) -> Result<CauchyTangentStiffness, ConstitutiveError> {
71        let (deformation_gradient_2_inverse, jacobian_2) =
72            deformation_gradient_2.inverse_and_determinant();
73        let deformation_gradient_1 = deformation_gradient * &deformation_gradient_2_inverse;
74        Ok(CauchyTangentStiffness1::from(
75            self.0
76                .cauchy_tangent_stiffness(&deformation_gradient_1.into())?,
77        ) * (deformation_gradient_2_inverse.transpose() / jacobian_2))
78    }
79    /// Calculates and returns the first Piola-Kirchhoff stress.
80    ///
81    /// ```math
82    /// \mathbf{P} = \mathbf{P}_1\cdot\mathbf{F}_2^{-T}
83    /// ```
84    fn first_piola_kirchhoff_stress(
85        &self,
86        deformation_gradient: &DeformationGradient,
87        deformation_gradient_2: &DeformationGradient2,
88    ) -> Result<FirstPiolaKirchhoffStress, ConstitutiveError> {
89        Ok(
90            self.cauchy_stress(deformation_gradient, deformation_gradient_2)?
91                * deformation_gradient.inverse_transpose()
92                * deformation_gradient.determinant(),
93        )
94    }
95    /// Calculates and returns the second Piola-Kirchhoff stress.
96    ///
97    /// ```math
98    /// \mathbf{S} = \mathbf{F}_2^{-1}\cdot\mathbf{S}_1\cdot\mathbf{F}_2^{-T}
99    /// ```
100    fn second_piola_kirchhoff_stress(
101        &self,
102        deformation_gradient: &DeformationGradient,
103        deformation_gradient_2: &DeformationGradient2,
104    ) -> Result<SecondPiolaKirchhoffStress, ConstitutiveError> {
105        Ok(deformation_gradient.inverse()
106            * self.first_piola_kirchhoff_stress(deformation_gradient, deformation_gradient_2)?)
107    }
108    fn internal_variables_initial(&self) -> DeformationGradient2 {
109        DeformationGradient2::identity()
110    }
111    /// Calculates and returns the residual associated with the second deformation gradient.
112    ///
113    /// ```math
114    /// \mathbf{R} = \mathbf{P}_2 - \mathbf{M}_1\cdot\mathbf{F}_2^{-T}
115    /// ```
116    fn internal_variables_residual(
117        &self,
118        deformation_gradient: &DeformationGradient,
119        deformation_gradient_2: &DeformationGradient2,
120    ) -> Result<DeformationGradient2, ConstitutiveError> {
121        let deformation_gradient_2_inverse = deformation_gradient_2.inverse();
122        let deformation_gradient_1 = deformation_gradient * &deformation_gradient_2_inverse;
123        Ok(FirstPiolaKirchhoffStress2::from(
124            self.1
125                .first_piola_kirchhoff_stress(deformation_gradient_2.into())?,
126        ) - deformation_gradient_1.transpose()
127            * FirstPiolaKirchhoffStress1::from(
128                self.0
129                    .first_piola_kirchhoff_stress(&deformation_gradient_1.into())?,
130            )
131            * deformation_gradient_2_inverse.transpose())
132    }
133    /// Calculates and returns the tangents associated with the internal variables.
134    ///
135    /// ```math
136    /// \frac{\partial P_{iJ}}{\partial F_{KL}^2} = -P_{iL}F_{KJ}^{2-T} - \mathcal{C}_{iJmL}F_{mK}^1
137    /// ```
138    /// ```math
139    /// \frac{\partial R_{IJ}}{\partial F_{kL}} = -F_{IL}^{2-T}P_{kJ} - F_{mI}^1\mathcal{C}_{mJkL}
140    /// ```
141    /// ```math
142    /// \frac{\partial R_{IJ}}{\partial F_{KL}^2} = \mathcal{C}_{IJKL}^2 + F_{IM}^1P_{ML}{F_{KJ}^{2-T}} - \frac{\partial R_{IJ}}{\partial F_{mL}}\,F_{mK}^1
143    /// ```
144    fn internal_variables_tangents(
145        &self,
146        deformation_gradient: &DeformationGradient,
147        deformation_gradient_2: &DeformationGradient2,
148    ) -> Result<
149        (
150            TensorRank4<3, 2, 0, 1, 0>,
151            TensorRank4<3, 1, 0, 2, 0>,
152            FirstPiolaKirchhoffTangentStiffness2,
153        ),
154        ConstitutiveError,
155    > {
156        //
157        // If hyperelastic, tangent_1 should equal a transpose of tangent_2.
158        // Could add a method to utilize that for minimize() in hyperelastic/multiplicative.
159        //
160        let deformation_gradient_2_inverse = deformation_gradient_2.inverse();
161        let deformation_gradient_2_inverse_transpose = deformation_gradient_2_inverse.transpose();
162        let deformation_gradient_1 = deformation_gradient * &deformation_gradient_2_inverse;
163        let deformation_gradient_1_transpose = deformation_gradient_1.transpose();
164        let first_piola_kirchhoff_stress =
165            self.first_piola_kirchhoff_stress(deformation_gradient, deformation_gradient_2)?;
166        let tangent_0 = self.first_piola_kirchhoff_tangent_stiffness(
167            deformation_gradient,
168            deformation_gradient_2,
169        )?;
170        let tangent_1 = TensorRank4::dyad_il_kj(
171            &(deformation_gradient_2_inverse_transpose * -1.0),
172            &first_piola_kirchhoff_stress,
173        ) - &deformation_gradient_1_transpose * &tangent_0;
174        let tangent_2 = TensorRank4::dyad_il_jk(
175            &first_piola_kirchhoff_stress,
176            &(&deformation_gradient_2_inverse * -1.0),
177        ) - tangent_0
178            .contract_third_index_with_first_index_of(&deformation_gradient_1);
179        let tangent_3 = FirstPiolaKirchhoffTangentStiffness2::from(
180            self.1
181                .first_piola_kirchhoff_tangent_stiffness(deformation_gradient_2.into())?,
182        ) - tangent_1
183            .contract_third_index_with_first_index_of(&deformation_gradient_1)
184            + TensorRank4::dyad_il_jk(
185                &(deformation_gradient_1_transpose * first_piola_kirchhoff_stress),
186                &deformation_gradient_2_inverse,
187            );
188        Ok((tangent_1, tangent_2, tangent_3))
189    }
190    fn internal_variables_constraints(&self) -> (&[usize], usize) {
191        (&[10, 11, 14], 9)
192    }
193}