#[allow(dead_code)]
pub(super) fn edge_matrix_raw(
p0: [f64; 3],
p1: [f64; 3],
p2: [f64; 3],
p3: [f64; 3],
) -> [[f64; 3]; 3] {
[
[p1[0] - p0[0], p2[0] - p0[0], p3[0] - p0[0]],
[p1[1] - p0[1], p2[1] - p0[1], p3[1] - p0[1]],
[p1[2] - p0[2], p2[2] - p0[2], p3[2] - p0[2]],
]
}
#[allow(dead_code)]
pub(super) fn det3x3(m: [[f64; 3]; 3]) -> f64 {
m[0][0] * (m[1][1] * m[2][2] - m[1][2] * m[2][1])
- m[0][1] * (m[1][0] * m[2][2] - m[1][2] * m[2][0])
+ m[0][2] * (m[1][0] * m[2][1] - m[1][1] * m[2][0])
}
#[allow(dead_code)]
pub(super) fn inv3x3(m: [[f64; 3]; 3]) -> [[f64; 3]; 3] {
let det = det3x3(m);
if det.abs() < 1e-30 {
return [[1.0, 0.0, 0.0], [0.0, 1.0, 0.0], [0.0, 0.0, 1.0]];
}
let inv_det = 1.0 / det;
[
[
(m[1][1] * m[2][2] - m[1][2] * m[2][1]) * inv_det,
(m[0][2] * m[2][1] - m[0][1] * m[2][2]) * inv_det,
(m[0][1] * m[1][2] - m[0][2] * m[1][1]) * inv_det,
],
[
(m[1][2] * m[2][0] - m[1][0] * m[2][2]) * inv_det,
(m[0][0] * m[2][2] - m[0][2] * m[2][0]) * inv_det,
(m[0][2] * m[1][0] - m[0][0] * m[1][2]) * inv_det,
],
[
(m[1][0] * m[2][1] - m[1][1] * m[2][0]) * inv_det,
(m[0][1] * m[2][0] - m[0][0] * m[2][1]) * inv_det,
(m[0][0] * m[1][1] - m[0][1] * m[1][0]) * inv_det,
],
]
}
#[allow(dead_code)]
pub(super) fn inv3x3_transpose(m: [[f64; 3]; 3]) -> [[f64; 3]; 3] {
transpose3x3(inv3x3(m))
}
#[allow(dead_code)]
pub(super) fn transpose3x3(m: [[f64; 3]; 3]) -> [[f64; 3]; 3] {
[
[m[0][0], m[1][0], m[2][0]],
[m[0][1], m[1][1], m[2][1]],
[m[0][2], m[1][2], m[2][2]],
]
}
#[allow(dead_code)]
pub(super) fn mul3x3(a: [[f64; 3]; 3], b: [[f64; 3]; 3]) -> [[f64; 3]; 3] {
let mut c = [[0.0_f64; 3]; 3];
for i in 0..3 {
for j in 0..3 {
for k in 0..3 {
c[i][j] += a[i][k] * b[k][j];
}
}
}
c
}
#[allow(dead_code)]
pub(super) fn isotropic_d_matrix(lambda: f64, mu: f64) -> [[f64; 6]; 6] {
let l2m = lambda + 2.0 * mu;
[
[l2m, lambda, lambda, 0.0, 0.0, 0.0],
[lambda, l2m, lambda, 0.0, 0.0, 0.0],
[lambda, lambda, l2m, 0.0, 0.0, 0.0],
[0.0, 0.0, 0.0, mu, 0.0, 0.0],
[0.0, 0.0, 0.0, 0.0, mu, 0.0],
[0.0, 0.0, 0.0, 0.0, 0.0, mu],
]
}