conspire/domain/vem/block/solid/elastic/
mod.rs

1use crate::{
2    constitutive::solid::elastic::Elastic,
3    math::{
4        Tensor,
5        optimize::{
6            EqualityConstraint, FirstOrderRootFinding, OptimizationError, ZerothOrderRootFinding,
7        },
8    },
9    vem::{
10        NodalCoordinates,
11        block::{
12            Block, FirstOrderRoot, VirtualElementBlockError, ZerothOrderRoot,
13            element::{VirtualElementError, solid::elastic::ElasticVirtualElement},
14            solid::{NodalForcesSolid, NodalStiffnessesSolid, SolidVirtualElementBlock},
15        },
16    },
17};
18
19pub trait ElasticVirtualElementBlock<C, F>
20where
21    C: Elastic,
22    F: ElasticVirtualElement<C>,
23{
24    fn nodal_forces(
25        &self,
26        nodal_coordinates: &NodalCoordinates,
27    ) -> Result<NodalForcesSolid, VirtualElementBlockError>;
28    fn nodal_stiffnesses(
29        &self,
30        nodal_coordinates: &NodalCoordinates,
31    ) -> Result<NodalStiffnessesSolid, VirtualElementBlockError>;
32}
33
34impl<C, F> ElasticVirtualElementBlock<C, F> for Block<C, F>
35where
36    C: Elastic,
37    F: ElasticVirtualElement<C>,
38    Self: SolidVirtualElementBlock<C, F>,
39{
40    fn nodal_forces(
41        &self,
42        nodal_coordinates: &NodalCoordinates,
43    ) -> Result<NodalForcesSolid, VirtualElementBlockError> {
44        let mut nodal_forces = NodalForcesSolid::zero(nodal_coordinates.len());
45        match self
46            .elements()
47            .iter()
48            .zip(self.elements_nodes())
49            .try_for_each(|(element, nodes)| {
50                element
51                    .nodal_forces(
52                        self.constitutive_model(),
53                        Self::element_coordinates(nodal_coordinates, nodes),
54                    )?
55                    .iter()
56                    .zip(nodes)
57                    .for_each(|(nodal_force, &node)| nodal_forces[node] += nodal_force);
58                Ok::<(), VirtualElementError>(())
59            }) {
60            Ok(()) => Ok(nodal_forces),
61            Err(error) => Err(VirtualElementBlockError::Upstream(
62                format!("{error}"),
63                format!("{self:?}"),
64            )),
65        }
66    }
67    fn nodal_stiffnesses(
68        &self,
69        nodal_coordinates: &NodalCoordinates,
70    ) -> Result<NodalStiffnessesSolid, VirtualElementBlockError> {
71        let mut nodal_stiffnesses = NodalStiffnessesSolid::zero(nodal_coordinates.len());
72        match self
73            .elements()
74            .iter()
75            .zip(self.elements_nodes())
76            .try_for_each(|(element, nodes)| {
77                element
78                    .nodal_stiffnesses(
79                        self.constitutive_model(),
80                        Self::element_coordinates(nodal_coordinates, nodes),
81                    )?
82                    .iter()
83                    .zip(nodes)
84                    .for_each(|(object, &node_a)| {
85                        object
86                            .iter()
87                            .zip(nodes)
88                            .for_each(|(nodal_stiffness, &node_b)| {
89                                nodal_stiffnesses[node_a][node_b] += nodal_stiffness
90                            })
91                    });
92                Ok::<(), VirtualElementError>(())
93            }) {
94            Ok(()) => Ok(nodal_stiffnesses),
95            Err(error) => Err(VirtualElementBlockError::Upstream(
96                format!("{error}"),
97                format!("{self:?}"),
98            )),
99        }
100    }
101}
102
103impl<C, F> ZerothOrderRoot<C, F, NodalCoordinates> for Block<C, F>
104where
105    C: Elastic,
106    F: ElasticVirtualElement<C>,
107{
108    fn root(
109        &self,
110        equality_constraint: EqualityConstraint,
111        solver: impl ZerothOrderRootFinding<NodalCoordinates>,
112    ) -> Result<NodalCoordinates, OptimizationError> {
113        solver.root(
114            |nodal_coordinates: &NodalCoordinates| Ok(self.nodal_forces(nodal_coordinates)?),
115            self.coordinates().clone().into(),
116            equality_constraint,
117        )
118    }
119}
120
121impl<C, F> FirstOrderRoot<C, F, NodalForcesSolid, NodalStiffnessesSolid, NodalCoordinates>
122    for Block<C, F>
123where
124    C: Elastic,
125    F: ElasticVirtualElement<C>,
126{
127    fn root(
128        &self,
129        equality_constraint: EqualityConstraint,
130        solver: impl FirstOrderRootFinding<NodalForcesSolid, NodalStiffnessesSolid, NodalCoordinates>,
131    ) -> Result<NodalCoordinates, OptimizationError> {
132        solver.root(
133            |nodal_coordinates: &NodalCoordinates| Ok(self.nodal_forces(nodal_coordinates)?),
134            |nodal_coordinates: &NodalCoordinates| Ok(self.nodal_stiffnesses(nodal_coordinates)?),
135            self.coordinates().clone().into(),
136            equality_constraint,
137        )
138    }
139}