diff -r 9738b51d90d7 -r 4f468d35fa29 src/forward_model.rs --- a/src/forward_model.rs Sun Apr 27 15:03:51 2025 -0500 +++ b/src/forward_model.rs Thu Feb 26 11:38:43 2026 -0500 @@ -2,14 +2,15 @@ Forward models from discrete measures to observations. */ -use alg_tools::error::DynError; -use alg_tools::euclidean::Euclidean; -use alg_tools::instance::Instance; +use crate::dataterm::QuadraticDataTerm; +use crate::measures::{Radon, RNDM}; +use crate::types::*; +use alg_tools::error::{DynError, DynResult}; +use alg_tools::euclidean::{ClosedEuclidean, Euclidean}; pub use alg_tools::linops::*; use alg_tools::norms::{Norm, NormExponent, L2}; +use serde::{Deserialize, Serialize}; -use crate::measures::Radon; -use crate::types::*; pub mod bias; pub mod sensor_grid; @@ -21,13 +22,12 @@ + GEMV + Preadjointable where - for<'a> Self::Observable: Instance, - Domain: Norm, + Domain: Norm, { /// The codomain or value space (of “observables”) for this operator. /// It is assumed to be a [`Euclidean`] space, and therefore also (identified with) /// the domain of the preadjoint. - type Observable: Euclidean + AXPY + Space + Clone; + type Observable: ClosedEuclidean + Clone; /// Write an observable into a file. fn write_observable(&self, b: &Self::Observable, prefix: String) -> DynError; @@ -36,48 +36,66 @@ fn zero_observable(&self) -> Self::Observable; } -/// Trait for operators $A$ for which $A_*A$ is bounded by some other operator. -pub trait AdjointProductBoundedBy: Linear { - type FloatType: Float; - /// Return $L$ such that $A_*A ≤ LD$. - fn adjoint_product_bound(&self, other: &D) -> Option; +/// Guess for [`BoundedCurvature`] calculations. +#[derive(Clone, Copy, Eq, PartialEq, Serialize, Deserialize, Debug)] +pub enum BoundedCurvatureGuess { + /// No iterate $μ^k$ is worse than $μ=0$. + BetterThanZero, } -/// Trait for operators $A$ for which $A_*A$ is bounded by a diagonal operator. -pub trait AdjointProductPairBoundedBy: Linear { - type FloatType: Float; - /// Return $(L, L_z)$ such that $A_*A ≤ (L_1 D_1, L_2 D_2)$. - fn adjoint_product_pair_bound( +/// Curvature error control: helper bounds for (4.2d), (5.2a), (5.2b), (5.15a), and (5.16a). +/// +/// Based on Lemma 5.11 and Example 5.12, the helper bound for (5.15a) and (5.16a) is (3.8). +/// Thus, subject to `guess` being correct, returns factor $(ℓ_F, Θ²)$ such that +/// $B_{F'(μ)} dγ ≤ ℓ_F c_2$ and $⟨F'(μ+Δ)-F'(μ)|Δ⟩ ≤ Θ²|γ|(c_2)$, where $Δ=(π_♯^1-π_♯^0)γ$. +/// +/// This trait is supposed to be implemented by the data term $F$, in the basic case a +/// [`Mapping`] from [`RNDM`] to a [`Float`] `F`. +/// The generic implementation for operators that satisfy [`BasicCurvatureBoundEstimates`] +/// uses Remark 5.15 and Example 5.16 for (4.2d) and (5.2a), and (5.2b); +/// and Lemma 3.8 for (3.8). +pub trait BoundedCurvature { + /// Returns $(ℓ_F, Θ²)$ or individual errors for each. + fn curvature_bound_components( &self, - other1: &D1, - other_2: &D2, - ) -> Option<(Self::FloatType, Self::FloatType)>; + guess: BoundedCurvatureGuess, + ) -> (DynResult, DynResult); } -/* -/// Trait for [`ForwardModel`]s whose preadjoint has Lipschitz values. -pub trait LipschitzValues { - type FloatType : Float; - /// Return (if one exists) a factor $L$ such that $A_*z$ is $L$-Lipschitz for all - /// $z$ in the unit ball. - fn value_unit_lipschitz_factor(&self) -> Option { - None - } +/// Curvature error control: helper bounds for (4.2d), (5.2a), (5.2b), (5.15a), and (5.16a) +/// for quadratic dataterms $F(μ) = \frac{1}{2}\|Aμ-b\|^2$. +/// +/// This trait is to be implemented by the [`Linear`] operator $A$, in the basic from +/// [`RNDM`] to a an Euclidean space. +/// It is used by implementations of [`BoundedCurvature`] for $F$. +/// +/// Based on Lemma 5.11 and Example 5.12, the helper bound for (5.15a) and (5.16a) is (3.8). +/// This trait provides the factor $θ²$ of (3.8) as determined by Lemma 3.8. +/// To aid in calculating (4.2d), (5.2a), (5.2b), motivated by Example 5.16, it also provides +/// $ℓ_F^0$ such that $∇v^k$ $ℓ_F^0 \|Aμ-b\|$-Lipschitz. Here $v^k := F'(∪^k)$. +pub trait BasicCurvatureBoundEstimates { + /// Returns $(ℓ_F^0, Θ²)$ or individual errors for each. + fn basic_curvature_bound_components(&self) -> (DynResult, DynResult); +} - /// Return (if one exists) a factor $L$ such that $∇A_*z$ is $L$-Lipschitz for all - /// $z$ in the unit ball. - fn value_diff_unit_lipschitz_factor(&self) -> Option { - None +impl BoundedCurvature for QuadraticDataTerm, A> +where + F: Float, + Z: Clone + Space + Euclidean, + A: Mapping, Codomain = Z>, + A: BasicCurvatureBoundEstimates, +{ + fn curvature_bound_components( + &self, + guess: BoundedCurvatureGuess, + ) -> (DynResult, DynResult) { + match guess { + BoundedCurvatureGuess::BetterThanZero => { + let opA = self.operator(); + let b = self.data(); + let (ℓ_F0, θ2) = opA.basic_curvature_bound_components(); + (ℓ_F0.map(|l| l * b.norm2()), θ2) + } + } } } -*/ - -/// Trait for [`ForwardModel`]s that satisfy bounds on curvature. -pub trait BoundedCurvature { - type FloatType: Float; - - /// Returns factor $ℓ_F$ and $ℓ_r$ such that - /// $B_{F'(μ)} dγ ≤ ℓ_F c_2$ and $⟨F'(μ)+F'(μ+Δ)|Δ⟩ ≤ ℓ_r|γ|(c_2)$, - /// where $Δ=(π_♯^1-π_♯^0)γ$. - fn curvature_bound_components(&self) -> (Option, Option); -}