src/mapping/dataterm.rs

changeset 198
3868555d135c
parent 194
a5ee4bfb0b87
--- /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<Domain>,
+    G: Mapping<A::Codomain, Codomain = F>,
+> {
+    // The operator A
+    opA: A,
+    // The data b
+    b: <A as Mapping<Domain>>::Codomain,
+    // The outer fidelity
+    g: G,
+}
+
+// Derive has troubles with `b`.
+impl<F, Domain, A, G> Clone for DataTerm<F, Domain, A, G>
+where
+    F: Float,
+    Domain: Space,
+    A: Mapping<Domain> + Clone,
+    <A as Mapping<Domain>>::Codomain: Clone,
+    G: Mapping<A::Codomain, Codomain = F> + Clone,
+{
+    fn clone(&self) -> Self {
+        DataTerm { opA: self.opA.clone(), b: self.b.clone(), g: self.g.clone() }
+    }
+}
+
+#[allow(non_snake_case)]
+impl<F: Float, Domain: Space, A: Mapping<Domain>, G: Mapping<A::Codomain, Codomain = F>>
+    DataTerm<F, Domain, A, G>
+{
+    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) -> &'_ <A as Mapping<Domain>>::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) -> <A as Mapping<Domain>>::Codomain
+    where
+        &'a Domain: Instance<Domain>,
+        <A as Mapping<Domain>>::Codomain:
+            Sub<&'b <A as Mapping<Domain>>::Codomain, Output = <A as Mapping<Domain>>::Codomain>,
+    {
+        self.opA.apply(x) - &self.b
+    }
+}
+
+//+ AdjointProductBoundedBy<RNDM<N, F>, P, FloatType = F>,
+
+impl<F, X, A, G> Mapping<X> for DataTerm<F, X, A, G>
+where
+    F: Float,
+    X: Space,
+    A: Mapping<X>,
+    G: Mapping<A::Codomain, Codomain = F>,
+    A::Codomain: ClosedSpace + for<'a> Sub<&'a A::Codomain, Output = A::Codomain>,
+{
+    type Codomain = F;
+
+    fn apply<I: Instance<X>>(&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<F, X, A, G> ConvexMapping<X, F> for DataTerm<F, X, A, G>
+where
+    F: Float,
+    X: Normed<F>,
+    A: Linear<X>,
+    G: ConvexMapping<A::Codomain, F>,
+    A::Codomain: ClosedSpace + Normed<F> + for<'a> Sub<&'a A::Codomain, Output = A::Codomain>,
+{
+}
+
+impl<F, X, Y, A, G> DifferentiableImpl<X> for DataTerm<F, X, A, G>
+where
+    F: Float,
+    X: Space,
+    Y: Space + Instance<Y> + for<'a> Sub<&'a Y, Output = Y>,
+    A: Linear<X, Codomain = Y> + Preadjointable<X, G::DerivativeDomain>,
+    G::DerivativeDomain: Instance<G::DerivativeDomain>,
+    A::PreadjointCodomain: ClosedSpace,
+    G: DifferentiableMapping<Y, Codomain = F>,
+    Self: Mapping<X, Codomain = F>,
+{
+    type Derivative = A::PreadjointCodomain;
+
+    fn differential_impl<I: Instance<X>>(&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<I: Instance<X>>(&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<X, X::NormExp> for DataTerm<F, X, A, G>
+where
+    F: Float,
+    X: Normed<F>,
+    Y: Normed<F>,
+    A: BoundedLinear<X, X::NormExp, L2, F, Codomain = Y>,
+    G: Mapping<Y, Codomain = F> + LipschitzDifferentiableImpl<Y, Y::NormExp>,
+    Self: DifferentiableImpl<X>,
+{
+    type FloatType = F;
+
+    fn diff_lipschitz_factor(&self, seminorm: X::NormExp) -> DynResult<F> {
+        Ok(self.opA.opnorm_bound(seminorm, L2)?.powi(2))
+    }
+}

mercurial