conspire/fem/block/element/linear/tetrahedron/
mod.rs

1#[cfg(test)]
2mod test;
3
4use crate::{
5    fem::{
6        GradientVectors, ReferenceNodalCoordinates, StandardGradientOperators,
7        block::element::{Element, FiniteElement},
8    },
9    math::{Scalar, Scalars, TensorRank1List, tensor_rank_1},
10};
11
12#[cfg(test)]
13use crate::fem::ShapeFunctionsAtIntegrationPoints;
14
15const G: usize = 1;
16const M: usize = 3;
17const N: usize = 4;
18const P: usize = G;
19
20#[cfg(test)]
21const Q: usize = N;
22
23pub type Tetrahedron = Element<G, N>;
24
25impl FiniteElement<G, N> for Tetrahedron {
26    fn new(reference_nodal_coordinates: ReferenceNodalCoordinates<N>) -> Self {
27        let (gradient_vectors, integration_weights) = Self::initialize(reference_nodal_coordinates);
28        Self {
29            gradient_vectors,
30            integration_weights,
31        }
32    }
33    fn reference() -> ReferenceNodalCoordinates<N> {
34        ReferenceNodalCoordinates::const_from([
35            tensor_rank_1([0.0, 0.0, 0.0]),
36            tensor_rank_1([1.0, 0.0, 0.0]),
37            tensor_rank_1([0.0, 1.0, 0.0]),
38            tensor_rank_1([0.0, 0.0, 1.0]),
39        ])
40    }
41    fn reset(&mut self) {
42        let (gradient_vectors, integration_weights) = Self::initialize(Self::reference());
43        self.gradient_vectors = gradient_vectors;
44        self.integration_weights = integration_weights;
45    }
46}
47
48impl Tetrahedron {
49    fn initialize(
50        reference_nodal_coordinates: ReferenceNodalCoordinates<N>,
51    ) -> (GradientVectors<G, N>, Scalars<G>) {
52        let standard_gradient_operator = &Self::standard_gradient_operators()[0];
53        let (operator, jacobian) = (reference_nodal_coordinates * standard_gradient_operator)
54            .inverse_transpose_and_determinant();
55        let gradient_vectors = GradientVectors::const_from([operator * standard_gradient_operator]);
56        let integration_weights = Scalars::const_from([jacobian * Self::integration_weight()]);
57        (gradient_vectors, integration_weights)
58    }
59    const fn integration_weight() -> Scalar {
60        1.0 / 6.0
61    }
62    #[cfg(test)]
63    const fn shape_functions_at_integration_points() -> ShapeFunctionsAtIntegrationPoints<G, Q> {
64        ShapeFunctionsAtIntegrationPoints::const_from([tensor_rank_1([0.25; Q])])
65    }
66    const fn standard_gradient_operators() -> StandardGradientOperators<M, N, P> {
67        StandardGradientOperators::const_from([TensorRank1List::const_from([
68            tensor_rank_1([-1.0, -1.0, -1.0]),
69            tensor_rank_1([1.0, 0.0, 0.0]),
70            tensor_rank_1([0.0, 1.0, 0.0]),
71            tensor_rank_1([0.0, 0.0, 1.0]),
72        ])])
73    }
74}