diff -r 1f19c6bbf07b -r 3868555d135c src/mapping/dataterm.rs --- /dev/null Thu Jan 01 00:00:00 1970 +0000 +++ b/src/mapping/dataterm.rs Fri May 15 14:46:30 2026 -0500 @@ -0,0 +1,155 @@ +/*! +General deata terms of the form $g(Ax-b)$ for an operator $A$ +to a [`Euclidean`] space, and a function g on that space. +*/ + +#![allow(non_snake_case)] + +use super::{DifferentiableImpl, DifferentiableMapping, LipschitzDifferentiableImpl, Mapping}; +use crate::convex::ConvexMapping; +use crate::error::DynResult; +use crate::instance::{ClosedSpace, Instance, Space}; +use crate::linops::{BoundedLinear, Linear, Preadjointable}; +use crate::norms::{Normed, L2}; +use crate::types::Float; +use std::ops::Sub; + +//use serde::{Deserialize, Serialize}; + +/// Functions of the form $g(Ax-b)$ for an operator $A$, data $b$, and fidelity $g$. +pub struct DataTerm< + F: Float, + Domain: Space, + A: Mapping, + G: Mapping, +> { + // The operator A + opA: A, + // The data b + b: >::Codomain, + // The outer fidelity + g: G, +} + +// Derive has troubles with `b`. +impl Clone for DataTerm +where + F: Float, + Domain: Space, + A: Mapping + Clone, + >::Codomain: Clone, + G: Mapping + Clone, +{ + fn clone(&self) -> Self { + DataTerm { opA: self.opA.clone(), b: self.b.clone(), g: self.g.clone() } + } +} + +#[allow(non_snake_case)] +impl, G: Mapping> + DataTerm +{ + pub fn new(opA: A, b: A::Codomain, g: G) -> Self { + DataTerm { opA, b, g } + } + + pub fn operator(&self) -> &'_ A { + &self.opA + } + + pub fn data(&self) -> &'_ >::Codomain { + &self.b + } + + pub fn fidelity(&self) -> &'_ G { + &self.g + } + + /// Returns the residual $Ax-b$. + pub fn residual<'a, 'b>(&'b self, x: &'a Domain) -> >::Codomain + where + &'a Domain: Instance, + >::Codomain: + Sub<&'b >::Codomain, Output = >::Codomain>, + { + self.opA.apply(x) - &self.b + } +} + +//+ AdjointProductBoundedBy, P, FloatType = F>, + +impl Mapping for DataTerm +where + F: Float, + X: Space, + A: Mapping, + G: Mapping, + A::Codomain: ClosedSpace + for<'a> Sub<&'a A::Codomain, Output = A::Codomain>, +{ + type Codomain = F; + + fn apply>(&self, x: I) -> F { + // TODO: possibly (if at all more effcient) use GEMV once generalised + // to not require preallocation. However, Rust should be pretty efficient + // at not doing preallocations or anything here, as the result of self.opA.apply() + // can be consumed, so maybe GEMV is no use. + self.g.apply(self.opA.apply(x) - &self.b) + } +} + +impl ConvexMapping for DataTerm +where + F: Float, + X: Normed, + A: Linear, + G: ConvexMapping, + A::Codomain: ClosedSpace + Normed + for<'a> Sub<&'a A::Codomain, Output = A::Codomain>, +{ +} + +impl DifferentiableImpl for DataTerm +where + F: Float, + X: Space, + Y: Space + Instance + for<'a> Sub<&'a Y, Output = Y>, + A: Linear + Preadjointable, + G::DerivativeDomain: Instance, + A::PreadjointCodomain: ClosedSpace, + G: DifferentiableMapping, + Self: Mapping, +{ + type Derivative = A::PreadjointCodomain; + + fn differential_impl>(&self, x: I) -> Self::Derivative { + // TODO: possibly (if at all more effcient) use GEMV once generalised + // to not require preallocation. However, Rust should be pretty efficient + // at not doing preallocations or anything here, as the result of self.opA.apply() + // can be consumed, so maybe GEMV is no use. + //self.opA.preadjoint().apply(self.opA.apply(x) - self.b) + self.opA + .preadjoint() + .apply(self.g.differential(self.opA.apply(x) - &self.b)) + } + + fn apply_and_differential_impl>(&self, x: I) -> (F, Self::Derivative) { + let j = self.opA.apply(x) - &self.b; + let (v, d) = self.g.apply_and_differential(j); + (v, self.opA.preadjoint().apply(d)) + } +} + +impl<'a, F, X, Y, A, G> LipschitzDifferentiableImpl for DataTerm +where + F: Float, + X: Normed, + Y: Normed, + A: BoundedLinear, + G: Mapping + LipschitzDifferentiableImpl, + Self: DifferentiableImpl, +{ + type FloatType = F; + + fn diff_lipschitz_factor(&self, seminorm: X::NormExp) -> DynResult { + Ok(self.opA.opnorm_bound(seminorm, L2)?.powi(2)) + } +}