conspire/domain/fem/solid/elastic/
mod.rs1use crate::{
2 fem::{
3 Blocks, ElementModel, ElementModelError, Elements, FirstOrderRoot, Model, NodalCoordinates,
4 ZerothOrderRoot,
5 block::{finalize_node_neighbors, solver_from_neighbors},
6 solid::{NodalForcesSolid, NodalStiffnessesSolid},
7 },
8 math::{
9 Tensor,
10 optimize::{
11 EqualityConstraint, FirstOrderRootFinding, OptimizationError, ZerothOrderRootFinding,
12 },
13 },
14};
15
16pub trait ElasticElements<const D: usize>
17where
18 Self: Elements,
19{
20 fn nodal_forces_into(
21 &self,
22 nodal_coordinates: &NodalCoordinates<D>,
23 nodal_forces: &mut NodalForcesSolid<D>,
24 ) -> Result<(), ElementModelError>;
25 fn nodal_forces(
26 &self,
27 nodal_coordinates: &NodalCoordinates<D>,
28 ) -> Result<NodalForcesSolid<D>, ElementModelError> {
29 let mut nodal_forces = NodalForcesSolid::zero(nodal_coordinates.len());
30 self.nodal_forces_into(nodal_coordinates, &mut nodal_forces)?;
31 Ok(nodal_forces)
32 }
33 fn nodal_stiffnesses_into(
34 &self,
35 nodal_coordinates: &NodalCoordinates<D>,
36 nodal_stiffnesses: &mut NodalStiffnessesSolid<D>,
37 ) -> Result<(), ElementModelError>;
38 fn nodal_stiffnesses(
39 &self,
40 nodal_coordinates: &NodalCoordinates<D>,
41 ) -> Result<NodalStiffnessesSolid<D>, ElementModelError> {
42 let mut nodal_stiffnesses = NodalStiffnessesSolid::zero(nodal_coordinates.len());
43 self.nodal_stiffnesses_into(nodal_coordinates, &mut nodal_stiffnesses)?;
44 Ok(nodal_stiffnesses)
45 }
46}
47
48impl<B, const D: usize> ElasticElements<D> for Model<B, D>
49where
50 B: ElasticElements<D>,
51{
52 fn nodal_forces_into(
53 &self,
54 nodal_coordinates: &NodalCoordinates<D>,
55 nodal_forces: &mut NodalForcesSolid<D>,
56 ) -> Result<(), ElementModelError> {
57 self.blocks
58 .nodal_forces_into(nodal_coordinates, nodal_forces)
59 }
60 fn nodal_stiffnesses_into(
61 &self,
62 nodal_coordinates: &NodalCoordinates<D>,
63 nodal_stiffnesses: &mut NodalStiffnessesSolid<D>,
64 ) -> Result<(), ElementModelError> {
65 self.blocks
66 .nodal_stiffnesses_into(nodal_coordinates, nodal_stiffnesses)
67 }
68}
69
70impl<B1, B2, const D: usize> ElasticElements<D> for Blocks<B1, B2>
71where
72 B1: ElasticElements<D>,
73 B2: ElasticElements<D>,
74{
75 fn nodal_forces_into(
76 &self,
77 nodal_coordinates: &NodalCoordinates<D>,
78 nodal_forces: &mut NodalForcesSolid<D>,
79 ) -> Result<(), ElementModelError> {
80 self.0.nodal_forces_into(nodal_coordinates, nodal_forces)?;
81 self.1.nodal_forces_into(nodal_coordinates, nodal_forces)
82 }
83 fn nodal_stiffnesses_into(
84 &self,
85 nodal_coordinates: &NodalCoordinates<D>,
86 nodal_stiffnesses: &mut NodalStiffnessesSolid<D>,
87 ) -> Result<(), ElementModelError> {
88 self.0
89 .nodal_stiffnesses_into(nodal_coordinates, nodal_stiffnesses)?;
90 self.1
91 .nodal_stiffnesses_into(nodal_coordinates, nodal_stiffnesses)
92 }
93}
94
95impl<B, const D: usize> ZerothOrderRoot<NodalCoordinates<D>> for Model<B, D>
96where
97 B: ElasticElements<D>,
98{
99 fn root(
100 &self,
101 equality_constraint: EqualityConstraint,
102 solver: impl ZerothOrderRootFinding<NodalCoordinates<D>>,
103 ) -> Result<NodalCoordinates<D>, OptimizationError> {
104 solver.root(
105 |nodal_coordinates: &NodalCoordinates<D>| Ok(self.nodal_forces(nodal_coordinates)?),
106 self.coordinates().clone().into(),
107 equality_constraint,
108 )
109 }
110}
111
112impl<B, const D: usize>
113 FirstOrderRoot<NodalForcesSolid<D>, NodalStiffnessesSolid<D>, NodalCoordinates<D>>
114 for Model<B, D>
115where
116 B: ElasticElements<D>,
117{
118 fn root(
119 &self,
120 equality_constraint: EqualityConstraint,
121 solver: impl FirstOrderRootFinding<
122 NodalForcesSolid<D>,
123 NodalStiffnessesSolid<D>,
124 NodalCoordinates<D>,
125 >,
126 ) -> Result<NodalCoordinates<D>, OptimizationError> {
127 let mut neighbors = vec![Vec::new(); self.coordinates().len()];
128 self.node_neighbors(&mut neighbors);
129 finalize_node_neighbors(&mut neighbors);
130 let sparse = solver_from_neighbors(&neighbors, &equality_constraint, D, false);
131 solver.root(
132 |nodal_coordinates: &NodalCoordinates<D>| Ok(self.nodal_forces(nodal_coordinates)?),
133 |nodal_coordinates: &NodalCoordinates<D>| {
134 Ok(self.nodal_stiffnesses(nodal_coordinates)?)
135 },
136 self.coordinates().clone().into(),
137 equality_constraint,
138 Some(sparse),
139 )
140 }
141}