conspire/constitutive/solid/elastic_hyperviscous/
mod.rs#[cfg(test)]
pub mod test;
mod almansi_hamel;
pub use almansi_hamel::AlmansiHamel;
use super::{super::fluid::viscous::Viscous, viscoelastic::Viscoelastic, *};
use crate::math::optimize::{NewtonRaphson, SecondOrder};
use std::fmt::Debug;
pub trait ElasticHyperviscous<'a>
where
Self: Viscoelastic<'a> + Debug,
{
fn dissipation_potential(
&self,
deformation_gradient: &DeformationGradient,
deformation_gradient_rate: &DeformationGradientRate,
) -> Result<Scalar, ConstitutiveError> {
Ok(self
.first_piola_kirchhoff_stress(deformation_gradient, &ZERO_10)?
.full_contraction(deformation_gradient_rate)
+ self.viscous_dissipation(deformation_gradient, deformation_gradient_rate)?)
}
fn viscous_dissipation(
&self,
deformation_gradient: &DeformationGradient,
deformation_gradient_rate: &DeformationGradientRate,
) -> Result<Scalar, ConstitutiveError>;
fn solve_uniaxial<const W: usize>(
&self,
deformation_gradient_rate_11: impl Fn(Scalar) -> Scalar,
evaluation_times: [Scalar; W],
) -> Result<(DeformationGradients<W>, CauchyStresses<W>), ConstitutiveError> {
let mut cauchy_stresses = CauchyStresses::zero();
let mut deformation_gradients = DeformationGradients::identity();
let time_steps = evaluation_times.windows(2).map(|time| time[1] - time[0]);
for ((index, time_step), time) in time_steps.enumerate().zip(evaluation_times.into_iter()) {
(deformation_gradients[index + 1], cauchy_stresses[index + 1]) = self
.solve_uniaxial_inner(
&deformation_gradients[index],
deformation_gradient_rate_11(time),
time_step,
)?;
}
Ok((deformation_gradients, cauchy_stresses))
}
#[doc(hidden)]
fn solve_uniaxial_inner(
&self,
deformation_gradient_previous: &DeformationGradient,
deformation_gradient_rate_11: Scalar,
time_step: Scalar,
) -> Result<(DeformationGradient, CauchyStress), ConstitutiveError> {
let optimization = NewtonRaphson {
..Default::default()
};
let deformation_gradient = optimization.minimize(
|deformation_gradient: &DeformationGradient| {
let (deformation_gradient_rate, _) = self.solve_uniaxial_inner_inner(
deformation_gradient,
&deformation_gradient_rate_11,
)?;
Ok(deformation_gradient.clone()
- deformation_gradient_previous
- &deformation_gradient_rate * time_step)
},
|deformation_gradient: &DeformationGradient| {
let (deformation_gradient_rate, _) = self.solve_uniaxial_inner_inner(
deformation_gradient,
&deformation_gradient_rate_11,
)?;
Ok(IDENTITY_1010
- TensorRank4::dyad_ik_jl(
&(&deformation_gradient_rate * deformation_gradient.inverse()),
&IDENTITY_00,
) * time_step)
},
IDENTITY_10,
None,
None,
)?;
let (_, cauchy_stress) =
self.solve_uniaxial_inner_inner(&deformation_gradient, &deformation_gradient_rate_11)?;
Ok((deformation_gradient, cauchy_stress))
}
#[doc(hidden)]
fn solve_uniaxial_inner_inner(
&self,
deformation_gradient: &DeformationGradient,
deformation_gradient_rate_11: &Scalar,
) -> Result<(DeformationGradientRate, CauchyStress), ConstitutiveError> {
let optimization = NewtonRaphson {
..Default::default()
};
let deformation_gradient_rate_33 = optimization.minimize(
|deformation_gradient_rate_33: &Scalar| {
Ok(self.cauchy_stress(
deformation_gradient,
&DeformationGradientRate::new([
[*deformation_gradient_rate_11, 0.0, 0.0],
[0.0, *deformation_gradient_rate_33, 0.0],
[0.0, 0.0, *deformation_gradient_rate_33],
]),
)?[2][2])
},
|deformation_gradient_rate_33: &Scalar| {
Ok(self.cauchy_rate_tangent_stiffness(
deformation_gradient,
&DeformationGradientRate::new([
[*deformation_gradient_rate_11, 0.0, 0.0],
[0.0, *deformation_gradient_rate_33, 0.0],
[0.0, 0.0, *deformation_gradient_rate_33],
]),
)?[2][2][2][2])
},
-deformation_gradient_rate_11 / deformation_gradient[0][0].powf(1.5),
None,
None,
)?;
let deformation_gradient_rate = DeformationGradientRate::new([
[*deformation_gradient_rate_11, 0.0, 0.0],
[0.0, deformation_gradient_rate_33, 0.0],
[0.0, 0.0, deformation_gradient_rate_33],
]);
let cauchy_stress = self.cauchy_stress(deformation_gradient, &deformation_gradient_rate)?;
Ok((deformation_gradient_rate, cauchy_stress))
}
}