conspire/constitutive/solid/hyperelastic/yeoh/
mod.rs1#[cfg(test)]
2mod test;
3
4use crate::{
5 constitutive::{
6 Constitutive, ConstitutiveError,
7 solid::{FIVE_THIRDS, Solid, TWO_THIRDS, elastic::Elastic, hyperelastic::Hyperelastic},
8 },
9 math::{IDENTITY, Rank2},
10 mechanics::{CauchyStress, CauchyTangentStiffness, Deformation, DeformationGradient, Scalar},
11};
12use std::iter::once;
13
14const SEVEN_THIRDS: Scalar = 7.0 / 3.0;
15
16#[doc = include_str!("doc.md")]
17#[derive(Debug)]
18pub struct Yeoh<const N: usize> {
19 pub bulk_modulus: Scalar,
21 pub shear_modulus: Scalar,
23 pub extra_moduli: [Scalar; N],
25}
26
27impl<const N: usize> Yeoh<N> {
28 pub fn extra_moduli(&self) -> &[Scalar] {
30 &self.extra_moduli
31 }
32}
33
34impl<const N: usize> Solid for Yeoh<N> {
35 fn bulk_modulus(&self) -> &Scalar {
36 &self.bulk_modulus
37 }
38 fn shear_modulus(&self) -> &Scalar {
39 &self.shear_modulus
40 }
41}
42
43impl<const N: usize> Elastic for Yeoh<N> {
44 #[doc = include_str!("cauchy_stress.md")]
45 fn cauchy_stress(
46 &self,
47 deformation_gradient: &DeformationGradient,
48 ) -> Result<CauchyStress, ConstitutiveError> {
49 let jacobian = self.jacobian(deformation_gradient)?;
50 let (deviatoric_left_cauchy_green_deformation, left_cauchy_green_deformation_trace) =
51 deformation_gradient
52 .left_cauchy_green()
53 .deviatoric_and_trace();
54 let scalar_term = left_cauchy_green_deformation_trace / jacobian.powf(TWO_THIRDS) - 3.0;
55 Ok(deviatoric_left_cauchy_green_deformation
56 * once(self.shear_modulus())
57 .chain(self.extra_moduli().iter())
58 .enumerate()
59 .map(|(n, modulus)| ((n as Scalar) + 1.0) * modulus * scalar_term.powi(n as i32))
60 .sum::<Scalar>()
61 / jacobian.powf(FIVE_THIRDS)
62 + IDENTITY * self.bulk_modulus() * 0.5 * (jacobian - 1.0 / jacobian))
63 }
64 #[doc = include_str!("cauchy_tangent_stiffness.md")]
65 fn cauchy_tangent_stiffness(
66 &self,
67 deformation_gradient: &DeformationGradient,
68 ) -> Result<CauchyTangentStiffness, ConstitutiveError> {
69 let jacobian = self.jacobian(deformation_gradient)?;
70 let inverse_transpose_deformation_gradient = deformation_gradient.inverse_transpose();
71 let left_cauchy_green_deformation = deformation_gradient.left_cauchy_green();
72 let scalar_term = left_cauchy_green_deformation.trace() / jacobian.powf(TWO_THIRDS) - 3.0;
73 let scaled_modulus = once(self.shear_modulus())
74 .chain(self.extra_moduli().iter())
75 .enumerate()
76 .map(|(n, modulus)| ((n as Scalar) + 1.0) * modulus * scalar_term.powi(n as i32))
77 .sum::<Scalar>()
78 / jacobian.powf(FIVE_THIRDS);
79 let deviatoric_left_cauchy_green_deformation = left_cauchy_green_deformation.deviatoric();
80 let last_term = CauchyTangentStiffness::dyad_ij_kl(
81 &deviatoric_left_cauchy_green_deformation,
82 &((left_cauchy_green_deformation.deviatoric()
83 * &inverse_transpose_deformation_gradient)
84 * (2.0
85 * self
86 .extra_moduli()
87 .iter()
88 .enumerate()
89 .map(|(n, modulus)| {
90 ((n as Scalar) + 2.0)
91 * ((n as Scalar) + 1.0)
92 * modulus
93 * scalar_term.powi(n as i32)
94 })
95 .sum::<Scalar>()
96 / jacobian.powf(SEVEN_THIRDS))),
97 );
98 Ok(
99 (CauchyTangentStiffness::dyad_ik_jl(&IDENTITY, deformation_gradient)
100 + CauchyTangentStiffness::dyad_il_jk(deformation_gradient, &IDENTITY)
101 - CauchyTangentStiffness::dyad_ij_kl(&IDENTITY, deformation_gradient)
102 * (TWO_THIRDS))
103 * scaled_modulus
104 + CauchyTangentStiffness::dyad_ij_kl(
105 &(IDENTITY * (0.5 * self.bulk_modulus() * (jacobian + 1.0 / jacobian))
106 - deviatoric_left_cauchy_green_deformation
107 * (scaled_modulus * FIVE_THIRDS)),
108 &inverse_transpose_deformation_gradient,
109 )
110 + last_term,
111 )
112 }
113}
114
115impl<const N: usize> Hyperelastic for Yeoh<N> {
116 #[doc = include_str!("helmholtz_free_energy_density.md")]
117 fn helmholtz_free_energy_density(
118 &self,
119 deformation_gradient: &DeformationGradient,
120 ) -> Result<Scalar, ConstitutiveError> {
121 let jacobian = self.jacobian(deformation_gradient)?;
122 let scalar_term =
123 deformation_gradient.left_cauchy_green().trace() / jacobian.powf(TWO_THIRDS) - 3.0;
124 Ok(0.5
125 * (once(self.shear_modulus())
126 .chain(self.extra_moduli().iter())
127 .enumerate()
128 .map(|(n, modulus)| modulus * scalar_term.powi((n + 1) as i32))
129 .sum::<Scalar>()
130 + self.bulk_modulus() * (0.5 * (jacobian.powi(2) - 1.0) - jacobian.ln())))
131 }
132}