Skip to main content

geop_core_math/solvers/
refine.rs

1//! Adaptive refinement: something built at a resolution — samples along a
2//! curve, cells across a face, stations along a blend — built again finer
3//! until what it measures of itself is good enough.
4//!
5//! The resolution is a free choice; whether the result at one is good
6//! enough is the caller's measurement, made in `step`. The solver promises
7//! that every step is taken from the resolution the one before asked for,
8//! at most `most` of them, and says which way it ended: with the result a
9//! step accepted ([`Refined::Done`]), or out of steps with the resolution
10//! the last one asked for, never accepted ([`Refined::Exhausted`]). A caller
11//! that has to stop at a coarse result says so by building it from the
12//! exhausted resolution; one that must not, errs.
13
14use crate::solvers::repeat::{Repeated, Step, repeat};
15
16/// What one step made of its resolution.
17pub enum Refine<X, T> {
18    /// Good enough: the result.
19    Done(T),
20    /// Not yet: build again at this resolution.
21    Finer(X),
22}
23
24/// How a refinement ended.
25#[derive(Debug, PartialEq, Eq)]
26pub enum Refined<X, T> {
27    /// A step accepted this.
28    Done(T),
29    /// Every step was taken; the last asked for this resolution.
30    Exhausted(X),
31}
32
33impl<X> Refined<X, X> {
34    /// What it ended with, accepted or not.
35    pub fn either(self) -> X {
36        match self {
37            Refined::Done(x) | Refined::Exhausted(x) => x,
38        }
39    }
40}
41
42/// Takes `step` from `start`, then from each resolution it asks for, at
43/// most `most` times, until one is [`Refine::Done`]. `step` may keep state
44/// of its own (a measure of the step before, to tell whether refining
45/// still helps), since every step is taken once, in order.
46pub fn refine<X, T, E>(
47    start: X,
48    most: usize,
49    mut step: impl FnMut(X) -> Result<Refine<X, T>, E>,
50) -> Result<Refined<X, T>, E> {
51    // The resolution the next step is taken from: there between steps.
52    let mut x = Some(start);
53    let run = repeat(most, || {
54        let at = x.take().expect("each step leaves the next resolution");
55        Ok(match step(at)? {
56            Refine::Done(done) => Step::Done(done),
57            Refine::Finer(finer) => {
58                x = Some(finer);
59                Step::Again
60            }
61        })
62    })?;
63    Ok(match run {
64        Repeated::Done { value, .. } => Refined::Done(value),
65        Repeated::Exhausted { .. } => {
66            Refined::Exhausted(x.expect("a step not done leaves the next resolution"))
67        }
68    })
69}
70
71#[cfg(test)]
72mod tests {
73    use super::*;
74
75    /// Midpoint sums of `x²` over `[0, 1]`, the intervals doubled until
76    /// two successive sums agree to `1e-6`: the sum it accepts is that
77    /// close to a third.
78    #[test]
79    fn doubling_until_successive_sums_agree() {
80        let sum = |n: usize| {
81            (0..n)
82                .map(|i| ((i as f64 + 0.5) / n as f64).powi(2) / n as f64)
83                .sum::<f64>()
84        };
85        let out = refine::<_, _, ()>((1usize, f64::INFINITY), 40, |(n, before)| {
86            let s = sum(n);
87            Ok(if (s - before).abs() < 1e-6 {
88                Refine::Done(s)
89            } else {
90                Refine::Finer((2 * n, s))
91            })
92        })
93        .unwrap();
94        let Refined::Done(s) = out else {
95            panic!("{out:?}")
96        };
97        assert!((s - 1.0 / 3.0).abs() < 1e-6);
98    }
99
100    /// A refinement that is never good enough reports that it ran out,
101    /// with the resolution the last step asked for — not a result.
102    #[test]
103    fn running_out_is_reported() {
104        let out = refine::<usize, (), ()>(1, 5, |n| Ok(Refine::Finer(2 * n))).unwrap();
105        assert_eq!(out, Refined::Exhausted(32));
106    }
107}