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
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
//! Property and layout tests for the two-electron (ERI) builder.
//!
//! Pure-Rust, no external library. These tests assert the
//! permutational symmetry and the documented block layout/strides — the latter
//! matters now that ERI blocks are *not* index-symmetric the way the Phase-1
//! operators were.
use integral::{Basis, Shell};
use integral_math::am::n_cart;
/// Four shells of different angular momenta on different centers, so every
/// (ij|kl) block is genuinely distinct (no accidental symmetry to hide a
/// transposition bug).
fn mixed_basis() -> Basis {
Basis::new(vec![
Shell::new(0, [0.0, 0.0, 0.0], vec![1.4, 0.5], vec![0.6, 0.5]).unwrap(), // s
Shell::new(1, [0.6, -0.3, 0.2], vec![0.9], vec![1.0]).unwrap(), // p
Shell::new(2, [-0.4, 0.7, -0.1], vec![1.1], vec![1.0]).unwrap(), // d
Shell::new(3, [0.2, 0.5, 0.8], vec![0.7], vec![1.0]).unwrap(), // f
])
}
fn nao(b: &Basis) -> usize {
b.nao_cart()
}
#[test]
fn ssss_is_positive_and_finite() {
// A single s shell: (00|00) must be a finite positive number.
let basis = Basis::new(vec![
Shell::new(0, [0.0, 0.0, 0.0], vec![0.8], vec![1.0]).unwrap()
]);
let eri = basis.eri();
assert_eq!(eri.len(), 1);
assert!(eri[0].is_finite() && eri[0] > 0.0, "(ss|ss) = {}", eri[0]);
}
#[test]
fn eight_fold_permutational_symmetry() {
// Genuine assertion: the tensor is built block-by-block with no
// symmetrization, so equality of the 8 permutations is a real check.
let basis = mixed_basis();
let n = nao(&basis);
let eri = basis.eri();
let at = |i: usize, j: usize, k: usize, l: usize| eri[((i * n + j) * n + k) * n + l];
// Tight *relative* residual on significant elements (|base| ≥ 1e-3·peak) plus an
// absolute floor for the rest (the `eri_maxima` convention). The two permutation
// paths of a structurally-near-zero element differ only by f64 round-off
// (~1e-16 abs); a floorless relative residual would over-weight those and reject
// any value-changing kernel optimisation, while the consumer reuses one value per
// unique quartet and never observes the difference.
let peak = eri.iter().fold(0.0f64, |m, &x| m.max(x.abs()));
let floor = 1e-3 * peak;
let mut worst_sig = 0.0_f64;
let mut worst_abs = 0.0_f64;
for i in 0..n {
for j in 0..n {
for k in 0..n {
for l in 0..n {
let base = at(i, j, k, l);
let perms = [
at(j, i, k, l),
at(i, j, l, k),
at(j, i, l, k),
at(k, l, i, j),
at(l, k, i, j),
at(k, l, j, i),
at(l, k, j, i),
];
for &p in &perms {
let d = (base - p).abs();
worst_abs = worst_abs.max(d);
if base.abs() >= floor {
worst_sig = worst_sig.max(d / base.abs());
}
}
}
}
}
}
assert!(
worst_sig < 1e-11,
"worst symmetry relative residual {worst_sig:e}"
);
assert!(
worst_abs < 1e-11 * peak.max(1.0) + 1e-12,
"worst symmetry absolute residual {worst_abs:e} (peak {peak:e})"
);
}
#[test]
fn block_layout_matches_dense_strides() {
// `eri_block(i,j,k,l)` must equal the dense tensor slice at the four shell
// offsets under the documented row-major (a,b,c,d) layout. Use all-distinct,
// different-L shells so a stride/transposition bug cannot stay hidden.
let basis = mixed_basis();
let n = nao(&basis);
let dense = basis.eri();
let shells = basis.shells();
let offs: Vec<usize> = {
let mut acc = 0;
let mut v = Vec::new();
for s in shells {
v.push(acc);
acc += s.n_cart();
}
v
};
for (qi, qj, qk, ql) in [(1usize, 0usize, 2usize, 3usize), (0, 1, 3, 2), (2, 3, 1, 0)] {
let block = basis.eri_block(qi, qj, qk, ql);
let (nb, nc, nd) = (
n_cart(shells[qj].l()),
n_cart(shells[qk].l()),
n_cart(shells[ql].l()),
);
let (na2,) = (n_cart(shells[qi].l()),);
assert_eq!(block.len(), na2 * nb * nc * nd);
for a in 0..na2 {
for b in 0..nb {
for c in 0..nc {
for d in 0..nd {
let src = ((a * nb + b) * nc + c) * nd + d;
let (i, j, k, l) = (offs[qi] + a, offs[qj] + b, offs[qk] + c, offs[ql] + d);
let dval = dense[((i * n + j) * n + k) * n + l];
assert!(
(block[src] - dval).abs() <= 1e-12 * dval.abs().max(1e-12),
"layout mismatch at block {src} ({i}{j}|{k}{l}): {} vs {}",
block[src],
dval
);
}
}
}
}
}
}