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}