conspire/domain/fem/block/element/surface/
mod.rs

1pub mod linear;
2
3use crate::{
4    fem::block::element::{
5        ElementNodalCoordinates, ElementNodalReferenceCoordinates, ElementNodalVelocities,
6        GradientVectors, StandardGradientOperators,
7    },
8    math::{IDENTITY, LEVI_CIVITA, Scalar, Scalars, Tensor, TensorArray, TensorRank2},
9    mechanics::{
10        Coordinates, Normal, NormalGradients, NormalRates, Normals, ReferenceNormals, SurfaceBases,
11    },
12};
13use std::fmt::{self, Debug, Formatter};
14
15pub struct SurfaceElement<const G: usize, const N: usize, const P: usize> {
16    gradient_vectors: GradientVectors<G, N>,
17    integration_weights: Scalars<G>,
18    reference_normals: ReferenceNormals<P>,
19}
20
21impl<const G: usize, const N: usize, const P: usize> Debug for SurfaceElement<G, N, P> {
22    fn fmt(&self, f: &mut Formatter<'_>) -> fmt::Result {
23        let element = match (G, N, P) {
24            (1, 3, 1) => "LinearTriangle",
25            _ => panic!(),
26        };
27        write!(f, "{element} {{ G: {G}, N: {N} , P: {P} }}",)
28    }
29}
30
31impl<const G: usize, const N: usize, const P: usize> SurfaceElement<G, N, P> {
32    fn gradient_vectors(&self) -> &GradientVectors<G, N> {
33        &self.gradient_vectors
34    }
35    fn integration_weights(&self) -> &Scalars<G> {
36        &self.integration_weights
37    }
38    fn reference_normals(&self) -> &ReferenceNormals<P> {
39        &self.reference_normals
40    }
41}
42
43pub trait SurfaceFiniteElement<const G: usize, const N: usize, const P: usize> {
44    fn new(
45        reference_nodal_coordinates: ElementNodalReferenceCoordinates<N>,
46        thickness: Scalar,
47    ) -> Self;
48}
49
50pub trait SurfaceFiniteElementMethods<
51    const G: usize,
52    const M: usize,
53    const N: usize,
54    const P: usize,
55> where
56    Self: SurfaceFiniteElementMethodsExtra<M, N, P>,
57{
58    fn bases<const I: usize>(nodal_coordinates: &Coordinates<I, N>) -> SurfaceBases<I, P>;
59    fn dual_bases<const I: usize>(nodal_coordinates: &Coordinates<I, N>) -> SurfaceBases<I, P>;
60    fn normals(nodal_coordinates: &ElementNodalCoordinates<N>) -> Normals<P>;
61    fn normal_gradients(nodal_coordinates: &ElementNodalCoordinates<N>) -> NormalGradients<N, P>;
62    fn normal_rates(
63        nodal_coordinates: &ElementNodalCoordinates<N>,
64        nodal_velocities: &ElementNodalVelocities<N>,
65    ) -> NormalRates<P>;
66}
67
68// make this a const fn and remove inherent impl of it once Rust stabilizes const fn trait methods
69pub trait SurfaceFiniteElementMethodsExtra<const M: usize, const N: usize, const P: usize> {
70    fn standard_gradient_operators() -> StandardGradientOperators<M, N, P>;
71}
72
73impl<const G: usize, const M: usize, const N: usize, const P: usize>
74    SurfaceFiniteElementMethods<G, M, N, P> for SurfaceElement<G, N, P>
75where
76    Self: SurfaceFiniteElementMethodsExtra<M, N, P>,
77{
78    fn bases<const I: usize>(nodal_coordinates: &Coordinates<I, N>) -> SurfaceBases<I, P> {
79        Self::standard_gradient_operators()
80            .iter()
81            .map(|standard_gradient_operator| {
82                standard_gradient_operator
83                    .iter()
84                    .zip(nodal_coordinates.iter())
85                    .map(|(standard_gradient_operator_a, nodal_coordinate_a)| {
86                        standard_gradient_operator_a
87                            .iter()
88                            .map(|standard_gradient_operator_a_m| {
89                                nodal_coordinate_a * standard_gradient_operator_a_m
90                            })
91                            .collect()
92                    })
93                    .sum()
94            })
95            .collect()
96    }
97    fn dual_bases<const I: usize>(nodal_coordinates: &Coordinates<I, N>) -> SurfaceBases<I, P> {
98        Self::bases(nodal_coordinates)
99            .iter()
100            .map(|basis_vectors| {
101                basis_vectors
102                    .iter()
103                    .map(|basis_vectors_m| {
104                        basis_vectors
105                            .iter()
106                            .map(|basis_vectors_n| basis_vectors_m * basis_vectors_n)
107                            .collect()
108                    })
109                    .collect::<TensorRank2<2, I, I>>()
110                    .inverse()
111                    .iter()
112                    .map(|metric_tensor_m| {
113                        metric_tensor_m
114                            .iter()
115                            .zip(basis_vectors.iter())
116                            .map(|(metric_tensor_mn, basis_vectors_n)| {
117                                basis_vectors_n * metric_tensor_mn
118                            })
119                            .sum()
120                    })
121                    .collect()
122            })
123            .collect()
124    }
125    fn normals(nodal_coordinates: &ElementNodalCoordinates<N>) -> Normals<P> {
126        Self::bases(nodal_coordinates)
127            .iter()
128            .map(|basis_vectors| basis_vectors[0].cross(&basis_vectors[1]).normalized())
129            .collect()
130    }
131    fn normal_gradients(nodal_coordinates: &ElementNodalCoordinates<N>) -> NormalGradients<N, P> {
132        let levi_civita_symbol = LEVI_CIVITA;
133        let mut normalization: Scalar = 0.0;
134        let mut normal_vector = Normal::zero();
135        Self::standard_gradient_operators().iter()
136        .zip(Self::bases(nodal_coordinates).iter())
137        .map(|(standard_gradient_operator, basis_vectors)|{
138            normalization = basis_vectors[0].cross(&basis_vectors[1]).norm();
139            normal_vector = basis_vectors[0].cross(&basis_vectors[1])/normalization;
140            standard_gradient_operator.iter()
141            .map(|standard_gradient_operator_a|
142                levi_civita_symbol.iter()
143                .map(|levi_civita_symbol_m|
144                    IDENTITY.iter()
145                    .zip(normal_vector.iter())
146                    .map(|(identity_i, normal_vector_i)|
147                        levi_civita_symbol_m.iter()
148                        .zip(basis_vectors[0].iter()
149                        .zip(basis_vectors[1].iter()))
150                        .map(|(levi_civita_symbol_mn, (basis_vector_0_n, basis_vector_1_n))|
151                            levi_civita_symbol_mn.iter()
152                            .zip(identity_i.iter()
153                            .zip(normal_vector.iter()))
154                            .map(|(levi_civita_symbol_mno, (identity_io, normal_vector_o))|
155                                levi_civita_symbol_mno * (identity_io - normal_vector_i * normal_vector_o)
156                            ).sum::<Scalar>() * (
157                                standard_gradient_operator_a[0] * basis_vector_1_n
158                              - standard_gradient_operator_a[1] * basis_vector_0_n
159                            )
160                        ).sum::<Scalar>() / normalization
161                    ).collect()
162                ).collect()
163            ).collect()
164        }).collect()
165    }
166    fn normal_rates(
167        nodal_coordinates: &ElementNodalCoordinates<N>,
168        nodal_velocities: &ElementNodalVelocities<N>,
169    ) -> NormalRates<P> {
170        let identity = IDENTITY;
171        let levi_civita_symbol = LEVI_CIVITA;
172        let mut normalization = 0.0;
173        Self::bases(nodal_coordinates)
174            .iter()
175            .zip(Self::normals(nodal_coordinates).iter()
176            .zip(Self::standard_gradient_operators().iter()))
177            .map(|(basis, (normal, standard_gradient_operator))| {
178                normalization = basis[0].cross(&basis[1]).norm();
179                identity.iter()
180                .zip(normal.iter())
181                .map(|(identity_i, normal_vector_i)|
182                    nodal_velocities.iter()
183                    .zip(standard_gradient_operator.iter())
184                    .map(|(nodal_velocity_a, standard_gradient_operator_a)|
185                        levi_civita_symbol.iter()
186                        .zip(nodal_velocity_a.iter())
187                        .map(|(levi_civita_symbol_m, nodal_velocity_a_m)|
188                            levi_civita_symbol_m.iter()
189                            .zip(basis[0].iter()
190                            .zip(basis[1].iter()))
191                            .map(|(levi_civita_symbol_mn, (basis_vector_0_n, basis_vector_1_n))|
192                                levi_civita_symbol_mn.iter()
193                                .zip(identity_i.iter()
194                                .zip(normal.iter()))
195                                .map(|(levi_civita_symbol_mno, (identity_io, normal_vector_o))|
196                                    levi_civita_symbol_mno * (identity_io - normal_vector_i * normal_vector_o)
197                                ).sum::<Scalar>() * (
198                                    standard_gradient_operator_a[0] * basis_vector_1_n
199                                - standard_gradient_operator_a[1] * basis_vector_0_n
200                                )
201                            ).sum::<Scalar>() * nodal_velocity_a_m
202                        ).sum::<Scalar>()
203                    ).sum::<Scalar>() / normalization
204                ).collect()
205        }).collect()
206    }
207}
conspire/domain/fem/block/element/surface/mod.rs

conspire/domain/fem/block/element/surface/
mod.rs
