conspire/fem/block/element/linear/tetrahedron/
mod.rs1#[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}