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
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
//! Chebyshev-Lobatto node and barycentric weight tables (ADR-0090, math.md §9.2.7).
//!
//! Provides `chebyshev_nodes` and `chebyshev_weights` for M ∈ {8, 16, 32, 64, 128, 256, 512}.
//!
//! Algorithm (NORMATIVE, Berrut-Trefethen 2004 §2):
//! Nodes: `x_k = cos(k·π / M)`, k = 0..M (Chebyshev-Lobatto on [−1, 1])
//! Weights: `w_k = (−1)^k · δ_k`, `δ_0` = `δ_M` = 1/2, `δ_k` = 1.0 else
//!
//! Generated by `scripts/generate_chebyshev_nodes.py`.
//! DO NOT hand-edit — re-run the script to regenerate.
//!
//! References:
//! Berrut & Trefethen 2004, SIAM Review 46:501 — barycentric formula (5.1).
//! Trefethen 2000 Spectral Methods in MATLAB — Table A.5 differentiation matrix.
// Grid/index/count values (usize) cast to f64 for coordinate and coefficient computations;
// all values are grid sizes or step counts ≪ 2^52, so precision loss is impossible in practice.
#![allow(clippy::cast_precision_loss)]
use crate::error::SemiflowError;
/// Supported M values for Chebyshev-Lobatto virtual node sets (ADR-0090 §"Supported M").
///
/// Used for documentation and `is_supported_m` validation.
pub(crate) const SUPPORTED_M: [usize; 7] = [8, 16, 32, 64, 128, 256, 512];
/// Returns `true` if M is in the supported set {8, 16, 32, 64, 128, 256, 512}.
#[inline]
pub(crate) fn is_supported_m(m: usize) -> bool {
SUPPORTED_M.contains(&m)
}
/// Return a static slice of M+1 Chebyshev-Lobatto nodes on [−1, 1] for the given M.
///
/// # Errors
/// Returns [`SemiflowError::DomainViolation`] if M ∉ {8, 16, 32, 64, 128, 256, 512}.
pub(crate) fn chebyshev_nodes(m: usize) -> Result<&'static [f64], SemiflowError> {
match m {
8 => Ok(&CHEB_NODES_8),
16 => Ok(&CHEB_NODES_16),
32 => Ok(&CHEB_NODES_32),
64 => Ok(&CHEB_NODES_64),
128 => Ok(&CHEB_NODES_128),
256 => Ok(&CHEB_NODES_256),
512 => Ok(&CHEB_NODES_512),
_ => Err(SemiflowError::DomainViolation {
what: "Chebyshev M must be in {8,16,32,64,128,256,512}",
value: m as f64,
}),
}
}
/// Return a static slice of M+1 barycentric weights for the given M.
///
/// # Errors
/// Returns [`SemiflowError::DomainViolation`] if M ∉ {8, 16, 32, 64, 128, 256, 512}.
pub(crate) fn chebyshev_weights(m: usize) -> Result<&'static [f64], SemiflowError> {
match m {
8 => Ok(&CHEB_WEIGHTS_8),
16 => Ok(&CHEB_WEIGHTS_16),
32 => Ok(&CHEB_WEIGHTS_32),
64 => Ok(&CHEB_WEIGHTS_64),
128 => Ok(&CHEB_WEIGHTS_128),
256 => Ok(&CHEB_WEIGHTS_256),
512 => Ok(&CHEB_WEIGHTS_512),
_ => Err(SemiflowError::DomainViolation {
what: "Chebyshev M must be in {8,16,32,64,128,256,512}",
value: m as f64,
}),
}
}
// ---------------------------------------------------------------------------
// AUTO-GENERATED node and weight tables
// Generated by: python3 scripts/generate_chebyshev_nodes.py
// Algorithm: x_k = cos(k*pi/M), w_k = (-1)^k * delta_k (delta_0=delta_M=0.5)
// ---------------------------------------------------------------------------
// cos(kπ/8) table: 7.071...e-1 = cos(π/4) = FRAC_1_SQRT_2, kept as literal for table consistency
#[allow(clippy::approx_constant)]
pub(crate) const CHEB_NODES_8: [f64; 9] = [
1.0,
9.238_795_325_112_867e-1,
7.071_067_811_865_476e-1,
3.826_834_323_650_898_4e-1,
6.123_233_995_736_766e-17,
-3.826_834_323_650_897e-1,
-7.071_067_811_865_475e-1,
-9.238_795_325_112_867e-1,
-1.0,
];
pub(crate) const CHEB_WEIGHTS_8: [f64; 9] = [0.5, -1.0, 1.0, -1.0, 1.0, -1.0, 1.0, -1.0, 0.5];
// cos(kπ/16) table: 7.071...e-1 = FRAC_1_SQRT_2, kept as literal for table consistency
#[allow(clippy::approx_constant)]
pub(crate) const CHEB_NODES_16: [f64; 17] = [
1.0,
9.807_852_804_032_304e-1,
9.238_795_325_112_867e-1,
8.314_696_123_025_452e-1,
7.071_067_811_865_476e-1,
5.555_702_330_196_022e-1,
3.826_834_323_650_898_4e-1,
1.950_903_220_161_283_6e-1,
6.123_233_995_736_766e-17,
-1.950_903_220_161_281_4e-1,
-3.826_834_323_650_897e-1,
-5.555_702_330_196_022e-1,
-7.071_067_811_865_475e-1,
-8.314_696_123_025_451e-1,
-9.238_795_325_112_867e-1,
-9.807_852_804_032_304e-1,
-1.0,
];
pub(crate) const CHEB_WEIGHTS_16: [f64; 17] = [
0.5, -1.0, 1.0, -1.0, 1.0, -1.0, 1.0, -1.0, 1.0, -1.0, 1.0, -1.0, 1.0, -1.0, 1.0, -1.0, 0.5,
];
// cos(kπ/32) table: 7.071...e-1 = FRAC_1_SQRT_2, kept as literal for table consistency
#[allow(clippy::approx_constant)]
pub(crate) const CHEB_NODES_32: [f64; 33] = [
1.0,
9.951_847_266_721_969e-1,
9.807_852_804_032_304e-1,
9.569_403_357_322_089e-1,
9.238_795_325_112_867e-1,
8.819_212_643_483_55e-1,
8.314_696_123_025_452e-1,
7.730_104_533_627_37e-1,
7.071_067_811_865_476e-1,
6.343_932_841_636_455e-1,
5.555_702_330_196_022e-1,
4.713_967_368_259_976_4e-1,
3.826_834_323_650_898_4e-1,
2.902_846_772_544_624e-1,
1.950_903_220_161_283_6e-1,
9.801_714_032_956_06e-2,
6.123_233_995_736_766e-17,
-9.801_714_032_956_045e-2,
-1.950_903_220_161_281_4e-1,
-2.902_846_772_544_623e-1,
-3.826_834_323_650_897e-1,
-4.713_967_368_259_976e-1,
-5.555_702_330_196_022e-1,
-6.343_932_841_636_455e-1,
-7.071_067_811_865_475e-1,
-7.730_104_533_627_37e-1,
-8.314_696_123_025_451e-1,
-8.819_212_643_483_55e-1,
-9.238_795_325_112_867e-1,
-9.569_403_357_322_089e-1,
-9.807_852_804_032_304e-1,
-9.951_847_266_721_969e-1,
-1.0,
];
pub(crate) const CHEB_WEIGHTS_32: [f64; 33] = [
0.5, -1.0, 1.0, -1.0, 1.0, -1.0, 1.0, -1.0, 1.0, -1.0, 1.0, -1.0, 1.0, -1.0, 1.0, -1.0, 1.0,
-1.0, 1.0, -1.0, 1.0, -1.0, 1.0, -1.0, 1.0, -1.0, 1.0, -1.0, 1.0, -1.0, 1.0, -1.0, 0.5,
];
// For M=64,128,256,512 the node tables are generated by scripts/generate_chebyshev_nodes.py.
// Stored as compact f64 literals; values are cos(k*pi/M) rounded to 18 significant digits.
include!(concat!(
env!("CARGO_MANIFEST_DIR"),
"/src/grid_chebyshev_nodes_large.rs"
));