semiflow 0.15.0-beta

no_std Rust PDE solver: evolution equations via Chernoff approximation of operator semigroups (heat, diffusion, Schrodinger)
Documentation 
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

// Tests for `resolvent_residual.rs` (LaplaceChernoffResolventResidual, ADR-0083).
//
// Properties asserted:
//   1. Construction succeeds (budget stored correctly).
//   2. budget() returns the configured value.
//   3. verify_residual returns Err for lambda <= 0.
//   4. verify_residual returns Err for lambda = NaN.
//   5. verify_residual returns Ok and finite residual for valid lambda.
//   6. fresh_from_fn creates a GridFn1D with correct values.
//   7. sample_at returns Err for empty coordinate slice.
//   8. sample_at returns Ok for a valid interior point.

use crate::{
    diffusion::DiffusionChernoff,
    grid::Grid1D,
    grid_fn::GridFn1D,
    resolvent::{LaplaceChernoffResolvent, LaplaceQuadrature},
    resolvent_residual::{LaplaceChernoffResolventResidual, Sampleable},
};

fn make_residual_wrapper(
    n: usize,
) -> LaplaceChernoffResolventResidual<DiffusionChernoff<f64>, f64> {
    let grid = Grid1D::new(0.0_f64, 1.0, n).unwrap();
    let kernel = DiffusionChernoff::new(|_| 1.0_f64, |_| 0.0, |_| 0.0, 1.0, grid);
    let resolvent = LaplaceChernoffResolvent::new(kernel, 4, LaplaceQuadrature::GaussLaguerre32)
        .expect("resolvent construction failed");
    LaplaceChernoffResolventResidual::new(resolvent, 1e-3_f64)
}

fn make_datum(n: usize) -> GridFn1D<f64> {
    let grid = Grid1D::new(0.0_f64, 1.0, n).unwrap();
    GridFn1D::from_fn(grid, |x| (-(x - 0.5) * (x - 0.5) * 50.0).exp())
}

// ── Construction ──────────────────────────────────────────────────────────────

#[test]
fn construction_succeeds() {
    let wrapper = make_residual_wrapper(32);
    assert!((wrapper.budget() - 1e-3).abs() < 1e-15, "budget mismatch");
}

// ── budget() ─────────────────────────────────────────────────────────────────

#[test]
fn budget_returns_configured_value() {
    let grid = Grid1D::new(0.0_f64, 1.0, 16).unwrap();
    let kernel = DiffusionChernoff::new(|_| 1.0_f64, |_| 0.0, |_| 0.0, 1.0, grid);
    let resolvent =
        LaplaceChernoffResolvent::new(kernel, 4, LaplaceQuadrature::GaussLaguerre32).unwrap();
    let wrapper = LaplaceChernoffResolventResidual::new(resolvent, 5e-4_f64);
    assert!((wrapper.budget() - 5e-4).abs() < 1e-15);
}

// ── verify_residual: Err for lambda <= 0 ─────────────────────────────────────

#[test]
fn verify_residual_rejects_nonpositive_lambda() {
    let wrapper = make_residual_wrapper(32);
    let f = make_datum(32);
    assert!(
        wrapper.verify_residual(0.0_f64, &f).is_err(),
        "expected Err for lambda=0"
    );
    assert!(
        wrapper.verify_residual(-1.0_f64, &f).is_err(),
        "expected Err for lambda=-1"
    );
}

// ── verify_residual: Err for lambda = NaN ────────────────────────────────────

#[test]
fn verify_residual_rejects_nan_lambda() {
    let wrapper = make_residual_wrapper(32);
    let f = make_datum(32);
    assert!(
        wrapper.verify_residual(f64::NAN, &f).is_err(),
        "expected Err for NaN lambda"
    );
}

// ── verify_residual: Ok for valid lambda ─────────────────────────────────────

#[test]
fn verify_residual_ok_for_positive_lambda() {
    let wrapper = make_residual_wrapper(32);
    let f = make_datum(32);
    let result = wrapper.verify_residual(1.0_f64, &f);
    assert!(result.is_ok(), "verify_residual failed: {result:?}");
    let err = result.unwrap();
    assert!(err.is_finite(), "residual is non-finite: {err}");
    assert!(err >= 0.0, "residual is negative: {err}");
}

// ── Sampleable: fresh_from_fn ─────────────────────────────────────────────────

#[test]
fn fresh_from_fn_gives_correct_values() {
    let grid = Grid1D::new(0.0_f64, 1.0, 16).unwrap();
    let template = GridFn1D::from_fn(grid, |_| 0.0_f64);
    let result = template.fresh_from_fn(&|x| x[0] * x[0]).unwrap();
    assert_eq!(result.values.len(), 16);
    for (i, &v) in result.values.iter().enumerate() {
        let x = grid.x_at(i);
        assert!((v - x * x).abs() < 1e-15, "at i={i}: expected {}, got {v}", x * x);
    }
}

// ── Sampleable: sample_at empty coords ───────────────────────────────────────

#[test]
fn sample_at_empty_coords_returns_err() {
    let grid = Grid1D::new(0.0_f64, 1.0, 16).unwrap();
    let gf = GridFn1D::from_fn(grid, f64::sin);
    let result = gf.sample_at(&[]);
    assert!(result.is_err(), "expected Err for empty coordinate slice");
}

// ── Sampleable: sample_at interior point ─────────────────────────────────────

#[test]
fn sample_at_valid_point_returns_ok() {
    let grid = Grid1D::new(0.0_f64, 1.0, 16).unwrap();
    let gf = GridFn1D::from_fn(grid, |x| x * 2.0);
    let result = gf.sample_at(&[0.5_f64]);
    assert!(result.is_ok(), "sample_at failed: {result:?}");
    let v = result.unwrap();
    assert!(v.is_finite(), "sample_at returned non-finite: {v}");
}