conspire/domain/fem/solid/elastic_hyperviscous/
mod.rs

1use crate::{
2    fem::{
3        Blocks, ElementModel, ElementModelError, Elements, Model, NodalCoordinates,
4        NodalCoordinatesHistory, NodalVelocities, NodalVelocitiesHistory,
5        block::{band_from_neighbors, finalize_node_neighbors},
6        solid::{NodalForcesSolid, NodalStiffnessesSolid, viscoelastic::ViscoelasticElements},
7    },
8    math::{
9        Scalar, Tensor,
10        integrate::{ImplicitDaeSecondOrderMinimize, IntegrationError},
11        optimize::{EqualityConstraint, SecondOrderOptimization},
12    },
13    mechanics::Times,
14};
15
16pub trait ElasticHyperviscousElements<const D: usize>
17where
18    Self: ViscoelasticElements<D>,
19{
20    fn viscous_dissipation(
21        &self,
22        nodal_coordinates: &NodalCoordinates<D>,
23        nodal_velocities: &NodalVelocities<D>,
24    ) -> Result<Scalar, ElementModelError>;
25    fn dissipation_potential(
26        &self,
27        nodal_coordinates: &NodalCoordinates<D>,
28        nodal_velocities: &NodalVelocities<D>,
29    ) -> Result<Scalar, ElementModelError>;
30}
31
32impl<B, const D: usize> ElasticHyperviscousElements<D> for Model<B, D>
33where
34    B: ElasticHyperviscousElements<D>,
35{
36    fn viscous_dissipation(
37        &self,
38        nodal_coordinates: &NodalCoordinates<D>,
39        nodal_velocities: &NodalVelocities<D>,
40    ) -> Result<Scalar, ElementModelError> {
41        self.blocks
42            .viscous_dissipation(nodal_coordinates, nodal_velocities)
43    }
44    fn dissipation_potential(
45        &self,
46        nodal_coordinates: &NodalCoordinates<D>,
47        nodal_velocities: &NodalVelocities<D>,
48    ) -> Result<Scalar, ElementModelError> {
49        self.blocks
50            .dissipation_potential(nodal_coordinates, nodal_velocities)
51    }
52}
53
54impl<B1, B2, const D: usize> ElasticHyperviscousElements<D> for Blocks<B1, B2>
55where
56    B1: ElasticHyperviscousElements<D>,
57    B2: ElasticHyperviscousElements<D>,
58{
59    fn viscous_dissipation(
60        &self,
61        nodal_coordinates: &NodalCoordinates<D>,
62        nodal_velocities: &NodalVelocities<D>,
63    ) -> Result<Scalar, ElementModelError> {
64        Ok(self
65            .0
66            .viscous_dissipation(nodal_coordinates, nodal_velocities)?
67            + self
68                .1
69                .viscous_dissipation(nodal_coordinates, nodal_velocities)?)
70    }
71    fn dissipation_potential(
72        &self,
73        nodal_coordinates: &NodalCoordinates<D>,
74        nodal_velocities: &NodalVelocities<D>,
75    ) -> Result<Scalar, ElementModelError> {
76        Ok(self
77            .0
78            .dissipation_potential(nodal_coordinates, nodal_velocities)?
79            + self
80                .1
81                .dissipation_potential(nodal_coordinates, nodal_velocities)?)
82    }
83}
84
85pub trait SecondOrderMinimize<const D: usize> {
86    fn minimize(
87        &self,
88        equality_constraint: EqualityConstraint,
89        integrator: impl ImplicitDaeSecondOrderMinimize<
90            Scalar,
91            NodalForcesSolid<D>,
92            NodalStiffnessesSolid<D>,
93            NodalVelocities<D>,
94            NodalVelocitiesHistory<D>,
95        >,
96        time: &[Scalar],
97        solver: impl SecondOrderOptimization<
98            Scalar,
99            NodalForcesSolid<D>,
100            NodalStiffnessesSolid<D>,
101            NodalCoordinates<D>,
102        >,
103    ) -> Result<(Times, NodalCoordinatesHistory<D>, NodalVelocitiesHistory<D>), IntegrationError>;
104}
105
106impl<B, const D: usize> SecondOrderMinimize<D> for Model<B, D>
107where
108    B: ElasticHyperviscousElements<D>,
109{
110    fn minimize(
111        &self,
112        equality_constraint: EqualityConstraint,
113        integrator: impl ImplicitDaeSecondOrderMinimize<
114            Scalar,
115            NodalForcesSolid<D>,
116            NodalStiffnessesSolid<D>,
117            NodalVelocities<D>,
118            NodalVelocitiesHistory<D>,
119        >,
120        time: &[Scalar],
121        solver: impl SecondOrderOptimization<
122            Scalar,
123            NodalForcesSolid<D>,
124            NodalStiffnessesSolid<D>,
125            NodalCoordinates<D>,
126        >,
127    ) -> Result<(Times, NodalCoordinatesHistory<D>, NodalVelocitiesHistory<D>), IntegrationError>
128    {
129        let mut neighbors = vec![Vec::new(); self.coordinates().len()];
130        self.node_neighbors(&mut neighbors);
131        finalize_node_neighbors(&mut neighbors);
132        let banded = band_from_neighbors(&neighbors, &equality_constraint, D);
133        integrator.integrate(
134            |_: Scalar,
135             nodal_coordinates: &NodalCoordinates<D>,
136             nodal_velocities: &NodalVelocities<D>| {
137                Ok(self.dissipation_potential(nodal_coordinates, nodal_velocities)?)
138            },
139            |_: Scalar,
140             nodal_coordinates: &NodalCoordinates<D>,
141             nodal_velocities: &NodalVelocities<D>| {
142                Ok(self.nodal_forces(nodal_coordinates, nodal_velocities)?)
143            },
144            |_: Scalar,
145             nodal_coordinates: &NodalCoordinates<D>,
146             nodal_velocities: &NodalVelocities<D>| {
147                Ok(self.nodal_stiffnesses(nodal_coordinates, nodal_velocities)?)
148            },
149            solver,
150            time,
151            self.coordinates().clone().into(),
152            |_: Scalar| equality_constraint.clone(),
153            Some(banded),
154        )
155    }
156}
conspire/domain/fem/solid/elastic_hyperviscous/mod.rs

conspire/domain/fem/solid/elastic_hyperviscous/
mod.rs
