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
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
//! #2339 — the ARD PSD majorizer's hard clamp `max(V'', 0)` is replaced by the
//! smooth homogeneity-preserving softplus envelope `α·softplus_{τ₀}(cos κt)`, so
//! the streaming `½log|B̃|` criterion becomes a composite-analytic estimand (no
//! kink at the clamp seam). These gates pin the four load-bearing properties:
//!
//! 1. SEAM IDENTITY — `psd_majorizer_hess + negative_hessian_remainder == V''`
//! (the majorizer + its complementary concave remainder reconstruct the exact
//! signed prior curvature), and the Euclidean path is a bit-identical no-op.
//! 2. HOMOGENEITY — the majorizer is exactly degree-one in `α = e^{ρ_ard}`, which
//! the `½log|B|` explicit-ρ θ-adjoint traces rely on (`∂/∂ρ = value`).
//! 3. POSITIVITY / FLOOR — `B ≻ 0` everywhere, and the deviation from the hard
//! clamp stays at/below the criterion's spectral-deflation floor `α·1e-8`.
//! 4. C¹ AT THE SEAM — the majorizer's coordinate derivative is the continuous
//! `α·logistic(cos/τ₀)·(−κ sin κt)` (the hard clamp's `1{cos>0}` step becomes
//! the logistic passing smoothly through `½` at the seam), matched by finite
//! differences of the value.
use crate::manifold::ArdAxisPrior;
use gam_solve::arrow_schur::SPECTRAL_DEFLATION_REL_FLOOR;
use std::f64::consts::{LN_2, TAU};
/// τ₀, the dimensionless softplus temperature. Mirrors the derived constant so
/// the test fails loudly if the derivation constant is ever silently retuned.
fn tau0() -> f64 {
ArdAxisPrior::CLAMP_TEMPERATURE
}
/// The maximum admissible deviation of the smooth clamp from the hard clamp, in
/// curvature units: `α·τ₀·ln2 = α·SPECTRAL_DEFLATION_REL_FLOOR` (the deflation
/// floor relative to the axis's own curvature scale `α`).
fn deviation_ceiling(alpha: f64) -> f64 {
alpha * tau0() * LN_2
}
#[test]
fn tau0_is_the_derived_deflation_floor_constant() {
// τ₀ = SPECTRAL_DEFLATION_REL_FLOOR / ln2, and hence α·τ₀·ln2 = α·floor.
assert_eq!(tau0(), SPECTRAL_DEFLATION_REL_FLOOR / LN_2);
assert!(
(deviation_ceiling(1.0) - SPECTRAL_DEFLATION_REL_FLOOR).abs()
<= 4.0 * f64::EPSILON * SPECTRAL_DEFLATION_REL_FLOOR,
"the per-axis deviation ceiling must equal the relative deflation floor"
);
}
#[test]
fn euclidean_clamp_is_a_bit_identical_no_op() {
// Non-periodic axes have constant curvature `V'' = α > 0`; the clamp must be a
// literal no-op (bit-for-bit), and the remainder exactly zero.
for &alpha in &[1e-3, 0.5, 1.0, 7.0, 250.0] {
for &t in &[-3.0, -0.25, 0.0, 0.8, 5.5] {
let prior = ArdAxisPrior::eval(alpha, t, None);
assert_eq!(prior.psd_majorizer_hess(), alpha);
assert_eq!(prior.negative_hessian_remainder(), 0.0);
assert_eq!(
prior.psd_majorizer_hess() + prior.negative_hessian_remainder(),
prior.hess
);
}
}
}
#[test]
fn seam_identity_reconstructs_exact_hess() {
// `B + R == V''`. Outside the ~1e-8-wide transition band one term underflows to
// exactly zero (softplus collapses to the hard clamp), so the reconstruction is
// bit-for-bit exact; the generic offset grid below never lands in that band.
let period = 1.0;
for &alpha in &[1e-2, 1.0, 13.0, 400.0] {
for i in 0..257 {
let t = 0.001_3 + i as f64 * (0.997 / 257.0);
let prior = ArdAxisPrior::eval(alpha, t, Some(period));
let sum = prior.psd_majorizer_hess() + prior.negative_hessian_remainder();
assert_eq!(
sum, prior.hess,
"seam identity must reconstruct V'' exactly at alpha={alpha}, t={t}"
);
}
}
}
#[test]
fn majorizer_is_positive_semidefinite() {
// `B ⪰ 0` everywhere (the majorizer never contributes negative curvature),
// matching the hard clamp it replaces. Mathematically softplus > 0, but on the
// deep concave half `e^{cos/τ₀}` underflows so the returned value is exactly 0
// (bit-identical to `max(cos,0)` there) — hence PSD, not strictly PD.
let period = 1.0;
let kappa = TAU / period;
for &alpha in &[1e-6, 1e-2, 1.0, 50.0, 1e4] {
for i in 0..512 {
let t = i as f64 * (period / 512.0);
let cos = (kappa * t).cos();
let b = ArdAxisPrior::eval(alpha, t, Some(period)).psd_majorizer_hess();
assert!(b >= 0.0, "majorizer must be non-negative; got {b} at t={t}");
// On the convex half the majorizer is strictly positive (it is at least
// the hard clamp `α·cos > 0`).
if cos > 1e-4 {
assert!(b > 0.0, "convex-half majorizer must be positive; got {b} at t={t}");
}
}
}
// Within the transition band around the seam the smoothing lifts the majorizer
// strictly above zero even where the exact curvature is zero/negative.
for &alpha in &[1.0, 100.0] {
assert!(ArdAxisPrior::smooth_clamp(alpha, 0.0) > 0.0);
assert!(ArdAxisPrior::smooth_clamp(alpha, -tau0()) > 0.0);
}
}
#[test]
fn concave_remainder_is_nonpositive() {
// `R = V'' − B ≤ 0` (since B = α·softplus(cos) ≥ α·max(cos,0) ≥ α cos = V''),
// i.e. the restored `E = −R ⪰ 0` — the exact `A = B − E` split is preserved.
let period = 1.0;
for &alpha in &[1e-2, 1.0, 30.0] {
for i in 0..512 {
let t = i as f64 * (period / 512.0);
let r = ArdAxisPrior::eval(alpha, t, Some(period)).negative_hessian_remainder();
assert!(r <= 0.0, "remainder must be non-positive; got {r} at t={t}");
}
}
}
#[test]
fn deviation_from_hard_clamp_within_floor() {
// |B − α·max(cos,0)| ≤ α·τ₀·ln2 everywhere, with equality (the kink is filled)
// at the seam cos κt = 0.
let period = 1.0;
let kappa = TAU / period;
for &alpha in &[1e-2, 1.0, 25.0, 1e3] {
let ceiling = deviation_ceiling(alpha);
let mut max_dev = 0.0_f64;
for i in 0..2048 {
let t = i as f64 * (period / 2048.0);
let cos = (kappa * t).cos();
let hard = alpha * cos.max(0.0);
let smooth = ArdAxisPrior::eval(alpha, t, Some(period)).psd_majorizer_hess();
let dev = (smooth - hard).abs();
assert!(
dev <= ceiling + 8.0 * f64::EPSILON * alpha,
"deviation {dev} exceeds floor ceiling {ceiling} at alpha={alpha}, t={t}"
);
max_dev = max_dev.max(dev);
}
// At the exact seam cos = 0 the kink is filled to the ceiling: the hard
// clamp gives 0, the smooth clamp gives α·softplus(0) = α·τ₀·ln2.
let seam = ArdAxisPrior::smooth_clamp(alpha, 0.0);
assert!(
(seam - ceiling).abs() <= 1e-12 * ceiling.max(f64::MIN_POSITIVE),
"seam value {seam} must equal the deviation ceiling {ceiling}"
);
assert!(max_dev > 0.5 * ceiling, "the smoothing must actually engage");
}
}
#[test]
fn homogeneity_degree_one_in_alpha() {
// ∂/∂ρ_ard[α·softplus(cos)] = α·softplus(cos) requires exact degree-one
// homogeneity in α (the explicit-ρ log-det traces depend on it). Doubling α is
// an exact power-of-two scaling, so the majorizer must double bit-for-bit.
let period = 1.0;
for &alpha in &[1e-3, 0.7, 3.0, 90.0] {
for i in 0..129 {
let t = 0.002 + i as f64 * (0.996 / 129.0);
let b1 = ArdAxisPrior::eval(alpha, t, Some(period)).psd_majorizer_hess();
let b2 = ArdAxisPrior::eval(2.0 * alpha, t, Some(period)).psd_majorizer_hess();
assert_eq!(2.0 * b1, b2, "majorizer must be exactly degree-1 in alpha");
// A general (non power-of-two) scale to tight relative tolerance.
let lambda = 3.5;
let bl = ArdAxisPrior::eval(lambda * alpha, t, Some(period)).psd_majorizer_hess();
assert!(
(bl - lambda * b1).abs() <= 8.0 * f64::EPSILON * bl.abs().max(f64::MIN_POSITIVE),
"degree-1 homogeneity must hold for a general scale"
);
}
}
}
#[test]
fn slope_is_continuous_through_the_seam() {
// The hard clamp's derivative is a step `1{cos>0}`; the smooth clamp's is the
// logistic `s'(c) = logistic(c/τ₀)`, which passes continuously through ½ at the
// seam (no kink). Monotone, symmetric about (0, ½), saturating to {0,1}.
assert_eq!(ArdAxisPrior::clamp_slope(0.0), 0.5);
assert_eq!(ArdAxisPrior::clamp_slope(1.0), 1.0); // saturates within f64
assert_eq!(ArdAxisPrior::clamp_slope(-1.0), 0.0);
let mut prev = 0.0_f64;
for i in 0..401 {
let c = -1.0 + i as f64 * (2.0 / 400.0);
let s = ArdAxisPrior::clamp_slope(c);
assert!((0.0..=1.0).contains(&s));
assert!(s >= prev - 1e-15, "slope must be non-decreasing in cos");
prev = s;
// Symmetry logistic(−z) = 1 − logistic(z).
let s_neg = ArdAxisPrior::clamp_slope(-c);
assert!((s + s_neg - 1.0).abs() <= 1e-12);
}
}
#[test]
fn value_slope_consistency_by_finite_difference() {
// The analytic coordinate derivative `α·s'(cos)·(−κ sin κt)` (the closed form
// `ard_majorized_hessian_derivative` assembles) must equal a central finite
// difference of the majorizer value. Away from the seam the function is the
// (smooth) hard clamp; at the seam a τ₀-scaled step resolves the transition and
// recovers the ½-slope, demonstrating C¹.
let period = 1.0;
let kappa = TAU / period;
let alpha = 4.0;
// (a) Away from the seam: h = 1e-6 is fine (function ≈ α·max(cos,0), smooth).
for &cos_target in &[0.9_f64, 0.5, -0.5, -0.9] {
// pick t with cos(κt) = cos_target on the falling branch (sin > 0).
let t = cos_target.acos() / kappa;
let h = 1e-6;
let f = |tt: f64| ArdAxisPrior::eval(alpha, tt, Some(period)).psd_majorizer_hess();
let fd = (f(t + h) - f(t - h)) / (2.0 * h);
let cos = (kappa * t).cos();
let sin = (kappa * t).sin();
let analytic = -alpha * kappa * sin * ArdAxisPrior::clamp_slope(cos);
assert!(
(fd - analytic).abs() <= 1e-4 * (1.0 + analytic.abs()),
"value/slope FD mismatch away from seam at cos={cos_target}: fd={fd}, analytic={analytic}"
);
}
// (b) At the seam, resolve the τ₀-wide band via smooth_clamp directly. Central
// FD in the dimensionless cosine with a τ₀-scaled step recovers exactly the
// logistic midpoint slope ½ — the kink is filled, not a step.
for &k in &[1.0_f64, 2.0, 5.0] {
let h = k * tau0();
let fd = (ArdAxisPrior::smooth_clamp(alpha, h) - ArdAxisPrior::smooth_clamp(alpha, -h))
/ (2.0 * h);
let analytic = alpha * ArdAxisPrior::clamp_slope(0.0); // = 0.5·α
assert!(
(fd - analytic).abs() <= 1e-9 * alpha,
"seam FD slope {fd} must equal the logistic midpoint {analytic}"
);
}
}