Rustb 0.7.0

A package for calculating band, angle state, linear and nonlinear conductivities based on tight-binding models
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83

//! Symmetric simplex quadrature rules and barycentric interpolation.

use ndarray::prelude::*;
use num_complex::Complex;

// ── Quadrature rules ────────────────────────────────────────────────────

/// Degree‑2 symmetric 3‑point rule for the reference triangle.
/// Barycentric coordinates; weights sum to 1.
pub(crate) const TRI_QUAD_PTS_3: [[f64; 3]; 3] = [
    [2.0 / 3.0, 1.0 / 6.0, 1.0 / 6.0],
    [1.0 / 6.0, 2.0 / 3.0, 1.0 / 6.0],
    [1.0 / 6.0, 1.0 / 6.0, 2.0 / 3.0],
];
pub(crate) const TRI_QUAD_WTS_3: [f64; 3] = [1.0 / 3.0, 1.0 / 3.0, 1.0 / 3.0];

/// Degree‑2 symmetric 4‑point rule for the reference tetrahedron.
/// Barycentric coordinates; weights sum to 1.
pub(crate) const TET_QUAD_PTS_4: [[f64; 4]; 4] = [
    [
        0.5854101966249685,
        0.1381966011250105,
        0.1381966011250105,
        0.1381966011250105,
    ],
    [
        0.1381966011250105,
        0.5854101966249685,
        0.1381966011250105,
        0.1381966011250105,
    ],
    [
        0.1381966011250105,
        0.1381966011250105,
        0.5854101966249685,
        0.1381966011250105,
    ],
    [
        0.1381966011250105,
        0.1381966011250105,
        0.1381966011250105,
        0.5854101966249685,
    ],
];
pub(crate) const TET_QUAD_WTS_4: [f64; 4] = [0.25, 0.25, 0.25, 0.25];

// ── Barycentric interpolation ───────────────────────────────────────────

/// Interpolate band energies at barycentric coords `lam`.
pub(crate) fn bary_interp_band(bands: &[Vec<f64>], lam: &[f64], nsta: usize) -> Vec<f64> {
    let mut out = vec![0.0; nsta];
    for v in 0..bands.len() {
        let lv = lam[v];
        if lv == 0.0 {
            continue;
        }
        for n in 0..nsta {
            out[n] += bands[v][n] * lv;
        }
    }
    out
}

/// Interpolate a matrix field `mats[d+1]` each `(nsta, nsta)` at
/// barycentric coords `lam`.
pub(crate) fn bary_interp_matrix(
    mats: &[Array2<Complex<f64>>],
    lam: &[f64],
) -> Array2<Complex<f64>> {
    let n = mats[0].nrows();
    let mut out = Array2::<Complex<f64>>::zeros((n, n));
    for (mat, &w) in mats.iter().zip(lam.iter()) {
        if w == 0.0 {
            continue;
        }
        for i in 0..n {
            for j in 0..n {
                out[[i, j]] += mat[[i, j]] * w;
            }
        }
    }
    out
}