--- a/src/convex.rs Mon Apr 28 23:16:56 2025 -0500
+++ b/src/convex.rs Tue Apr 29 00:03:12 2025 -0500
@@ -2,27 +2,29 @@
Some convex analysis basics
*/
-use std::marker::PhantomData;
-use crate::types::*;
+use crate::euclidean::Euclidean;
+use crate::instance::{DecompositionMut, Instance, InstanceMut};
+use crate::linops::{IdOp, Scaled};
use crate::mapping::{Mapping, Space};
-use crate::linops::IdOp;
-use crate::instance::{Instance, InstanceMut, DecompositionMut};
+use crate::norms::*;
use crate::operator_arithmetic::{Constant, Weighted};
-use crate::norms::*;
+use crate::types::*;
+use std::marker::PhantomData;
/// Trait for convex mappings. Has no features, just serves as a constraint
///
/// TODO: should constrain `Mapping::Codomain` to implement a partial order,
/// but this makes everything complicated with little benefit.
-pub trait ConvexMapping<Domain : Space, F : Num = f64> : Mapping<Domain, Codomain = F>
-{}
+pub trait ConvexMapping<Domain: Space, F: Num = f64>: Mapping<Domain, Codomain = F> {}
/// Trait for mappings with a Fenchel conjugate
///
/// The conjugate type has to implement [`ConvexMapping`], but a `Conjugable` mapping need
/// not be convex.
-pub trait Conjugable<Domain : HasDual<F>, F : Num = f64> : Mapping<Domain> {
- type Conjugate<'a> : ConvexMapping<Domain::DualSpace, F> where Self : 'a;
+pub trait Conjugable<Domain: HasDual<F>, F: Num = f64>: Mapping<Domain> {
+ type Conjugate<'a>: ConvexMapping<Domain::DualSpace, F>
+ where
+ Self: 'a;
fn conjugate(&self) -> Self::Conjugate<'_>;
}
@@ -31,12 +33,14 @@
///
/// In contrast to [`Conjugable`], the preconjugate need not implement [`ConvexMapping`],
/// but a `Preconjugable` mapping has to be convex.
-pub trait Preconjugable<Domain, Predual, F : Num = f64> : ConvexMapping<Domain, F>
+pub trait Preconjugable<Domain, Predual, F: Num = f64>: ConvexMapping<Domain, F>
where
- Domain : Space,
- Predual : HasDual<F>
+ Domain: Space,
+ Predual: HasDual<F>,
{
- type Preconjugate<'a> : Mapping<Predual> where Self : 'a;
+ type Preconjugate<'a>: Mapping<Predual>
+ where
+ Self: 'a;
fn preconjugate(&self) -> Self::Preconjugate<'_>;
}
@@ -45,53 +49,54 @@
///
/// The conjugate type has to implement [`ConvexMapping`], but a `Conjugable` mapping need
/// not be convex.
-pub trait Prox<Domain : Space> : Mapping<Domain> {
- type Prox<'a> : Mapping<Domain, Codomain=Domain> where Self : 'a;
+pub trait Prox<Domain: Space>: Mapping<Domain> {
+ type Prox<'a>: Mapping<Domain, Codomain = Domain>
+ where
+ Self: 'a;
/// Returns a proximal mapping with weight τ
- fn prox_mapping(&self, τ : Self::Codomain) -> Self::Prox<'_>;
+ fn prox_mapping(&self, τ: Self::Codomain) -> Self::Prox<'_>;
/// Calculate the proximal mapping with weight τ
- fn prox<I : Instance<Domain>>(&self, τ : Self::Codomain, z : I) -> Domain {
+ fn prox<I: Instance<Domain>>(&self, τ: Self::Codomain, z: I) -> Domain {
self.prox_mapping(τ).apply(z)
}
/// Calculate the proximal mapping with weight τ in-place
- fn prox_mut<'b>(&self, τ : Self::Codomain, y : &'b mut Domain)
+ fn prox_mut<'b>(&self, τ: Self::Codomain, y: &'b mut Domain)
where
- &'b mut Domain : InstanceMut<Domain>,
- Domain:: Decomp : DecompositionMut<Domain>,
- for<'a> &'a Domain : Instance<Domain>,
+ &'b mut Domain: InstanceMut<Domain>,
+ Domain::Decomp: DecompositionMut<Domain>,
+ for<'a> &'a Domain: Instance<Domain>,
{
*y = self.prox(τ, &*y);
}
}
-
/// Constraint to the unit ball of the norm described by `E`.
-pub struct NormConstraint<F : Float, E : NormExponent> {
- radius : F,
- norm : NormMapping<F, E>,
+pub struct NormConstraint<F: Float, E: NormExponent> {
+ radius: F,
+ norm: NormMapping<F, E>,
}
impl<Domain, E, F> ConvexMapping<Domain, F> for NormMapping<F, E>
where
- Domain : Space,
- E : NormExponent,
- F : Float,
- Self : Mapping<Domain, Codomain=F>
-{}
-
+ Domain: Space,
+ E: NormExponent,
+ F: Float,
+ Self: Mapping<Domain, Codomain = F>,
+{
+}
impl<F, E, Domain> Mapping<Domain> for NormConstraint<F, E>
where
- Domain : Space + Norm<F, E>,
- F : Float,
- E : NormExponent,
+ Domain: Space + Norm<F, E>,
+ F: Float,
+ E: NormExponent,
{
type Codomain = F;
- fn apply<I : Instance<Domain>>(&self, d : I) -> F {
+ fn apply<I: Instance<Domain>>(&self, d: I) -> F {
if d.eval(|x| x.norm(self.norm.exponent)) <= self.radius {
F::ZERO
} else {
@@ -102,68 +107,80 @@
impl<Domain, E, F> ConvexMapping<Domain, F> for NormConstraint<F, E>
where
- Domain : Space,
- E : NormExponent,
- F : Float,
- Self : Mapping<Domain, Codomain=F>
-{}
-
+ Domain: Space,
+ E: NormExponent,
+ F: Float,
+ Self: Mapping<Domain, Codomain = F>,
+{
+}
impl<E, F, Domain> Conjugable<Domain, F> for NormMapping<F, E>
where
- E : HasDualExponent,
- F : Float,
- Domain : HasDual<F> + Norm<F, E> + Space,
- <Domain as HasDual<F>>::DualSpace : Norm<F, E::DualExp>
+ E: HasDualExponent,
+ F: Float,
+ Domain: HasDual<F> + Norm<F, E> + Space,
+ <Domain as HasDual<F>>::DualSpace: Norm<F, E::DualExp>,
{
- type Conjugate<'a> = NormConstraint<F, E::DualExp> where Self : 'a;
+ type Conjugate<'a>
+ = NormConstraint<F, E::DualExp>
+ where
+ Self: 'a;
fn conjugate(&self) -> Self::Conjugate<'_> {
NormConstraint {
- radius : F::ONE,
- norm : self.exponent.dual_exponent().as_mapping()
+ radius: F::ONE,
+ norm: self.exponent.dual_exponent().as_mapping(),
}
}
}
impl<C, E, F, Domain> Conjugable<Domain, F> for Weighted<NormMapping<F, E>, C>
where
- C : Constant<Type = F>,
- E : HasDualExponent,
- F : Float,
- Domain : HasDual<F> + Norm<F, E> + Space,
- <Domain as HasDual<F>>::DualSpace : Norm<F, E::DualExp>
+ C: Constant<Type = F>,
+ E: HasDualExponent,
+ F: Float,
+ Domain: HasDual<F> + Norm<F, E> + Space,
+ <Domain as HasDual<F>>::DualSpace: Norm<F, E::DualExp>,
{
- type Conjugate<'a> = NormConstraint<F, E::DualExp> where Self : 'a;
+ type Conjugate<'a>
+ = NormConstraint<F, E::DualExp>
+ where
+ Self: 'a;
fn conjugate(&self) -> Self::Conjugate<'_> {
NormConstraint {
- radius : self.weight.value(),
- norm : self.base_fn.exponent.dual_exponent().as_mapping()
+ radius: self.weight.value(),
+ norm: self.base_fn.exponent.dual_exponent().as_mapping(),
}
}
}
impl<Domain, E, F> Prox<Domain> for NormConstraint<F, E>
where
- Domain : Space + Norm<F, E>,
- E : NormExponent,
- F : Float,
- NormProjection<F, E> : Mapping<Domain, Codomain=Domain>,
+ Domain: Space + Norm<F, E>,
+ E: NormExponent,
+ F: Float,
+ NormProjection<F, E>: Mapping<Domain, Codomain = Domain>,
{
- type Prox<'a> = NormProjection<F, E> where Self : 'a;
+ type Prox<'a>
+ = NormProjection<F, E>
+ where
+ Self: 'a;
#[inline]
- fn prox_mapping(&self, _τ : Self::Codomain) -> Self::Prox<'_> {
+ fn prox_mapping(&self, _τ: Self::Codomain) -> Self::Prox<'_> {
assert!(self.radius >= F::ZERO);
- NormProjection{ radius : self.radius, exponent : self.norm.exponent }
+ NormProjection {
+ radius: self.radius,
+ exponent: self.norm.exponent,
+ }
}
}
/// Projection to the unit ball of the norm described by `E`.
-pub struct NormProjection<F : Float, E : NormExponent> {
- radius : F,
- exponent : E,
+pub struct NormProjection<F: Float, E: NormExponent> {
+ radius: F,
+ exponent: E,
}
/*
@@ -182,41 +199,42 @@
impl<F, E, Domain> Mapping<Domain> for NormProjection<F, E>
where
- Domain : Space + Projection<F, E>,
- F : Float,
- E : NormExponent,
+ Domain: Space + Projection<F, E>,
+ F: Float,
+ E: NormExponent,
{
type Codomain = Domain;
- fn apply<I : Instance<Domain>>(&self, d : I) -> Domain {
+ fn apply<I: Instance<Domain>>(&self, d: I) -> Domain {
d.own().proj_ball(self.radius, self.exponent)
}
}
+/// The zero mapping
+pub struct Zero<Domain: Space, F: Num>(PhantomData<(Domain, F)>);
-/// The zero mapping
-pub struct Zero<Domain : Space, F : Num>(PhantomData<(Domain, F)>);
-
-impl<Domain : Space, F : Num> Zero<Domain, F> {
+impl<Domain: Space, F: Num> Zero<Domain, F> {
pub fn new() -> Self {
Zero(PhantomData)
}
}
-impl<Domain : Space, F : Num> Mapping<Domain> for Zero<Domain, F> {
+impl<Domain: Space, F: Num> Mapping<Domain> for Zero<Domain, F> {
type Codomain = F;
/// Compute the value of `self` at `x`.
- fn apply<I : Instance<Domain>>(&self, _x : I) -> Self::Codomain {
+ fn apply<I: Instance<Domain>>(&self, _x: I) -> Self::Codomain {
F::ZERO
}
}
-impl<Domain : Space, F : Num> ConvexMapping<Domain, F> for Zero<Domain, F> { }
-
+impl<Domain: Space, F: Num> ConvexMapping<Domain, F> for Zero<Domain, F> {}
-impl<Domain : HasDual<F>, F : Float> Conjugable<Domain, F> for Zero<Domain, F> {
- type Conjugate<'a> = ZeroIndicator<Domain::DualSpace, F> where Self : 'a;
+impl<Domain: HasDual<F>, F: Float> Conjugable<Domain, F> for Zero<Domain, F> {
+ type Conjugate<'a>
+ = ZeroIndicator<Domain::DualSpace, F>
+ where
+ Self: 'a;
#[inline]
fn conjugate(&self) -> Self::Conjugate<'_> {
@@ -224,12 +242,15 @@
}
}
-impl<Domain, Predual, F : Float> Preconjugable<Domain, Predual, F> for Zero<Domain, F>
+impl<Domain, Predual, F: Float> Preconjugable<Domain, Predual, F> for Zero<Domain, F>
where
- Domain : Space,
- Predual : HasDual<F>
+ Domain: Space,
+ Predual: HasDual<F>,
{
- type Preconjugate<'a> = ZeroIndicator<Predual, F> where Self : 'a;
+ type Preconjugate<'a>
+ = ZeroIndicator<Predual, F>
+ where
+ Self: 'a;
#[inline]
fn preconjugate(&self) -> Self::Preconjugate<'_> {
@@ -237,38 +258,43 @@
}
}
-impl<Domain : Space + Clone, F : Num> Prox<Domain> for Zero<Domain, F> {
- type Prox<'a> = IdOp<Domain> where Self : 'a;
+impl<Domain: Space + Clone, F: Num> Prox<Domain> for Zero<Domain, F> {
+ type Prox<'a>
+ = IdOp<Domain>
+ where
+ Self: 'a;
#[inline]
- fn prox_mapping(&self, _τ : Self::Codomain) -> Self::Prox<'_> {
+ fn prox_mapping(&self, _τ: Self::Codomain) -> Self::Prox<'_> {
IdOp::new()
}
}
+/// The zero indicator
+pub struct ZeroIndicator<Domain: Space, F: Num>(PhantomData<(Domain, F)>);
-/// The zero indicator
-pub struct ZeroIndicator<Domain : Space, F : Num>(PhantomData<(Domain, F)>);
-
-impl<Domain : Space, F : Num> ZeroIndicator<Domain, F> {
+impl<Domain: Space, F: Num> ZeroIndicator<Domain, F> {
pub fn new() -> Self {
ZeroIndicator(PhantomData)
}
}
-impl<Domain : Normed<F>, F : Float> Mapping<Domain> for ZeroIndicator<Domain, F> {
+impl<Domain: Normed<F>, F: Float> Mapping<Domain> for ZeroIndicator<Domain, F> {
type Codomain = F;
/// Compute the value of `self` at `x`.
- fn apply<I : Instance<Domain>>(&self, x : I) -> Self::Codomain {
+ fn apply<I: Instance<Domain>>(&self, x: I) -> Self::Codomain {
x.eval(|x̃| if x̃.is_zero() { F::ZERO } else { F::INFINITY })
}
}
-impl<Domain : Normed<F>, F : Float> ConvexMapping<Domain, F> for ZeroIndicator<Domain, F> { }
+impl<Domain: Normed<F>, F: Float> ConvexMapping<Domain, F> for ZeroIndicator<Domain, F> {}
-impl<Domain : HasDual<F>, F : Float> Conjugable<Domain, F> for ZeroIndicator<Domain, F> {
- type Conjugate<'a> = Zero<Domain::DualSpace, F> where Self : 'a;
+impl<Domain: HasDual<F>, F: Float> Conjugable<Domain, F> for ZeroIndicator<Domain, F> {
+ type Conjugate<'a>
+ = Zero<Domain::DualSpace, F>
+ where
+ Self: 'a;
#[inline]
fn conjugate(&self) -> Self::Conjugate<'_> {
@@ -276,15 +302,77 @@
}
}
-impl<Domain, Predual, F : Float> Preconjugable<Domain, Predual, F> for ZeroIndicator<Domain, F>
+impl<Domain, Predual, F: Float> Preconjugable<Domain, Predual, F> for ZeroIndicator<Domain, F>
where
- Domain : Normed<F>,
- Predual : HasDual<F>
+ Domain: Normed<F>,
+ Predual: HasDual<F>,
{
- type Preconjugate<'a> = Zero<Predual, F> where Self : 'a;
+ type Preconjugate<'a>
+ = Zero<Predual, F>
+ where
+ Self: 'a;
#[inline]
fn preconjugate(&self) -> Self::Preconjugate<'_> {
Zero::new()
}
}
+
+/// The squared Euclidean norm divided by two
+pub struct Norm222<Domain: Space, F: Float>(PhantomData<(Domain, F)>);
+
+impl<Domain: Euclidean<F>, F: Float> Norm222<Domain, F> {
+ pub fn new() -> Self {
+ Norm222(PhantomData)
+ }
+}
+
+impl<Domain: Euclidean<F>, F: Float> Mapping<Domain> for Norm222<Domain, F> {
+ type Codomain = F;
+
+ /// Compute the value of `self` at `x`.
+ fn apply<I: Instance<Domain>>(&self, x: I) -> Self::Codomain {
+ x.eval(|z| z.norm2_squared() / F::TWO)
+ }
+}
+
+impl<Domain: Euclidean<F>, F: Float> ConvexMapping<Domain, F> for Norm222<Domain, F> {}
+
+impl<Domain: Euclidean<F>, F: Float> Conjugable<Domain, F> for Norm222<Domain, F> {
+ type Conjugate<'a>
+ = Self
+ where
+ Self: 'a;
+
+ #[inline]
+ fn conjugate(&self) -> Self::Conjugate<'_> {
+ Self::new()
+ }
+}
+
+impl<Domain: Euclidean<F>, F: Float> Preconjugable<Domain, Domain, F> for Norm222<Domain, F> {
+ type Preconjugate<'a>
+ = Self
+ where
+ Self: 'a;
+
+ #[inline]
+ fn preconjugate(&self) -> Self::Preconjugate<'_> {
+ Self::new()
+ }
+}
+
+impl<Domain, F> Prox<Domain> for Norm222<Domain, F>
+where
+ F: Float,
+ Domain: Euclidean<F, Output = Domain>,
+{
+ type Prox<'a>
+ = Scaled<F>
+ where
+ Self: 'a;
+
+ fn prox_mapping(&self, τ: F) -> Self::Prox<'_> {
+ Scaled(F::ONE / (F::ONE + τ))
+ }
+}