Mercurial > repos > pointsource_algs / file comparison
comparison: src/forward_model.rs
src/forward_model.rs
- changeset 72
- e9a460a0e638
- parent 63
- 7a8a55fd41c0
equal
deleted
inserted
replaced
| 71:e2953ffd4e0b | 72:e9a460a0e638 |
|---|---|
| 41 pub enum BoundedCurvatureGuess { | 41 pub enum BoundedCurvatureGuess { |
| 42 /// No iterate $μ^k$ is worse than $μ=0$. | 42 /// No iterate $μ^k$ is worse than $μ=0$. |
| 43 BetterThanZero, | 43 BetterThanZero, |
| 44 } | 44 } |
| 45 | 45 |
| 46 /// Curvature error control: helper bounds for (4.2d), (5.2a), (5.2b), (5.15a), and (5.16a). | 46 /// Curvature error control. |
| 47 /// | |
| 48 /// Based on Lemma 5.11 and Example 5.12, the helper bound for (5.15a) and (5.16a) is (3.8). | |
| 49 /// Thus, subject to `guess` being correct, returns factor $(ℓ_F, Θ²)$ such that | |
| 50 /// $B_{F'(μ)} dγ ≤ ℓ_F c_2$ and $⟨F'(μ+Δ)-F'(μ)|Δ⟩ ≤ Θ²|γ|(c_2)$, where $Δ=(π_♯^1-π_♯^0)γ$. | |
| 51 /// The latter is the firm transport Lipshitz property of (3.8b) and Lemma 5.11. | |
| 52 /// | |
| 53 /// This trait is supposed to be implemented by the data term $F$, in the basic case a | |
| 54 /// [`Mapping`] from [`RNDM<N, F>`] to a [`Float`] `F`. | |
| 55 /// The generic implementation for operators that satisfy [`BasicCurvatureBoundEstimates`] | |
| 56 /// uses Remark 5.15 and Example 5.16 for (4.2d) and (5.2a), and (5.2b); | |
| 57 /// and Lemma 3.8 for (3.8). | |
| 58 pub trait BoundedCurvature<F: Float = f64> { | 47 pub trait BoundedCurvature<F: Float = f64> { |
| 59 /// Returns $(ℓ_F, Θ²)$ or individual errors for each. | 48 /// Returns an estimate of $(ℓ_F, ℓ_∇v, Θ²)$ or individual errors for each. |
| 60 fn curvature_bound_components( | 49 fn curvature_bound_components( |
| 61 &self, | 50 &self, |
| 62 guess: BoundedCurvatureGuess, | 51 guess: BoundedCurvatureGuess, |
| 63 ) -> (DynResult<F>, DynResult<F>); | 52 ) -> (DynResult<F>, DynResult<F>, DynResult<F>); |
| 64 } | 53 } |
| 65 | 54 |
| 66 /// Curvature error control: helper bounds for (4.2d), (5.2a), (5.2b), (5.15a), and (5.16a) | 55 /// Curvature error control: helper Lipschitz-like bounds for the quadratic dataterms |
| 67 /// for quadratic dataterms $F(μ) = \frac{1}{2}\|Aμ-b\|^2$. | 56 /// $F(μ) = \frac{1}{2}\|Aμ-b\|^2$. |
| 68 /// | |
| 69 /// This trait is to be implemented by the [`Linear`] operator $A$, in the basic from | |
| 70 /// [`RNDM<N, F>`] to a an Euclidean space. | |
| 71 /// It is used by implementations of [`BoundedCurvature`] for $F$. | |
| 72 /// | |
| 73 /// Based on Lemma 5.11 and Example 5.12, the helper bound for (5.15a) and (5.16a) is (3.8). | |
| 74 /// This trait provides the factor $θ²$ of (3.8) as determined by Lemma 3.8. | |
| 75 /// To aid in calculating (4.2d), (5.2a), (5.2b), motivated by Example 5.16, it also provides | |
| 76 /// $ℓ_F^0$ such that $∇v^k$ $ℓ_F^0 \|Aμ-b\|$-Lipschitz. Here $v^k := F'(∪^k)$. | |
| 77 pub trait BasicCurvatureBoundEstimates<F: Float = f64> { | 57 pub trait BasicCurvatureBoundEstimates<F: Float = f64> { |
| 78 /// Returns $(ℓ_F^0, Θ²)$ or individual errors for each. | 58 /// Returns $(ℓ_F, ℓ_{∇v}^0, Θ²)$ or individual errors for each. |
| 79 fn basic_curvature_bound_components(&self) -> (DynResult<F>, DynResult<F>); | 59 fn basic_curvature_bound_components(&self) -> (DynResult<F>, DynResult<F>, DynResult<F>); |
| 80 } | 60 } |
| 81 | 61 |
| 82 impl<F, A, Z, const N: usize> BoundedCurvature<F> for QuadraticDataTerm<F, RNDM<N, F>, A> | 62 impl<F, A, Z, const N: usize> BoundedCurvature<F> for QuadraticDataTerm<F, RNDM<N, F>, A> |
| 83 where | 63 where |
| 84 F: Float, | 64 F: Float, |
| 87 A: BasicCurvatureBoundEstimates<F>, | 67 A: BasicCurvatureBoundEstimates<F>, |
| 88 { | 68 { |
| 89 fn curvature_bound_components( | 69 fn curvature_bound_components( |
| 90 &self, | 70 &self, |
| 91 guess: BoundedCurvatureGuess, | 71 guess: BoundedCurvatureGuess, |
| 92 ) -> (DynResult<F>, DynResult<F>) { | 72 ) -> (DynResult<F>, DynResult<F>, DynResult<F>) { |
| 93 match guess { | 73 match guess { |
| 94 BoundedCurvatureGuess::BetterThanZero => { | 74 BoundedCurvatureGuess::BetterThanZero => { |
| 95 let opA = self.operator(); | 75 let opA = self.operator(); |
| 96 let b = self.data(); | 76 let b = self.data(); |
| 97 let (ℓ_F0, θ2) = opA.basic_curvature_bound_components(); | 77 let (ℓ_F, ℓ_gradv_0, θ2) = opA.basic_curvature_bound_components(); |
| 98 (ℓ_F0.map(|l| l * b.norm2()), θ2) | 78 (ℓ_F, ℓ_gradv_0.map(|l| l * b.norm2()), θ2) |
| 99 } | 79 } |
| 100 } | 80 } |
| 101 } | 81 } |
| 102 } | 82 } |
