conspire/domain/fem/block/element/solid/elastic/
mod.rs1pub mod internal_variables;
2
3use crate::{
4 constitutive::solid::elastic::Elastic,
5 fem::block::element::{
6 Element, ElementNodalCoordinates, FiniteElement, FiniteElementError, GradientVectors,
7 solid::{ElementNodalForcesSolid, ElementNodalStiffnessesSolid, SolidFiniteElement},
8 surface::{SurfaceElement, SurfaceFiniteElement},
9 },
10 math::{ContractSecondFourthWithFirst, Current, IDENTITY, Quantity, Tensor, TensorRank2},
11 mechanics::{FirstPiolaKirchhoffStressList, FirstPiolaKirchhoffTangentStiffnessList},
12 units::StressPerArea,
13};
14
15pub trait ElasticFiniteElement<C, const G: usize, const M: usize, const N: usize, const P: usize>
16where
17 C: Elastic,
18 Self: SolidFiniteElement<G, M, N, P>,
19{
20 fn nodal_forces(
21 &self,
22 constitutive_model: &C,
23 nodal_coordinates: &ElementNodalCoordinates<N>,
24 ) -> Result<ElementNodalForcesSolid<N>, FiniteElementError>;
25 fn nodal_stiffnesses(
26 &self,
27 constitutive_model: &C,
28 nodal_coordinates: &ElementNodalCoordinates<N>,
29 ) -> Result<ElementNodalStiffnessesSolid<N>, FiniteElementError>;
30}
31
32impl<C, const G: usize, const N: usize, const O: usize, const P: usize>
33 ElasticFiniteElement<C, G, 3, N, P> for Element<3, G, N, O>
34where
35 C: Elastic,
36 Self: SolidFiniteElement<G, 3, N, P>,
37{
38 fn nodal_forces(
39 &self,
40 constitutive_model: &C,
41 nodal_coordinates: &ElementNodalCoordinates<N>,
42 ) -> Result<ElementNodalForcesSolid<N>, FiniteElementError> {
43 nodal_forces::<_, _, _, _, _, O, _>(
44 self,
45 constitutive_model,
46 self.gradient_vectors(),
47 nodal_coordinates,
48 )
49 }
50 fn nodal_stiffnesses(
51 &self,
52 constitutive_model: &C,
53 nodal_coordinates: &ElementNodalCoordinates<N>,
54 ) -> Result<ElementNodalStiffnessesSolid<N>, FiniteElementError> {
55 let first_piola_kirchhoff_tangent_stiffnesses = self
56 .deformation_gradients(nodal_coordinates)
57 .iter()
58 .map(|deformation_gradient| {
59 constitutive_model.first_piola_kirchhoff_tangent_stiffness(deformation_gradient)
60 })
61 .collect::<Result<FirstPiolaKirchhoffTangentStiffnessList<G>, _>>()
62 .map_err(|error| FiniteElementError::upstream(error, self))?;
63 Ok(first_piola_kirchhoff_tangent_stiffnesses
64 .iter()
65 .zip(
66 self.gradient_vectors()
67 .iter()
68 .zip(self.integration_weights()),
69 )
70 .map(
71 |(
72 first_piola_kirchhoff_tangent_stiffness,
73 (gradient_vectors, integration_weight),
74 )| {
75 gradient_vectors
76 .iter()
77 .map(|gradient_vector_a| {
78 gradient_vectors
79 .iter()
80 .map(|gradient_vector_b| {
81 first_piola_kirchhoff_tangent_stiffness
82 .contract_second_fourth_with_first(
83 gradient_vector_a,
84 gradient_vector_b,
85 )
86 * integration_weight
87 })
88 .collect()
89 })
90 .collect()
91 },
92 )
93 .sum())
94 }
95}
96
97impl<C, const G: usize, const N: usize, const O: usize> ElasticFiniteElement<C, G, 2, N, N>
98 for SurfaceElement<G, N, O>
99where
100 C: Elastic,
101 Self: SolidFiniteElement<G, 2, N, N>,
102{
103 fn nodal_forces(
104 &self,
105 constitutive_model: &C,
106 nodal_coordinates: &ElementNodalCoordinates<N>,
107 ) -> Result<ElementNodalForcesSolid<N>, FiniteElementError> {
108 nodal_forces::<_, _, _, _, _, O, _>(
109 self,
110 constitutive_model,
111 self.gradient_vectors(),
112 nodal_coordinates,
113 )
114 }
115 fn nodal_stiffnesses(
116 &self,
117 constitutive_model: &C,
118 nodal_coordinates: &ElementNodalCoordinates<N>,
119 ) -> Result<ElementNodalStiffnessesSolid<N>, FiniteElementError> {
120 let first_piola_kirchhoff_tangent_stiffnesses = self
121 .deformation_gradients(nodal_coordinates)
122 .iter()
123 .map(|deformation_gradient| {
124 constitutive_model.first_piola_kirchhoff_tangent_stiffness(deformation_gradient)
125 })
126 .collect::<Result<FirstPiolaKirchhoffTangentStiffnessList<G>, _>>()
127 .map_err(|error| FiniteElementError::upstream(error, self))?;
128 Ok(first_piola_kirchhoff_tangent_stiffnesses
129 .iter()
130 .zip(
131 self.gradient_vectors()
132 .iter()
133 .zip(self.integration_weights().iter()
134 .zip(self.reference_normals().iter()
135 .zip(Self::normal_gradients(nodal_coordinates))
136 )
137 ),
138 )
139 .map(
140 |(
141 first_piola_kirchhoff_tangent_stiffness,
142 (gradient_vectors, (integration_weight, (reference_normal, normal_gradients))),
143 )| {
144 gradient_vectors.iter()
145 .map(|gradient_vector_a|
146 gradient_vectors.iter()
147 .zip(normal_gradients.iter())
148 .map(|(gradient_vector_b, normal_gradient_b)|
149 first_piola_kirchhoff_tangent_stiffness.iter()
150 .map(|first_piola_kirchhoff_tangent_stiffness_m|
151 IDENTITY.iter()
152 .zip(normal_gradient_b.iter())
153 .map(|(identity_n, normal_gradient_b_n)|
154 first_piola_kirchhoff_tangent_stiffness_m.iter()
155 .zip(gradient_vector_a.iter())
156 .map(|(first_piola_kirchhoff_tangent_stiffness_mj, gradient_vector_a_j)|
157 first_piola_kirchhoff_tangent_stiffness_mj.iter()
158 .zip(identity_n.iter()
159 .zip(normal_gradient_b_n.iter()))
160 .map(|(first_piola_kirchhoff_tangent_stiffness_mjk, (identity_nk, normal_gradient_b_n_k))|
161 first_piola_kirchhoff_tangent_stiffness_mjk.iter()
162 .zip(gradient_vector_b.iter()
163 .zip(reference_normal.iter()))
164 .map(|(first_piola_kirchhoff_tangent_stiffness_mjkl, (gradient_vector_b_l, reference_normal_l))|
165 first_piola_kirchhoff_tangent_stiffness_mjkl * gradient_vector_a_j * (
166 identity_nk * gradient_vector_b_l + normal_gradient_b_n_k * reference_normal_l
167 )
168 ).sum::<Quantity<StressPerArea>>()
169 ).sum::<Quantity<StressPerArea>>()
170 ).sum::<Quantity<StressPerArea>>()
171 ).collect()
172 ).collect::<TensorRank2<3, Current, Current, StressPerArea>>() * integration_weight
173 ).collect()
174 ).collect()
175 }
176 )
177 .sum())
178 }
179}
180
181fn nodal_forces<
182 C,
183 F,
184 const G: usize,
185 const M: usize,
186 const N: usize,
187 const O: usize,
188 const P: usize,
189>(
190 element: &F,
191 constitutive_model: &C,
192 gradient_vectors: &GradientVectors<3, G, N>,
193 nodal_coordinates: &ElementNodalCoordinates<N>,
194) -> Result<ElementNodalForcesSolid<N>, FiniteElementError>
195where
196 C: Elastic,
197 F: SolidFiniteElement<G, M, N, P>,
198{
199 let first_piola_kirchhoff_stresses = element
200 .deformation_gradients(nodal_coordinates)
201 .iter()
202 .map(|deformation_gradient| {
203 constitutive_model.first_piola_kirchhoff_stress(deformation_gradient)
204 })
205 .collect::<Result<FirstPiolaKirchhoffStressList<G>, _>>()
206 .map_err(|error| FiniteElementError::upstream(error, element))?;
207 Ok(first_piola_kirchhoff_stresses
208 .iter()
209 .zip(gradient_vectors.iter().zip(element.integration_weights()))
210 .map(
211 |(first_piola_kirchhoff_stress, (gradient_vectors, integration_weight))| {
212 gradient_vectors
213 .iter()
214 .map(|gradient_vector| {
215 (first_piola_kirchhoff_stress * gradient_vector) * integration_weight
216 })
217 .collect()
218 },
219 )
220 .sum())
221}