conspire/constitutive/solid/elastic/seth_hill_lagrangian/
mod.rs1#[cfg(test)]
2mod test;
3
4use super::hencky::Hencky;
5use crate::{
6 constitutive::{
7 ConstitutiveError,
8 solid::{Solid, TWO_THIRDS, elastic::Elastic},
9 },
10 math::{
11 ContractThirdFourthWithFirstSecond, IDENTITY_00, Quantity, Rank2, Spectrum, TensorRank4,
12 },
13 mechanics::{
14 Deformation, DeformationGradient, Scalar, SecondPiolaKirchhoffStress,
15 SecondPiolaKirchhoffTangentStiffness,
16 },
17 units::Stress,
18};
19
20#[doc = include_str!("doc.md")]
21#[derive(Clone, Debug)]
22pub struct SethHillLagrangian {
23 pub bulk_modulus: Quantity<Stress>,
25 pub shear_modulus: Quantity<Stress>,
27 pub exponent: Scalar,
29}
30
31impl SethHillLagrangian {
32 pub fn exponent(&self) -> Scalar {
34 self.exponent
35 }
36 fn hencky(&self) -> Hencky {
37 Hencky {
38 bulk_modulus: self.bulk_modulus,
39 shear_modulus: self.shear_modulus,
40 }
41 }
42}
43
44impl Solid for SethHillLagrangian {
45 fn bulk_modulus(&self) -> Quantity<Stress> {
46 self.bulk_modulus
47 }
48 fn shear_modulus(&self) -> Quantity<Stress> {
49 self.shear_modulus
50 }
51}
52
53impl Elastic for SethHillLagrangian {
54 #[doc = include_str!("second_piola_kirchhoff_stress.md")]
55 fn second_piola_kirchhoff_stress(
56 &self,
57 deformation_gradient: &DeformationGradient,
58 ) -> Result<SecondPiolaKirchhoffStress, ConstitutiveError> {
59 if self.exponent() == 0.0 {
60 return self
61 .hencky()
62 .second_piola_kirchhoff_stress(deformation_gradient);
63 }
64 let _jacobian = self.jacobian(deformation_gradient)?;
65 let right_cauchy_green = deformation_gradient.right_cauchy_green();
66 let (deviatoric_strain, strain_trace) = ((right_cauchy_green
67 .powm(0.5 * self.exponent())
68 .map_err(|error| ConstitutiveError::upstream(error, self))?
69 - IDENTITY_00)
70 / self.exponent())
71 .deviatoric_and_trace();
72 Ok(deviatoric_strain * (2.0 * self.shear_modulus())
73 + IDENTITY_00 * (self.bulk_modulus() * strain_trace))
74 }
75 #[doc = include_str!("second_piola_kirchhoff_tangent_stiffness.md")]
76 fn second_piola_kirchhoff_tangent_stiffness(
77 &self,
78 deformation_gradient: &DeformationGradient,
79 ) -> Result<SecondPiolaKirchhoffTangentStiffness, ConstitutiveError> {
80 if self.exponent() == 0.0 {
81 return self
82 .hencky()
83 .second_piola_kirchhoff_tangent_stiffness(deformation_gradient);
84 }
85 let _jacobian = self.jacobian(deformation_gradient)?;
86 let right_cauchy_green = deformation_gradient.right_cauchy_green();
87 let half_exponent = 0.5 * self.exponent();
88 let spectrum = Spectrum::new(&right_cauchy_green)
89 .map_err(|error| ConstitutiveError::upstream(error, self))?;
90 let deformation_gradient_transpose = deformation_gradient.transpose();
91 let scaled_deformation_gradient_transpose =
92 &deformation_gradient_transpose * (2.0 * self.shear_modulus() / self.exponent());
93 let trace_term = deformation_gradient
94 * spectrum
95 .powm(half_exponent - 1.0)
96 .map_err(|error| ConstitutiveError::upstream(error, self))?;
97 Ok(spectrum
98 .dpowm(half_exponent)
99 .map_err(|error| ConstitutiveError::upstream(error, self))?
100 .contract_third_fourth_with_first_second(
101 &(TensorRank4::dyad_il_jk(&IDENTITY_00, &scaled_deformation_gradient_transpose)
102 + TensorRank4::dyad_ik_jl(
103 &scaled_deformation_gradient_transpose,
104 &IDENTITY_00,
105 )),
106 )
107 + TensorRank4::dyad_ij_kl(
108 &(IDENTITY_00 * (self.bulk_modulus() - TWO_THIRDS * self.shear_modulus())),
109 &trace_term,
110 ))
111 }
112}