Mercurial > repos > pointsource_algs / changeset
changeset
Generalise algorithms to arbitrary domains
Thu, 06 Aug 2026 11:14:15 -0500
- author
- Tuomo Valkonen <tuomov@iki.fi>
- date
- Thu, 06 Aug 2026 11:14:15 -0500
- changeset 93
- 8f1e57821d76
- parent 92
- 759634edb536
- child 94
- 88d994f2d63d
Generalise algorithms to arbitrary domains
| src/domain.rs | file | annotate | diff | comparison | revisions | |
| src/fb.rs | file | annotate | diff | comparison | revisions | |
| src/forward_pdps.rs | file | annotate | diff | comparison | revisions | |
| src/frank_wolfe.rs | file | annotate | diff | comparison | revisions | |
| src/pdps.rs | file | annotate | diff | comparison | revisions | |
| src/run.rs | file | annotate | diff | comparison | revisions | |
| src/sliding_fb.rs | file | annotate | diff | comparison | revisions | |
| src/sliding_pdps.rs | file | annotate | diff | comparison | revisions |
--- a/src/domain.rs Tue Aug 04 20:11:10 2026 -0500 +++ b/src/domain.rs Thu Aug 06 11:14:15 2026 -0500 @@ -28,7 +28,7 @@ } } -impl<F: Num + PartialOrd + Sized, const N: usize> AdmissibleDomain<Loc<N, F>> for FullSpace {} +impl<U: Space> AdmissibleDomain<U> for FullSpace {} pub struct UniformRandom<S, R: Rng = StdRng>(pub S, pub R);
--- a/src/fb.rs Tue Aug 04 20:11:10 2026 -0500 +++ b/src/fb.rs Thu Aug 06 11:14:15 2026 -0500 @@ -79,7 +79,7 @@ use crate::domain::{AdmissibleDomain, FullSpace}; use crate::measures::merging::SpikeMerging; -use crate::measures::{DiscreteMeasure, RNDM}; +use crate::measures::DiscreteMeasure; use crate::plot::Plotter; pub use crate::prox_penalty::{InsertionConfig, ProxPenalty, StepLengthBound}; use crate::regularisation::RegTerm; @@ -87,6 +87,7 @@ use alg_tools::error::DynResult; use alg_tools::instance::Instance; use alg_tools::iterate::AlgIteratorFactory; +use alg_tools::linops::Space; use alg_tools::mapping::DifferentiableMapping; use alg_tools::nalgebra_support::ToNalgebraRealField; use colored::Colorize; @@ -113,8 +114,8 @@ } } -pub(crate) fn prune_with_stats<F: Float, const N: usize>( - μ: &mut RNDM<N, F>, +pub(crate) fn prune_with_stats<U: Space, F: Float>( + μ: &mut DiscreteMeasure<U, F>, until: usize, ) -> usize { let n_before_prune = μ.len(); @@ -124,14 +125,16 @@ } #[replace_float_literals(F::cast_from(literal))] -pub(crate) fn postprocess<F: Float, Dat: Fn(&RNDM<N, F>) -> F, const N: usize>( - mut μ: RNDM<N, F>, +pub(crate) fn postprocess<F: Float, U, Dat>( + mut μ: DiscreteMeasure<U, F>, config: &InsertionConfig<F>, f: Dat, -) -> DynResult<RNDM<N, F>> +) -> DynResult<DiscreteMeasure<U, F>> where - RNDM<N, F>: SpikeMerging<F>, - for<'a> &'a RNDM<N, F>: Instance<RNDM<N, F>>, + Dat: Fn(&DiscreteMeasure<U, F>) -> F, + U: Space + Clone, + DiscreteMeasure<U, F>: SpikeMerging<F>, + for<'a> &'a DiscreteMeasure<U, F>: Instance<DiscreteMeasure<U, F>>, { //μ.merge_spikes_fitness(config.final_merging_method(), |μ̃| f.apply(μ̃), |&v| v); μ.merge_spikes_fitness(config.final_merging_method(), f, |&v| v); @@ -151,25 +154,26 @@ /// /// Returns the final iterate. #[replace_float_literals(F::cast_from(literal))] -pub fn pointsource_fb_dom<F, I, Dat, Reg, P, Plot, const N: usize>( +pub fn pointsource_fb_dom<F, I, Dat, Reg, P, Plot, U>( f: &Dat, reg: &Reg, prox_penalty: &P, fbconfig: &FBConfig<F>, iterator: I, mut plotter: Plot, - μ0: Option<RNDM<N, F>>, - domain: &mut impl AdmissibleDomain<Loc<N, F>>, -) -> DynResult<RNDM<N, F>> + μ0: Option<DiscreteMeasure<U, F>>, + domain: &mut impl AdmissibleDomain<U>, +) -> DynResult<DiscreteMeasure<U, F>> where + U: Space + Clone, F: Float + ToNalgebraRealField, I: AlgIteratorFactory<IterInfo<F>>, - RNDM<N, F>: SpikeMerging<F>, - Dat: DifferentiableMapping<RNDM<N, F>, Codomain = F>, + DiscreteMeasure<U, F>: SpikeMerging<F>, + Dat: DifferentiableMapping<DiscreteMeasure<U, F>, Codomain = F>, Dat::DerivativeDomain: ClosedMul<F>, - Reg: RegTerm<Loc<N, F>, F>, - P: ProxPenalty<Loc<N, F>, Dat::DerivativeDomain, Reg, F> + StepLengthBound<F, Dat>, - Plot: Plotter<P::ReturnMapping, Dat::DerivativeDomain, RNDM<N, F>>, + Reg: RegTerm<U, F>, + P: ProxPenalty<U, Dat::DerivativeDomain, Reg, F> + StepLengthBound<F, Dat>, + Plot: Plotter<P::ReturnMapping, Dat::DerivativeDomain, DiscreteMeasure<U, F>>, { // Set up parameters let config = &fbconfig.insertion; @@ -183,7 +187,7 @@ let mut μ = μ0.unwrap_or_else(|| DiscreteMeasure::new()); // Statistics - let full_stats = |μ: &RNDM<N, F>, ε, stats| IterInfo { + let full_stats = |μ: &DiscreteMeasure<U, F>, ε, stats| IterInfo { value: f.apply(μ) + reg.apply(μ), n_spikes: μ.len(), ε, @@ -219,7 +223,7 @@ ε, ins, ®, - Some(|μ̃: &RNDM<N, F>| f.apply(μ̃) + reg.apply(μ̃)), + Some(|μ̃: &DiscreteMeasure<U, F>| f.apply(μ̃) + reg.apply(μ̃)), ); // Prune. Do not prune new insertions on first round. @@ -267,24 +271,25 @@ } #[deprecated] -pub fn pointsource_fb_reg<F, I, Dat, Reg, P, Plot, const N: usize>( +pub fn pointsource_fb_reg<F, U, I, Dat, Reg, P, Plot>( f: &Dat, reg: &Reg, prox_penalty: &P, fbconfig: &FBConfig<F>, iterator: I, plotter: Plot, - μ0: Option<RNDM<N, F>>, -) -> DynResult<RNDM<N, F>> + μ0: Option<DiscreteMeasure<U, F>>, +) -> DynResult<DiscreteMeasure<U, F>> where + U: Space + Clone, F: Float + ToNalgebraRealField, I: AlgIteratorFactory<IterInfo<F>>, - RNDM<N, F>: SpikeMerging<F>, - Dat: DifferentiableMapping<RNDM<N, F>, Codomain = F>, + DiscreteMeasure<U, F>: SpikeMerging<F>, + Dat: DifferentiableMapping<DiscreteMeasure<U, F>, Codomain = F>, Dat::DerivativeDomain: ClosedMul<F>, - Reg: RegTerm<Loc<N, F>, F>, - P: ProxPenalty<Loc<N, F>, Dat::DerivativeDomain, Reg, F> + StepLengthBound<F, Dat>, - Plot: Plotter<P::ReturnMapping, Dat::DerivativeDomain, RNDM<N, F>>, + Reg: RegTerm<U, F>, + P: ProxPenalty<U, Dat::DerivativeDomain, Reg, F> + StepLengthBound<F, Dat>, + Plot: Plotter<P::ReturnMapping, Dat::DerivativeDomain, DiscreteMeasure<U, F>>, { pointsource_fb_dom( f, @@ -310,24 +315,25 @@ /// /// Returns the final iterate. #[replace_float_literals(F::cast_from(literal))] -pub fn pointsource_fista_reg<F, I, Dat, Reg, P, Plot, const N: usize>( +pub fn pointsource_fista_reg<F, U, I, Dat, Reg, P, Plot>( f: &Dat, reg: &Reg, prox_penalty: &P, fbconfig: &FBConfig<F>, iterator: I, mut plotter: Plot, - μ0: Option<RNDM<N, F>>, -) -> DynResult<RNDM<N, F>> + μ0: Option<DiscreteMeasure<U, F>>, +) -> DynResult<DiscreteMeasure<U, F>> where + U: Space + Clone, F: Float + ToNalgebraRealField, I: AlgIteratorFactory<IterInfo<F>>, - RNDM<N, F>: SpikeMerging<F>, - Dat: DifferentiableMapping<RNDM<N, F>, Codomain = F>, + DiscreteMeasure<U, F>: SpikeMerging<F>, + Dat: DifferentiableMapping<DiscreteMeasure<U, F>, Codomain = F>, Dat::DerivativeDomain: ClosedMul<F>, - Reg: RegTerm<Loc<N, F>, F>, - P: ProxPenalty<Loc<N, F>, Dat::DerivativeDomain, Reg, F> + StepLengthBound<F, Dat>, - Plot: Plotter<P::ReturnMapping, Dat::DerivativeDomain, RNDM<N, F>>, + Reg: RegTerm<U, F>, + P: ProxPenalty<U, Dat::DerivativeDomain, Reg, F> + StepLengthBound<F, Dat>, + Plot: Plotter<P::ReturnMapping, Dat::DerivativeDomain, DiscreteMeasure<U, F>>, { // Set up parameters let config = &fbconfig.insertion; @@ -344,7 +350,7 @@ let mut warned_merging = false; // Statistics - let full_stats = |ν: &RNDM<N, F>, ε, stats| IterInfo { + let full_stats = |ν: &DiscreteMeasure<U, F>, ε, stats| IterInfo { value: f.apply(ν) + reg.apply(ν), n_spikes: ν.len(), ε,
--- a/src/forward_pdps.rs Tue Aug 04 20:11:10 2026 -0500 +++ b/src/forward_pdps.rs Thu Aug 06 11:14:15 2026 -0500 @@ -6,7 +6,7 @@ use crate::domain::{AdmissibleDomain, FullSpace}; use crate::fb::*; use crate::measures::merging::SpikeMerging; -use crate::measures::{DiscreteMeasure, RNDM}; +use crate::measures::DiscreteMeasure; use crate::plot::Plotter; use crate::prox_penalty::{ProxPenalty, StepLengthBoundPair}; use crate::regularisation::RegTerm; @@ -16,8 +16,8 @@ use alg_tools::error::DynResult; use alg_tools::euclidean::ClosedEuclidean; use alg_tools::iterate::AlgIteratorFactory; -use alg_tools::linops::{BoundedLinear, IdOp, SimplyAdjointable, ZeroOp, AXPY, GEMV}; -use alg_tools::mapping::{DifferentiableMapping, DifferentiableRealMapping, Instance}; +use alg_tools::linops::{BoundedLinear, IdOp, SimplyAdjointable, Space, ZeroOp, AXPY, GEMV}; +use alg_tools::mapping::{DifferentiableMapping, Instance}; use alg_tools::nalgebra_support::ToNalgebraRealField; use alg_tools::norms::L2; use anyhow::ensure; @@ -48,7 +48,7 @@ } } -type MeasureZ<F, Z, const N: usize> = Pair<RNDM<N, F>, Z>; +type MeasureZ<F, Z, U> = Pair<DiscreteMeasure<U, F>, Z>; /// Iteratively solve the pointsource localisation with an additional variable /// using primal-dual proximal splitting with a forward step. @@ -89,6 +89,7 @@ #[replace_float_literals(F::cast_from(literal))] pub fn pointsource_forward_pdps_pair_dom< F, + U, I, S, Dat, @@ -100,7 +101,6 @@ /*KOpM, */ KOpZ, H, Plot, - const N: usize, >( f: &Dat, reg: &Reg, @@ -108,22 +108,23 @@ config: &ForwardPDPSConfig<F>, iterator: I, mut plotter: Plot, - (μ0, mut z, mut y): (Option<RNDM<N, F>>, Z, Y), + (μ0, mut z, mut y): (Option<DiscreteMeasure<U, F>>, Z, Y), //opKμ : KOpM, opKz: &KOpZ, fnR: &R, fnH: &H, - domain: &mut impl AdmissibleDomain<Loc<N, F>>, -) -> DynResult<MeasureZ<F, Z, N>> + domain: &mut impl AdmissibleDomain<U>, +) -> DynResult<MeasureZ<F, Z, U>> where + U: Space + Clone, F: Float + ToNalgebraRealField, I: AlgIteratorFactory<IterInfo<F>>, - Dat: DifferentiableMapping<MeasureZ<F, Z, N>, Codomain = F, DerivativeDomain = Pair<S, Z>>, + Dat: DifferentiableMapping<MeasureZ<F, Z, U>, Codomain = F, DerivativeDomain = Pair<S, Z>>, //Pair<S, Z>: ClosedMul<F>, // Doesn't really need to be closed, if make this signature more complex… - S: DifferentiableRealMapping<N, F> + ClosedMul<F>, - RNDM<N, F>: SpikeMerging<F>, - Reg: RegTerm<Loc<N, F>, F>, - P: ProxPenalty<Loc<N, F>, S, Reg, F>, + S: ClosedMul<F>, + DiscreteMeasure<U, F>: SpikeMerging<F>, + Reg: RegTerm<U, F>, + P: ProxPenalty<U, S, Reg, F>, for<'a> Pair<&'a P, &'a IdOp<Z>>: StepLengthBoundPair<F, Dat>, KOpZ: BoundedLinear<Z, L2, L2, F, Codomain = Y> + GEMV<F, Z, Y> @@ -136,7 +137,7 @@ R: Prox<Z, Codomain = F>, H: Conjugable<Y, F, Codomain = F>, for<'b> H::Conjugate<'b>: Prox<Y>, - Plot: Plotter<P::ReturnMapping, S, RNDM<N, F>>, + Plot: Plotter<P::ReturnMapping, S, DiscreteMeasure<U, F>>, { // Check parameters ensure!( @@ -188,7 +189,7 @@ let starH = fnH.conjugate(); // Statistics - let full_stats = |μ: &RNDM<N, F>, z: &Z, ε, stats| IterInfo { + let full_stats = |μ: &DiscreteMeasure<U, F>, z: &Z, ε, stats| IterInfo { value: f.apply(Pair(μ, z)) + fnR.apply(z) + reg.apply(μ) @@ -250,7 +251,7 @@ ε, ins, ®, - Some(|μ̃: &RNDM<N, F>| f.apply(Pair(μ̃, &z)) + reg.apply(μ̃)), + Some(|μ̃: &DiscreteMeasure<U, F>| f.apply(Pair(μ̃, &z)) + reg.apply(μ̃)), ); // Prune. Do not prune new insertions on first round. stats.pruned += prune_with_stats(&mut μ, μ_base_len); @@ -310,7 +311,7 @@ ε = tolerance.update(ε, iter); } - let fit = |μ̃: &RNDM<N, F>| { + let fit = |μ̃: &DiscreteMeasure<U, F>| { f.apply(Pair(μ̃, &z)) /*+ fnR.apply(z) + reg.apply(μ)*/ + fnH.apply(/* opKμ.apply(&μ̃) + */ opKz.apply(&z)) + reg.apply(μ̃) @@ -324,6 +325,7 @@ #[deprecated] pub fn pointsource_forward_pdps_pair< F, + U, I, S, Dat, @@ -335,7 +337,6 @@ /*KOpM, */ KOpZ, H, Plot, - const N: usize, >( f: &Dat, reg: &Reg, @@ -343,21 +344,22 @@ config: &ForwardPDPSConfig<F>, iterator: I, plotter: Plot, - init: (Option<RNDM<N, F>>, Z, Y), + init: (Option<DiscreteMeasure<U, F>>, Z, Y), //opKμ : KOpM, opKz: &KOpZ, fnR: &R, fnH: &H, -) -> DynResult<MeasureZ<F, Z, N>> +) -> DynResult<MeasureZ<F, Z, U>> where + U: Space + Clone, F: Float + ToNalgebraRealField, I: AlgIteratorFactory<IterInfo<F>>, - Dat: DifferentiableMapping<MeasureZ<F, Z, N>, Codomain = F, DerivativeDomain = Pair<S, Z>>, + Dat: DifferentiableMapping<MeasureZ<F, Z, U>, Codomain = F, DerivativeDomain = Pair<S, Z>>, //Pair<S, Z>: ClosedMul<F>, // Doesn't really need to be closed, if make this signature more complex… - S: DifferentiableRealMapping<N, F> + ClosedMul<F>, - RNDM<N, F>: SpikeMerging<F>, - Reg: RegTerm<Loc<N, F>, F>, - P: ProxPenalty<Loc<N, F>, S, Reg, F>, + S: ClosedMul<F>, + DiscreteMeasure<U, F>: SpikeMerging<F>, + Reg: RegTerm<U, F>, + P: ProxPenalty<U, S, Reg, F>, for<'a> Pair<&'a P, &'a IdOp<Z>>: StepLengthBoundPair<F, Dat>, KOpZ: BoundedLinear<Z, L2, L2, F, Codomain = Y> + GEMV<F, Z, Y> @@ -370,7 +372,7 @@ R: Prox<Z, Codomain = F>, H: Conjugable<Y, F, Codomain = F>, for<'b> H::Conjugate<'b>: Prox<Y>, - Plot: Plotter<P::ReturnMapping, S, RNDM<N, F>>, + Plot: Plotter<P::ReturnMapping, S, DiscreteMeasure<U, F>>, { pointsource_forward_pdps_pair_dom( f, @@ -393,31 +395,32 @@ /// The implementation uses [`pointsource_forward_pdps_pair`] with appropriate dummy /// variables, operators, and functions. #[replace_float_literals(F::cast_from(literal))] -pub fn pointsource_fb_pair_dom<F, I, S, Dat, Reg, P, Z, R, Plot, const N: usize>( +pub fn pointsource_fb_pair_dom<F, U, I, S, Dat, Reg, P, Z, R, Plot>( f: &Dat, reg: &Reg, prox_penalty: &P, config: &FBConfig<F>, iterator: I, plotter: Plot, - (μ0, z): (Option<RNDM<N, F>>, Z), + (μ0, z): (Option<DiscreteMeasure<U, F>>, Z), //opKμ : KOpM, fnR: &R, - domain: &mut impl AdmissibleDomain<Loc<N, F>>, -) -> DynResult<MeasureZ<F, Z, N>> + domain: &mut impl AdmissibleDomain<U>, +) -> DynResult<MeasureZ<F, Z, U>> where + U: Space + Clone, F: Float + ToNalgebraRealField, I: AlgIteratorFactory<IterInfo<F>>, - Dat: DifferentiableMapping<MeasureZ<F, Z, N>, Codomain = F, DerivativeDomain = Pair<S, Z>>, - S: DifferentiableRealMapping<N, F> + ClosedMul<F>, - RNDM<N, F>: SpikeMerging<F>, - Reg: RegTerm<Loc<N, F>, F>, - P: ProxPenalty<Loc<N, F>, S, Reg, F>, + Dat: DifferentiableMapping<MeasureZ<F, Z, U>, Codomain = F, DerivativeDomain = Pair<S, Z>>, + S: ClosedMul<F>, + DiscreteMeasure<U, F>: SpikeMerging<F>, + Reg: RegTerm<U, F>, + P: ProxPenalty<U, S, Reg, F>, for<'a> Pair<&'a P, &'a IdOp<Z>>: StepLengthBoundPair<F, Dat>, Z: ClosedEuclidean<F> + AXPY<Field = F> + Clone, for<'b> &'b Z: Instance<Z>, R: Prox<Z, Codomain = F>, - Plot: Plotter<P::ReturnMapping, S, RNDM<N, F>>, + Plot: Plotter<P::ReturnMapping, S, DiscreteMeasure<U, F>>, // We should not need to explicitly require this: for<'b> &'b Loc<0, F>: Instance<Loc<0, F>>, { @@ -445,30 +448,31 @@ #[deprecated] #[replace_float_literals(F::cast_from(literal))] -pub fn pointsource_fb_pair<F, I, S, Dat, Reg, P, Z, R, Plot, const N: usize>( +pub fn pointsource_fb_pair<F, U, I, S, Dat, Reg, P, Z, R, Plot>( f: &Dat, reg: &Reg, prox_penalty: &P, config: &FBConfig<F>, iterator: I, plotter: Plot, - init: (Option<RNDM<N, F>>, Z), + init: (Option<DiscreteMeasure<U, F>>, Z), //opKμ : KOpM, fnR: &R, -) -> DynResult<MeasureZ<F, Z, N>> +) -> DynResult<MeasureZ<F, Z, U>> where + U: Space + Clone, F: Float + ToNalgebraRealField, I: AlgIteratorFactory<IterInfo<F>>, - Dat: DifferentiableMapping<MeasureZ<F, Z, N>, Codomain = F, DerivativeDomain = Pair<S, Z>>, - S: DifferentiableRealMapping<N, F> + ClosedMul<F>, - RNDM<N, F>: SpikeMerging<F>, - Reg: RegTerm<Loc<N, F>, F>, - P: ProxPenalty<Loc<N, F>, S, Reg, F>, + Dat: DifferentiableMapping<MeasureZ<F, Z, U>, Codomain = F, DerivativeDomain = Pair<S, Z>>, + S: DifferentiableMapping<U, Codomain = F> + ClosedMul<F>, + DiscreteMeasure<U, F>: SpikeMerging<F>, + Reg: RegTerm<U, F>, + P: ProxPenalty<U, S, Reg, F>, for<'a> Pair<&'a P, &'a IdOp<Z>>: StepLengthBoundPair<F, Dat>, Z: ClosedEuclidean<F> + AXPY<Field = F> + Clone, for<'b> &'b Z: Instance<Z>, R: Prox<Z, Codomain = F>, - Plot: Plotter<P::ReturnMapping, S, RNDM<N, F>>, + Plot: Plotter<P::ReturnMapping, S, DiscreteMeasure<U, F>>, // We should not need to explicitly require this: for<'b> &'b Loc<0, F>: Instance<Loc<0, F>>, {
--- a/src/frank_wolfe.rs Tue Aug 04 20:11:10 2026 -0500 +++ b/src/frank_wolfe.rs Thu Aug 06 11:14:15 2026 -0500 @@ -20,7 +20,7 @@ use crate::dataterm::QuadraticDataTerm; use crate::forward_model::ForwardModel; use crate::measures::merging::{SpikeMerging, SpikeMergingMethod}; -use crate::measures::{DeltaMeasure, DiscreteMeasure, Radon, RNDM}; +use crate::measures::{DeltaMeasure, DiscreteMeasure, Radon}; use crate::plot::Plotter; use crate::regularisation::{NonnegRadonRegTerm, RadonRegTerm, RegTerm}; #[allow(unused_imports)] // Used in documentation @@ -29,18 +29,15 @@ }; use crate::tolerance::Tolerance; use crate::types::*; -use alg_tools::bisection_tree::P2Minimise; use alg_tools::bounds::MinMaxMapping; use alg_tools::error::DynResult; use alg_tools::euclidean::Euclidean; use alg_tools::instance::Instance; use alg_tools::iterate::{AlgIteratorFactory, AlgIteratorOptions, ValueIteratorFactory}; -use alg_tools::linops::Mapping; -use alg_tools::loc::Loc; +use alg_tools::linops::{Mapping, Space}; use alg_tools::nalgebra_support::ToNalgebraRealField; use alg_tools::norms::Norm; use alg_tools::norms::L2; -use alg_tools::sets::Cube; /// Settings for [`pointsource_fw_reg`]. #[derive(Clone, Copy, Eq, PartialEq, Serialize, Deserialize, Debug)] @@ -106,9 +103,9 @@ /// Trait for finite dimensional weight optimisation. pub trait WeightOptim< F: Float + ToNalgebraRealField, - A: ForwardModel<RNDM<N, F>, F>, + U: Space + Clone + PartialEq, + A: ForwardModel<DiscreteMeasure<U, F>, F>, I: AlgIteratorFactory<F>, - const N: usize, > { /// Return a pre-initialisation struct for [`Self::optimise_weights`]. @@ -135,7 +132,7 @@ /// Returns the number of iterations taken by the method configured in `inner`. fn optimise_weights<'a>( &self, - μ: &mut RNDM<N, F>, + μ: &mut DiscreteMeasure<U, F>, opA: &'a A, b: &A::Observable, findim_data: &FindimData<F>, @@ -147,10 +144,10 @@ /// Trait for regularisation terms supported by [`pointsource_fw_reg`]. pub trait RegTermFW< F: Float + ToNalgebraRealField, - A: ForwardModel<RNDM<N, F>, F>, + U: Space + Clone + PartialEq, + A: ForwardModel<DiscreteMeasure<U, F>, F>, I: AlgIteratorFactory<F>, - const N: usize, ->: RegTerm<Loc<N, F>, F> + WeightOptim<F, A, I, N> + Mapping<RNDM<N, F>, Codomain = F> +>: RegTerm<U, F> + WeightOptim<F, U, A, I> + Mapping<DiscreteMeasure<U, F>, Codomain = F> { /// With $g = A\_\*(Aμ-b)$, returns $(x, g(x))$ for $x$ a new point to be inserted /// into $μ$, as determined by the regulariser. @@ -162,26 +159,26 @@ g: &mut A::PreadjointCodomain, refinement_tolerance: F, max_steps: usize, - ) -> (Loc<N, F>, F); + ) -> (U, F); /// Insert point `ξ` into `μ` for the relaxed algorithm from Bredies–Pikkarainen. fn relaxed_insert<'a>( &self, - μ: &mut RNDM<N, F>, + μ: &mut DiscreteMeasure<U, F>, g: &A::PreadjointCodomain, opA: &'a A, - ξ: Loc<N, F>, + ξ: U, v_ξ: F, findim_data: &FindimData<F>, ); } #[replace_float_literals(F::cast_from(literal))] -impl<F: Float + ToNalgebraRealField, A, I, const N: usize> WeightOptim<F, A, I, N> - for RadonRegTerm<F> +impl<F: Float + ToNalgebraRealField, U, A, I> WeightOptim<F, U, A, I> for RadonRegTerm<F> where + U: Space + Clone + PartialEq, I: AlgIteratorFactory<F>, - A: FindimQuadraticModel<Loc<N, F>, F>, + A: FindimQuadraticModel<U, F>, { fn prepare_optimise_weights(&self, opA: &A, b: &A::Observable) -> DynResult<FindimData<F>> { Ok(FindimData { @@ -192,7 +189,7 @@ fn optimise_weights<'a>( &self, - μ: &mut RNDM<N, F>, + μ: &mut DiscreteMeasure<U, F>, opA: &'a A, b: &A::Observable, findim_data: &FindimData<F>, @@ -221,23 +218,24 @@ } #[replace_float_literals(F::cast_from(literal))] -impl<F: Float + ToNalgebraRealField, A, I, const N: usize> RegTermFW<F, A, I, N> for RadonRegTerm<F> +impl<F: Float + ToNalgebraRealField, U, A, I> RegTermFW<F, U, A, I> for RadonRegTerm<F> where - Cube<N, F>: P2Minimise<Loc<N, F>, F>, + U: Space + Clone + PartialEq, + RadonRegTerm<F>: RegTerm<U, F>, I: AlgIteratorFactory<F>, - A: FindimQuadraticModel<Loc<N, F>, F>, - A::PreadjointCodomain: MinMaxMapping<Loc<N, F>, F>, + A: FindimQuadraticModel<U, F>, + A::PreadjointCodomain: MinMaxMapping<U, F>, for<'a> &'a A::PreadjointCodomain: Instance<A::PreadjointCodomain>, // FIXME: the following *should not* be needed, they are already implied - RNDM<N, F>: Mapping<A::PreadjointCodomain, Codomain = F>, - DeltaMeasure<Loc<N, F>, F>: Mapping<A::PreadjointCodomain, Codomain = F>, + DiscreteMeasure<U, F>: Mapping<A::PreadjointCodomain, Codomain = F>, + DeltaMeasure<U, F>: Copy + Mapping<A::PreadjointCodomain, Codomain = F>, { fn find_insertion( &self, g: &mut A::PreadjointCodomain, refinement_tolerance: F, max_steps: usize, - ) -> (Loc<N, F>, F) { + ) -> (U, F) { let (ξmax, v_ξmax) = g.maximise(refinement_tolerance, max_steps); let (ξmin, v_ξmin) = g.minimise(refinement_tolerance, max_steps); if v_ξmin < 0.0 && -v_ξmin > v_ξmax { @@ -249,10 +247,10 @@ fn relaxed_insert<'a>( &self, - μ: &mut RNDM<N, F>, + μ: &mut DiscreteMeasure<U, F>, g: &A::PreadjointCodomain, opA: &'a A, - ξ: Loc<N, F>, + ξ: U, v_ξ: F, findim_data: &FindimData<F>, ) { @@ -285,11 +283,11 @@ } #[replace_float_literals(F::cast_from(literal))] -impl<F: Float + ToNalgebraRealField, A, I, const N: usize> WeightOptim<F, A, I, N> - for NonnegRadonRegTerm<F> +impl<F: Float + ToNalgebraRealField, A, U, I> WeightOptim<F, U, A, I> for NonnegRadonRegTerm<F> where + U: Space + Clone + PartialEq, I: AlgIteratorFactory<F>, - A: FindimQuadraticModel<Loc<N, F>, F>, + A: FindimQuadraticModel<U, F>, { fn prepare_optimise_weights(&self, opA: &A, b: &A::Observable) -> DynResult<FindimData<F>> { Ok(FindimData { @@ -300,7 +298,7 @@ fn optimise_weights<'a>( &self, - μ: &mut RNDM<N, F>, + μ: &mut DiscreteMeasure<U, F>, opA: &'a A, b: &A::Observable, findim_data: &FindimData<F>, @@ -329,33 +327,33 @@ } #[replace_float_literals(F::cast_from(literal))] -impl<F: Float + ToNalgebraRealField, A, I, const N: usize> RegTermFW<F, A, I, N> - for NonnegRadonRegTerm<F> +impl<F: Float + ToNalgebraRealField, U, A, I> RegTermFW<F, U, A, I> for NonnegRadonRegTerm<F> where - Cube<N, F>: P2Minimise<Loc<N, F>, F>, + U: Space + Clone + PartialEq, + NonnegRadonRegTerm<F>: RegTerm<U, F>, I: AlgIteratorFactory<F>, - A: FindimQuadraticModel<Loc<N, F>, F>, - A::PreadjointCodomain: MinMaxMapping<Loc<N, F>, F>, + A: FindimQuadraticModel<U, F>, + A::PreadjointCodomain: MinMaxMapping<U, F>, for<'a> &'a A::PreadjointCodomain: Instance<A::PreadjointCodomain>, // FIXME: the following *should not* be needed, they are already implied - RNDM<N, F>: Mapping<A::PreadjointCodomain, Codomain = F>, - DeltaMeasure<Loc<N, F>, F>: Mapping<A::PreadjointCodomain, Codomain = F>, + DiscreteMeasure<U, F>: Mapping<A::PreadjointCodomain, Codomain = F>, + DeltaMeasure<U, F>: Mapping<A::PreadjointCodomain, Codomain = F>, { fn find_insertion( &self, g: &mut A::PreadjointCodomain, refinement_tolerance: F, max_steps: usize, - ) -> (Loc<N, F>, F) { + ) -> (U, F) { g.maximise(refinement_tolerance, max_steps) } fn relaxed_insert<'a>( &self, - μ: &mut RNDM<N, F>, + μ: &mut DiscreteMeasure<U, F>, g: &A::PreadjointCodomain, opA: &'a A, - ξ: Loc<N, F>, + ξ: U, v_ξ: F, findim_data: &FindimData<F>, ) { @@ -401,25 +399,25 @@ /// `iterator` is used to iterate the steps of the method, and `plotter` may be used to /// save intermediate iteration states as images. #[replace_float_literals(F::cast_from(literal))] -pub fn pointsource_fw_reg<'a, F, I, A, Reg, Plot, const N: usize>( - f: &'a QuadraticDataTerm<F, RNDM<N, F>, A>, +pub fn pointsource_fw_reg<'a, F, U, I, A, Reg, Plot>( + f: &'a QuadraticDataTerm<F, DiscreteMeasure<U, F>, A>, reg: &Reg, - //domain : Cube<N, F>, config: &FWConfig<F>, iterator: I, mut plotter: Plot, - μ0: Option<RNDM<N, F>>, -) -> DynResult<RNDM<N, F>> + μ0: Option<DiscreteMeasure<U, F>>, +) -> DynResult<DiscreteMeasure<U, F>> where F: Float + ToNalgebraRealField, + U: Space + Clone + PartialEq, I: AlgIteratorFactory<IterInfo<F>>, - A: ForwardModel<RNDM<N, F>, F>, - A::PreadjointCodomain: MinMaxMapping<Loc<N, F>, F>, + A: ForwardModel<DiscreteMeasure<U, F>, F>, + A::PreadjointCodomain: MinMaxMapping<U, F>, &'a A::PreadjointCodomain: Instance<A::PreadjointCodomain>, - Cube<N, F>: P2Minimise<Loc<N, F>, F>, - RNDM<N, F>: SpikeMerging<F>, - Reg: RegTermFW<F, A, ValueIteratorFactory<F, AlgIteratorOptions>, N>, - Plot: Plotter<A::PreadjointCodomain, A::PreadjointCodomain, RNDM<N, F>>, + //Cube<N, F>: P2Minimise<U, F>, + DiscreteMeasure<U, F>: SpikeMerging<F>, + Reg: RegTermFW<F, U, A, ValueIteratorFactory<F, AlgIteratorOptions>>, + Plot: Plotter<A::PreadjointCodomain, A::PreadjointCodomain, DiscreteMeasure<U, F>>, { let opA = f.operator(); let b = f.data(); @@ -438,7 +436,7 @@ let mut residual = f.residual(&μ); // Statistics - let full_stats = |residual: &A::Observable, ν: &RNDM<N, F>, ε, stats| IterInfo { + let full_stats = |residual: &A::Observable, ν: &DiscreteMeasure<U, F>, ε, stats| IterInfo { value: residual.norm2_squared_div2() + reg.apply(ν), n_spikes: ν.len(), ε,
--- a/src/pdps.rs Tue Aug 04 20:11:10 2026 -0500 +++ b/src/pdps.rs Thu Aug 06 11:14:15 2026 -0500 @@ -42,7 +42,7 @@ use crate::forward_model::ForwardModel; use crate::measures::merging::SpikeMerging; use crate::measures::merging::SpikeMergingMethod; -use crate::measures::{DiscreteMeasure, RNDM}; +use crate::measures::DiscreteMeasure; use crate::plot::Plotter; pub use crate::prox_penalty::{InsertionConfig, ProxPenalty, StepLengthBoundPD}; use crate::regularisation::RegTerm; @@ -50,7 +50,7 @@ use alg_tools::convex::{Conjugable, ConvexMapping, Prox}; use alg_tools::error::DynResult; use alg_tools::iterate::AlgIteratorFactory; -use alg_tools::linops::{Mapping, AXPY}; +use alg_tools::linops::{Mapping, Space, AXPY}; use alg_tools::mapping::{DataTerm, Instance}; use alg_tools::nalgebra_support::ToNalgebraRealField; use anyhow::ensure; @@ -140,27 +140,29 @@ /// /// Returns the final iterate. #[replace_float_literals(F::cast_from(literal))] -pub fn pointsource_pdps_reg<'a, F, I, A, Phi, Reg, Plot, P, const N: usize>( - f: &'a DataTerm<F, RNDM<N, F>, A, Phi>, +pub fn pointsource_pdps_reg<'a, F, U, I, A, Phi, Reg, Plot, P>( + f: &'a DataTerm<F, DiscreteMeasure<U, F>, A, Phi>, reg: &Reg, prox_penalty: &P, pdpsconfig: &PDPSConfig<F>, iterator: I, mut plotter: Plot, - μ0: Option<RNDM<N, F>>, -) -> DynResult<RNDM<N, F>> + μ0: Option<DiscreteMeasure<U, F>>, +) -> DynResult<DiscreteMeasure<U, F>> where + U: Space + Clone + PartialEq, F: Float + ToNalgebraRealField, I: AlgIteratorFactory<IterInfo<F>>, - A: ForwardModel<RNDM<N, F>, F>, + A: ForwardModel<DiscreteMeasure<U, F>, F>, for<'b> &'b A::Observable: Instance<A::Observable>, A::Observable: AXPY<Field = F>, - RNDM<N, F>: SpikeMerging<F>, - Reg: RegTerm<Loc<N, F>, F>, + DiscreteMeasure<U, F>: SpikeMerging<F>, + Reg: RegTerm<U, F>, Phi: Conjugable<A::Observable, F>, for<'b> Phi::Conjugate<'b>: Prox<A::Observable>, - P: ProxPenalty<Loc<N, F>, A::PreadjointCodomain, Reg, F> + StepLengthBoundPD<F, A, RNDM<N, F>>, - Plot: Plotter<P::ReturnMapping, A::PreadjointCodomain, RNDM<N, F>>, + P: ProxPenalty<U, A::PreadjointCodomain, Reg, F> + + StepLengthBoundPD<F, A, DiscreteMeasure<U, F>>, + Plot: Plotter<P::ReturnMapping, A::PreadjointCodomain, DiscreteMeasure<U, F>>, { // Check parameters ensure!( @@ -187,7 +189,7 @@ // Initialise iterates let mut μ = μ0.unwrap_or_else(|| DiscreteMeasure::new()); let mut y = f.residual(&μ); - let full_stats = |μ: &RNDM<N, F>, ε, stats| IterInfo { + let full_stats = |μ: &DiscreteMeasure<U, F>, ε, stats| IterInfo { value: f.apply(μ) + reg.apply(μ), n_spikes: μ.len(), ε,
--- a/src/run.rs Tue Aug 04 20:11:10 2026 -0500 +++ b/src/run.rs Thu Aug 06 11:14:15 2026 -0500 @@ -55,7 +55,7 @@ }; use alg_tools::maputil::map3; use alg_tools::nalgebra_support::ToNalgebraRealField; -use alg_tools::norms::{NormExponent, L1, L2}; +use alg_tools::norms::{Dist, NormExponent, L1, L2}; use alg_tools::operator_arithmetic::Weighted; use alg_tools::sets::Cube; use alg_tools::sets::SetOrd; @@ -1086,30 +1086,36 @@ /// $L^2$-squared data term is enabled through the type signatures. /// /// `cont` gives a continuation attempt to find the algorithm matching the description `alg`. -pub fn run_fw<'a, F, A, Reg, I, Plot, const N: usize>( - f: &'a DataTerm<F, RNDM<N, F>, A, Norm222<F>>, +pub fn run_fw<'a, F, U, A, Reg, I, Plot>( + f: &'a DataTerm<F, DiscreteMeasure<U, F>, A, Norm222<F>>, reg: &Reg, (alg, iterator, plotter, μ0, running): ( &AlgorithmConfig<F>, I, Plot, - Option<RNDM<N, F>>, + Option<DiscreteMeasure<U, F>>, String, ), cont: impl FnOnce( - (&AlgorithmConfig<F>, I, Plot, Option<RNDM<N, F>>, String), - ) -> DynResult<RNDM<N, F>>, -) -> DynResult<RNDM<N, F>> + ( + &AlgorithmConfig<F>, + I, + Plot, + Option<DiscreteMeasure<U, F>>, + String, + ), + ) -> DynResult<DiscreteMeasure<U, F>>, +) -> DynResult<DiscreteMeasure<U, F>> where F: Float + ToNalgebraRealField, + U: Space + Clone + PartialEq, I: AlgIteratorFactory<IterInfo<F>>, - A: ForwardModel<RNDM<N, F>, F>, - A::PreadjointCodomain: MinMaxMapping<Loc<N, F>, F>, + A: ForwardModel<DiscreteMeasure<U, F>, F>, + A::PreadjointCodomain: MinMaxMapping<U, F>, for<'b> &'b A::PreadjointCodomain: Instance<A::PreadjointCodomain>, - Cube<N, F>: P2Minimise<Loc<N, F>, F>, - RNDM<N, F>: SpikeMerging<F>, - Reg: RegTermFW<F, A, ValueIteratorFactory<F, AlgIteratorOptions>, N>, - Plot: Plotter<A::PreadjointCodomain, A::PreadjointCodomain, RNDM<N, F>>, + DiscreteMeasure<U, F>: SpikeMerging<F>, + Reg: RegTermFW<F, U, A, ValueIteratorFactory<F, AlgIteratorOptions>>, + Plot: Plotter<A::PreadjointCodomain, A::PreadjointCodomain, DiscreteMeasure<U, F>>, { match alg { &AlgorithmConfig::FW(ref algconfig) => { @@ -1381,10 +1387,10 @@ } } -type MeasureZ<F, Z, const N: usize> = Pair<RNDM<N, F>, Z>; +type MeasureZ<F, Z, U> = Pair<DiscreteMeasure<U, F>, Z>; #[deprecated] -pub fn run_pdps_pair<F, S, Dat, Reg, Z, R, Y, KOpZ, H, P, I, Plot, const N: usize>( +pub fn run_pdps_pair<F, S, Dat, Reg, Z, R, Y, KOpZ, H, P, I, Plot, U>( f: &Dat, reg: &Reg, prox_penalty: &P, @@ -1395,7 +1401,7 @@ &AlgorithmConfig<F>, I, Plot, - (Option<RNDM<N, F>>, Z, Y), + (Option<DiscreteMeasure<U, F>>, Z, Y), String, ), cont: impl FnOnce( @@ -1403,22 +1409,24 @@ &AlgorithmConfig<F>, I, Plot, - (Option<RNDM<N, F>>, Z, Y), + (Option<DiscreteMeasure<U, F>>, Z, Y), String, ), - ) -> DynResult<Pair<RNDM<N, F>, Z>>, -) -> DynResult<Pair<RNDM<N, F>, Z>> + ) -> DynResult<Pair<DiscreteMeasure<U, F>, Z>>, +) -> DynResult<Pair<DiscreteMeasure<U, F>, Z>> where F: Float + ToNalgebraRealField, + U: Space + Clone + Dist<L2, F>, + for<'a> &'a U: Instance<U>, I: AlgIteratorFactory<IterInfo<F>>, - Dat: DifferentiableMapping<MeasureZ<F, Z, N>, Codomain = F, DerivativeDomain = Pair<S, Z>> + Dat: DifferentiableMapping<MeasureZ<F, Z, U>, Codomain = F, DerivativeDomain = Pair<S, Z>> + BoundedCurvature<F>, - S: DifferentiableRealMapping<N, F> + ClosedMul<F> + MinMaxMapping<Loc<N, F>, F>, + S: DifferentiableMapping<U, Codomain = F> + ClosedMul<F> + MinMaxMapping<U, F>, for<'a> Pair<&'a P, &'a IdOp<Z>>: StepLengthBoundPair<F, Dat>, //Pair<S, Z>: ClosedMul<F>, - RNDM<N, F>: SpikeMerging<F>, - Reg: SlidingRegTerm<Loc<N, F>, F>, - P: TransportProxPenalty<Loc<N, F>, S, Reg, F>, + DiscreteMeasure<U, F>: SpikeMerging<F>, + Reg: SlidingRegTerm<U, F>, + P: TransportProxPenalty<U, S, Reg, F>, KOpZ: BoundedLinear<Z, L2, L2, F, Codomain = Y> + GEMV<F, Z> + SimplyAdjointable<Z, Y, AdjointCodomain = Z>, @@ -1430,7 +1438,7 @@ R: Prox<Z, Codomain = F>, H: Conjugable<Y, F, Codomain = F>, for<'b> H::Conjugate<'b>: Prox<Y>, - Plot: Plotter<P::ReturnMapping, S, RNDM<N, F>>, + Plot: Plotter<P::ReturnMapping, S, DiscreteMeasure<U, F>>, { run_pdps_pair_dom( f, @@ -1445,7 +1453,7 @@ ) } -pub fn run_pdps_pair_dom<F, S, Dat, Reg, Z, R, Y, KOpZ, H, P, I, Plot, const N: usize>( +pub fn run_pdps_pair_dom<F, S, Dat, Reg, Z, R, Y, KOpZ, H, P, I, Plot, U>( f: &Dat, reg: &Reg, prox_penalty: &P, @@ -1456,31 +1464,33 @@ &AlgorithmConfig<F>, I, Plot, - (Option<RNDM<N, F>>, Z, Y), + (Option<DiscreteMeasure<U, F>>, Z, Y), String, ), - domain: &mut impl AdmissibleDomain<Loc<N, F>>, + domain: &mut impl AdmissibleDomain<U>, cont: impl FnOnce( ( &AlgorithmConfig<F>, I, Plot, - (Option<RNDM<N, F>>, Z, Y), + (Option<DiscreteMeasure<U, F>>, Z, Y), String, ), - ) -> DynResult<Pair<RNDM<N, F>, Z>>, -) -> DynResult<Pair<RNDM<N, F>, Z>> + ) -> DynResult<Pair<DiscreteMeasure<U, F>, Z>>, +) -> DynResult<Pair<DiscreteMeasure<U, F>, Z>> where F: Float + ToNalgebraRealField, + U: Space + Clone + Dist<L2, F>, + for<'a> &'a U: Instance<U>, I: AlgIteratorFactory<IterInfo<F>>, - Dat: DifferentiableMapping<MeasureZ<F, Z, N>, Codomain = F, DerivativeDomain = Pair<S, Z>> + Dat: DifferentiableMapping<MeasureZ<F, Z, U>, Codomain = F, DerivativeDomain = Pair<S, Z>> + BoundedCurvature<F>, - S: DifferentiableRealMapping<N, F> + ClosedMul<F> + MinMaxMapping<Loc<N, F>, F>, + S: DifferentiableMapping<U, Codomain = F> + ClosedMul<F> + MinMaxMapping<U, F>, for<'a> Pair<&'a P, &'a IdOp<Z>>: StepLengthBoundPair<F, Dat>, //Pair<S, Z>: ClosedMul<F>, - RNDM<N, F>: SpikeMerging<F>, - Reg: SlidingRegTerm<Loc<N, F>, F>, - P: TransportProxPenalty<Loc<N, F>, S, Reg, F>, + DiscreteMeasure<U, F>: SpikeMerging<F>, + Reg: SlidingRegTerm<U, F>, + P: TransportProxPenalty<U, S, Reg, F>, KOpZ: BoundedLinear<Z, L2, L2, F, Codomain = Y> + GEMV<F, Z> + SimplyAdjointable<Z, Y, AdjointCodomain = Z>, @@ -1492,7 +1502,7 @@ R: Prox<Z, Codomain = F>, H: Conjugable<Y, F, Codomain = F>, for<'b> H::Conjugate<'b>: Prox<Y>, - Plot: Plotter<P::ReturnMapping, S, RNDM<N, F>>, + Plot: Plotter<P::ReturnMapping, S, DiscreteMeasure<U, F>>, { let pt = P::prox_type(); @@ -1533,7 +1543,7 @@ } } -pub fn run_fb_pair<F, S, Dat, Reg, Z, R, P, I, Plot, const N: usize>( +pub fn run_fb_pair<F, S, Dat, Reg, Z, R, P, I, Plot, U>( f: &Dat, reg: &Reg, prox_penalty: &P, @@ -1542,7 +1552,7 @@ &AlgorithmConfig<F>, I, Plot, - (Option<RNDM<N, F>>, Z), + (Option<DiscreteMeasure<U, F>>, Z), String, ), cont: impl FnOnce( @@ -1550,32 +1560,34 @@ &AlgorithmConfig<F>, I, Plot, - (Option<RNDM<N, F>>, Z), + (Option<DiscreteMeasure<U, F>>, Z), String, ), - ) -> DynResult<Pair<RNDM<N, F>, Z>>, -) -> DynResult<Pair<RNDM<N, F>, Z>> + ) -> DynResult<Pair<DiscreteMeasure<U, F>, Z>>, +) -> DynResult<Pair<DiscreteMeasure<U, F>, Z>> where F: Float + ToNalgebraRealField, + U: Space + Clone + Dist<L2, F>, + for<'a> &'a U: Instance<U>, I: AlgIteratorFactory<IterInfo<F>>, - Dat: DifferentiableMapping<MeasureZ<F, Z, N>, Codomain = F, DerivativeDomain = Pair<S, Z>> + Dat: DifferentiableMapping<MeasureZ<F, Z, U>, Codomain = F, DerivativeDomain = Pair<S, Z>> + BoundedCurvature<F>, - S: DifferentiableRealMapping<N, F> + ClosedMul<F> + MinMaxMapping<Loc<N, F>, F>, - RNDM<N, F>: SpikeMerging<F>, - Reg: SlidingRegTerm<Loc<N, F>, F>, - P: TransportProxPenalty<Loc<N, F>, S, Reg, F>, + S: DifferentiableMapping<U, Codomain = F> + ClosedMul<F> + MinMaxMapping<U, F>, + DiscreteMeasure<U, F>: SpikeMerging<F>, + Reg: SlidingRegTerm<U, F>, + P: TransportProxPenalty<U, S, Reg, F>, for<'a> Pair<&'a P, &'a IdOp<Z>>: StepLengthBoundPair<F, Dat>, Z: ClosedEuclidean<F> + AXPY + Clone, for<'b> &'b Z: Instance<Z>, R: Prox<Z, Codomain = F>, - Plot: Plotter<P::ReturnMapping, S, RNDM<N, F>>, + Plot: Plotter<P::ReturnMapping, S, DiscreteMeasure<U, F>>, // We should not need to explicitly require this: for<'b> &'b Loc<0, F>: Instance<Loc<0, F>>, { run_fb_pair_dom(f, reg, prox_penalty, fnR, init, &mut FullSpace, cont) } -pub fn run_fb_pair_dom<F, S, Dat, Reg, Z, R, P, I, Plot, const N: usize>( +pub fn run_fb_pair_dom<F, S, Dat, Reg, Z, R, P, I, Plot, U>( f: &Dat, reg: &Reg, prox_penalty: &P, @@ -1584,34 +1596,36 @@ &AlgorithmConfig<F>, I, Plot, - (Option<RNDM<N, F>>, Z), + (Option<DiscreteMeasure<U, F>>, Z), String, ), - domain: &mut impl AdmissibleDomain<Loc<N, F>>, + domain: &mut impl AdmissibleDomain<U>, cont: impl FnOnce( ( &AlgorithmConfig<F>, I, Plot, - (Option<RNDM<N, F>>, Z), + (Option<DiscreteMeasure<U, F>>, Z), String, ), - ) -> DynResult<Pair<RNDM<N, F>, Z>>, -) -> DynResult<Pair<RNDM<N, F>, Z>> + ) -> DynResult<Pair<DiscreteMeasure<U, F>, Z>>, +) -> DynResult<Pair<DiscreteMeasure<U, F>, Z>> where F: Float + ToNalgebraRealField, + U: Space + Clone + Dist<L2, F>, + for<'a> &'a U: Instance<U>, I: AlgIteratorFactory<IterInfo<F>>, - Dat: DifferentiableMapping<MeasureZ<F, Z, N>, Codomain = F, DerivativeDomain = Pair<S, Z>> + Dat: DifferentiableMapping<MeasureZ<F, Z, U>, Codomain = F, DerivativeDomain = Pair<S, Z>> + BoundedCurvature<F>, - S: DifferentiableRealMapping<N, F> + ClosedMul<F> + MinMaxMapping<Loc<N, F>, F>, - RNDM<N, F>: SpikeMerging<F>, - Reg: SlidingRegTerm<Loc<N, F>, F>, - P: TransportProxPenalty<Loc<N, F>, S, Reg, F>, + S: DifferentiableMapping<U, Codomain = F> + ClosedMul<F> + MinMaxMapping<U, F>, + DiscreteMeasure<U, F>: SpikeMerging<F>, + Reg: SlidingRegTerm<U, F>, + P: TransportProxPenalty<U, S, Reg, F>, for<'a> Pair<&'a P, &'a IdOp<Z>>: StepLengthBoundPair<F, Dat>, Z: ClosedEuclidean<F> + AXPY + Clone, for<'b> &'b Z: Instance<Z>, R: Prox<Z, Codomain = F>, - Plot: Plotter<P::ReturnMapping, S, RNDM<N, F>>, + Plot: Plotter<P::ReturnMapping, S, DiscreteMeasure<U, F>>, // We should not need to explicitly require this: for<'b> &'b Loc<0, F>: Instance<Loc<0, F>>, {
--- a/src/sliding_fb.rs Tue Aug 04 20:11:10 2026 -0500 +++ b/src/sliding_fb.rs Thu Aug 06 11:14:15 2026 -0500 @@ -24,11 +24,11 @@ use alg_tools::bounds::{Bounds, MinMaxMapping}; use alg_tools::error::DynResult; use alg_tools::euclidean::Euclidean; -use alg_tools::instance::Space; +use alg_tools::instance::{Instance, Space}; use alg_tools::iterate::AlgIteratorFactory; use alg_tools::mapping::{DifferentiableMapping, DifferentiableRealMapping, Mapping, RealMapping}; use alg_tools::nalgebra_support::ToNalgebraRealField; -use alg_tools::norms::Norm; +use alg_tools::norms::{Dist, Norm, L2}; use anyhow::ensure; use std::ops::ControlFlow; @@ -193,11 +193,11 @@ } #[derive(Clone, Debug, Serialize)] -pub struct SingleTransport<Domain, F: Float> { +pub struct SingleTransport<U, F: Float> { /// Source point - x: Domain, + x: U, /// Target point - y: Domain, + y: U, /// Original mass α_μ_orig: F, /// Transported mass @@ -211,8 +211,8 @@ } #[derive(Clone, Debug, Serialize)] -pub struct Transport<Domain, F: Float> { - vec: Vec<SingleTransport<Domain, F>>, +pub struct Transport<U, F: Float> { + vec: Vec<SingleTransport<U, F>>, } /// Whether partially transported points are allowed. @@ -223,21 +223,21 @@ const MINIMAL_PARTIAL_TRANSPORT: bool = true; const NEW_APPROACH: bool = true; -pub trait TransportProxPenalty<Point, PreadjointCodomain, Reg, F = f64>: - ProxPenalty<Point, PreadjointCodomain, Reg, F> +pub trait TransportProxPenalty<U, PreadjointCodomain, Reg, F = f64>: + ProxPenalty<U, PreadjointCodomain, Reg, F> where F: Float + ToNalgebraRealField, - Reg: SlidingRegTerm<Point, F>, - Point: Space + Clone, + Reg: SlidingRegTerm<U, F>, + U: Space + Clone, { type TransportStepLength; /// Constrution of initial transport `γ1` from initial measure `μ` and `v=F'(μ)` /// with step lengh τ and transport step length `θ_or_adaptive`. - fn initial_transport<Dom: AdmissibleDomain<Point>>( + fn initial_transport<Dom: AdmissibleDomain<U>>( &self, - γ: &mut Transport<Point, F>, - μ: &DiscreteMeasure<Point, F>, + γ: &mut Transport<U, F>, + μ: &DiscreteMeasure<U, F>, ε: F, τ: F, τθ_or_adaptive: &mut Self::TransportStepLength, @@ -249,9 +249,9 @@ /// A posteriori transport adaptation. fn aposteriori_transport( &self, - γ: &mut Transport<Point, F>, - μ: &DiscreteMeasure<Point, F>, - μ̆: &DiscreteMeasure<Point, F>, + γ: &mut Transport<U, F>, + μ: &DiscreteMeasure<U, F>, + μ̆: &DiscreteMeasure<U, F>, τv̆: &mut PreadjointCodomain, v: &mut PreadjointCodomain, extra: Option<F>, @@ -377,14 +377,12 @@ } #[replace_float_literals(F::cast_from(literal))] -impl<F, M, Reg, 𝒟, O, const N: usize> TransportProxPenalty<Loc<N, F>, M, Reg, F> for 𝒟 +impl<F, M, Reg, 𝒟, const N: usize> TransportProxPenalty<Loc<N, F>, M, Reg, F> for 𝒟 where F: Float + ToNalgebraRealField, 𝒟: DiscreteMeasureOp<Loc<N, F>, F>, 𝒟::Codomain: RealMapping<N, F>, M: MinMaxMapping<Loc<N, F>, F> + DifferentiableRealMapping<N, F>, - for<'a> &'a M: std::ops::Add<𝒟::PreCodomain, Output = O>, - O: MinMaxMapping<Loc<N, F>, F>, Reg: SlidingRegTerm<Loc<N, F>, F>, Self: ProxPenalty<Loc<N, F>, M, Reg, F>, RNDM<N, F>: SpikeMerging<F>, @@ -608,25 +606,37 @@ } } +impl<U, F: Float> Norm<Radon, F> for Transport<U, F> { + fn norm(&self, _: Radon) -> F { + self.iter().map(|ρ| ρ.α_γ.abs()).sum() + } +} + +impl<U, F: Float> MulAssign<F> for Transport<U, F> { + fn mul_assign(&mut self, factor: F) { + for ρ in self.iter_mut() { + ρ.α_γ *= factor; + } + } +} + #[replace_float_literals(F::cast_from(literal))] -impl<const N: usize, F: Float> Transport<Loc<N, F>, F> { +impl<U, F: Float> Transport<U, F> { pub(crate) fn new() -> Self { Transport { vec: Vec::new() } } - pub(crate) fn iter(&self) -> impl Iterator<Item = &'_ SingleTransport<Loc<N, F>, F>> { + pub(crate) fn iter(&self) -> impl Iterator<Item = &'_ SingleTransport<U, F>> { self.vec.iter() } - pub(crate) fn iter_mut( - &mut self, - ) -> impl Iterator<Item = &'_ mut SingleTransport<Loc<N, F>, F>> { + pub(crate) fn iter_mut(&mut self) -> impl Iterator<Item = &'_ mut SingleTransport<U, F>> { self.vec.iter_mut() } pub(crate) fn extend<I>(&mut self, it: I) where - I: IntoIterator<Item = SingleTransport<Loc<N, F>, F>>, + I: IntoIterator<Item = SingleTransport<U, F>>, { self.vec.extend(it) } @@ -635,62 +645,44 @@ self.vec.len() } - // pub(crate) fn dist_matching(&self, μ: &RNDM<N, F>) -> F { - // self.iter() - // .zip(μ.iter_spikes()) - // .map(|(ρ, δ)| (ρ.α_γ - δ.α).abs()) - // .sum() - // } - - /// Construct `μ̆`, replacing the contents of `μ`. - pub(crate) fn μ̆_into(&self, μ: &mut RNDM<N, F>) { - assert!(self.len() <= μ.len()); - - // First transported points - for (δ, ρ) in izip!(μ.iter_spikes_mut(), self.iter()) { - if ρ.α_γ.abs() > 0.0 { - // Transport – transported point - δ.α = ρ.α_γ; - δ.x = ρ.y; - } else { - // No transport – original point - δ.α = ρ.α_μ_orig; - δ.x = ρ.x; - } - } - - // Then source points with partial transport - let mut i = self.len(); - // This can cause the number of points to explode, so cannot have partial transport. - for ρ in self.iter() { - let α = ρ.α_μ_orig - ρ.α_γ; - if ρ.α_γ.abs() > F::EPSILON && α.abs() > F::EPSILON { - let δ = DeltaMeasure { α, x: ρ.x }; - if i < μ.len() { - μ[i] = δ; - } else { - μ.push(δ) - } - i += 1; - } - } - μ.truncate(i); - } - /// Returns $‖μ\^k - π\_♯\^0γ\^{k+1}‖$ pub(crate) fn μ0_minus_γ0_radon(&self) -> F { self.vec.iter().map(|ρ| (ρ.α_μ_orig - ρ.α_γ).abs()).sum() } + /// Prune spikes with zero weight. To maintain correct ordering between μ and γ, also the + /// latter needs to be pruned when μ is. + pub(crate) fn prune_compat(&mut self, μ: &mut DiscreteMeasure<U, F>, stats: &mut IterInfo<F>) { + assert!(self.vec.len() <= μ.len()); + let old_len = μ.len(); + for (ρ, δ) in self.vec.iter_mut().zip(μ.iter_spikes()) { + ρ.retain = δ.α != F::ZERO; + } + // Do not prune new insertions on first round. + μ.prune_until(self.vec.len()); + stats.pruned += old_len - μ.len(); + self.vec.retain(|ρ| ρ.retain); + assert!(self.vec.len() <= μ.len()); + } +} + +#[replace_float_literals(F::cast_from(literal))] +impl<U, F> Transport<U, F> +where + U: Dist<L2, F>, + F: Float, + for<'a> &'a U: Instance<U>, +{ /// Returns $∫ c_2 d|γ|$ pub(crate) fn c2integral(&self) -> F { + let sq2 = |x| (x * x) / 2.0; self.vec .iter() - .map(|ρ| ρ.y.dist2_squared(&ρ.x) / 2.0 * ρ.α_γ.abs()) + .map(|ρ| sq2(ρ.y.dist(&ρ.x, L2)) * ρ.α_γ.abs()) .sum() } - pub(crate) fn get_transport_stats(&self, stats: &mut IterInfo<F>, μ: &RNDM<N, F>) { + pub(crate) fn get_transport_stats(&self, stats: &mut IterInfo<F>, μ: &DiscreteMeasure<U, F>) { // TODO: This doesn't take into account μ[i].α becoming zero in the latest tranport // attempt, for i < self.len(), when a corresponding source term also exists with index // j ≥ self.len(). For now, we let that be reflected in the prune count. @@ -726,22 +718,55 @@ }) }; } +} - /// Prune spikes with zero weight. To maintain correct ordering between μ and γ, also the - /// latter needs to be pruned when μ is. - pub(crate) fn prune_compat(&mut self, μ: &mut RNDM<N, F>, stats: &mut IterInfo<F>) { - assert!(self.vec.len() <= μ.len()); - let old_len = μ.len(); - for (ρ, δ) in self.vec.iter_mut().zip(μ.iter_spikes()) { - ρ.retain = δ.α != F::ZERO; +#[replace_float_literals(F::cast_from(literal))] +impl<U: Clone, F: Float> Transport<U, F> { + // pub(crate) fn dist_matching(&self, μ: &RNDM<N, F>) -> F { + // self.iter() + // .zip(μ.iter_spikes()) + // .map(|(ρ, δ)| (ρ.α_γ - δ.α).abs()) + // .sum() + // } + + /// Construct `μ̆`, replacing the contents of `μ`. + pub(crate) fn μ̆_into(&self, μ: &mut DiscreteMeasure<U, F>) { + assert!(self.len() <= μ.len()); + + // First transported points + for (δ, ρ) in izip!(μ.iter_spikes_mut(), self.iter()) { + if ρ.α_γ.abs() > 0.0 { + // Transport – transported point + δ.α = ρ.α_γ; + δ.x = ρ.y.clone(); + } else { + // No transport – original point + δ.α = ρ.α_μ_orig; + δ.x = ρ.x.clone(); + } } - // Do not prune new insertions on first round. - μ.prune_until(self.vec.len()); - stats.pruned += old_len - μ.len(); - self.vec.retain(|ρ| ρ.retain); - assert!(self.vec.len() <= μ.len()); + + // Then source points with partial transport + let mut i = self.len(); + // This can cause the number of points to explode, so cannot have partial transport. + for ρ in self.iter() { + let α = ρ.α_μ_orig - ρ.α_γ; + if ρ.α_γ.abs() > F::EPSILON && α.abs() > F::EPSILON { + let δ = DeltaMeasure { α, x: ρ.x.clone() }; + if i < μ.len() { + μ[i] = δ; + } else { + μ.push(δ) + } + i += 1; + } + } + μ.truncate(i); } +} +#[replace_float_literals(F::cast_from(literal))] +impl<const N: usize, F: Float> Transport<Loc<N, F>, F> { /// Helper for initial transport. Called from [`TransportProxPenalty::initial_transport`]. pub(super) fn do_init_transport<D, Dom: AdmissibleDomain<Loc<N, F>>>( &mut self, @@ -955,20 +980,6 @@ } } -impl<const N: usize, F: Float> Norm<Radon, F> for Transport<Loc<N, F>, F> { - fn norm(&self, _: Radon) -> F { - self.iter().map(|ρ| ρ.α_γ.abs()).sum() - } -} - -impl<const N: usize, F: Float> MulAssign<F> for Transport<Loc<N, F>, F> { - fn mul_assign(&mut self, factor: F) { - for ρ in self.iter_mut() { - ρ.α_γ *= factor; - } - } -} - /// Iteratively solve the pointsource localisation problem using sliding forward-backward /// splitting /// @@ -976,27 +987,29 @@ /// may be given. /// Inertia is currently not supported. #[replace_float_literals(F::cast_from(literal))] -pub fn pointsource_sliding_fb<F, I, Dat, Reg, Plot, P, const N: usize>( +pub fn pointsource_sliding_fb<F, U, I, Dat, Reg, Plot, P>( f: &Dat, reg: &Reg, prox_penalty: &P, config: &SlidingFBConfig<F>, iterator: I, mut plotter: Plot, - μ0: Option<RNDM<N, F>>, - domain: &mut impl AdmissibleDomain<Loc<N, F>>, -) -> DynResult<RNDM<N, F>> + μ0: Option<DiscreteMeasure<U, F>>, + domain: &mut impl AdmissibleDomain<U>, +) -> DynResult<DiscreteMeasure<U, F>> where F: Float + ToNalgebraRealField, + U: Space + Clone + Dist<L2, F>, + for<'a> &'a U: Instance<U>, I: AlgIteratorFactory<IterInfo<F>>, - Dat: DifferentiableMapping<RNDM<N, F>, Codomain = F> + BoundedCurvature<F>, + Dat: DifferentiableMapping<DiscreteMeasure<U, F>, Codomain = F> + BoundedCurvature<F>, Dat::DerivativeDomain: - DifferentiableRealMapping<N, F> + ClosedMul<F> + MinMaxMapping<Loc<N, F>, F>, + DifferentiableMapping<U, Codomain = F> + ClosedMul<F> + MinMaxMapping<U, F>, //for<'a> Dat::Differential<'a>: Lipschitz<&'a P, FloatType = F>, - RNDM<N, F>: SpikeMerging<F>, - Reg: SlidingRegTerm<Loc<N, F>, F>, - P: TransportProxPenalty<Loc<N, F>, Dat::DerivativeDomain, Reg, F> + StepLengthBound<F, Dat>, - Plot: Plotter<P::ReturnMapping, Dat::DerivativeDomain, RNDM<N, F>>, + DiscreteMeasure<U, F>: SpikeMerging<F>, + Reg: SlidingRegTerm<U, F>, + P: TransportProxPenalty<U, Dat::DerivativeDomain, Reg, F> + StepLengthBound<F, Dat>, + Plot: Plotter<P::ReturnMapping, Dat::DerivativeDomain, DiscreteMeasure<U, F>>, { // Check parameters ensure!(config.τ0 > 0.0, "Invalid step length parameter"); @@ -1025,7 +1038,7 @@ let mut ε = tolerance.initial(); // Statistics - let full_stats = |μ: &RNDM<N, F>, ε, stats| IterInfo { + let full_stats = |μ: &DiscreteMeasure<U, F>, ε, stats| IterInfo { value: f.apply(μ) + reg.apply(μ), n_spikes: μ.len(), ε, @@ -1114,7 +1127,7 @@ ε, ins, ®, - Some(|μ̃: &RNDM<N, F>| f.apply(μ̃) + reg.apply(μ̃)), + Some(|μ̃: &DiscreteMeasure<U, F>| f.apply(μ̃) + reg.apply(μ̃)), ); // Prune. Do not prune new insertions on first round. γ.prune_compat(&mut μ, &mut stats); @@ -1167,26 +1180,28 @@ #[deprecated] /// Version of `pointsource_sliding_fb` without domain restriction -pub fn pointsource_sliding_fb_reg<F, I, Dat, Reg, Plot, P, const N: usize>( +pub fn pointsource_sliding_fb_reg<F, U, I, Dat, Reg, Plot, P>( f: &Dat, reg: &Reg, prox_penalty: &P, config: &SlidingFBConfig<F>, iterator: I, plotter: Plot, - μ0: Option<RNDM<N, F>>, -) -> DynResult<RNDM<N, F>> + μ0: Option<DiscreteMeasure<U, F>>, +) -> DynResult<DiscreteMeasure<U, F>> where F: Float + ToNalgebraRealField, + U: Space + Clone + Dist<L2, F>, + for<'a> &'a U: Instance<U>, I: AlgIteratorFactory<IterInfo<F>>, - Dat: DifferentiableMapping<RNDM<N, F>, Codomain = F> + BoundedCurvature<F>, + Dat: DifferentiableMapping<DiscreteMeasure<U, F>, Codomain = F> + BoundedCurvature<F>, Dat::DerivativeDomain: - DifferentiableRealMapping<N, F> + ClosedMul<F> + MinMaxMapping<Loc<N, F>, F>, + DifferentiableMapping<U, Codomain = F> + ClosedMul<F> + MinMaxMapping<U, F>, //for<'a> Dat::Differential<'a>: Lipschitz<&'a P, FloatType = F>, - RNDM<N, F>: SpikeMerging<F>, - Reg: SlidingRegTerm<Loc<N, F>, F>, - P: TransportProxPenalty<Loc<N, F>, Dat::DerivativeDomain, Reg, F> + StepLengthBound<F, Dat>, - Plot: Plotter<P::ReturnMapping, Dat::DerivativeDomain, RNDM<N, F>>, + DiscreteMeasure<U, F>: SpikeMerging<F>, + Reg: SlidingRegTerm<U, F>, + P: TransportProxPenalty<U, Dat::DerivativeDomain, Reg, F> + StepLengthBound<F, Dat>, + Plot: Plotter<P::ReturnMapping, Dat::DerivativeDomain, DiscreteMeasure<U, F>>, { pointsource_sliding_fb( f,
--- a/src/sliding_pdps.rs Tue Aug 04 20:11:10 2026 -0500 +++ b/src/sliding_pdps.rs Thu Aug 06 11:14:15 2026 -0500 @@ -7,7 +7,7 @@ use crate::fb::*; use crate::forward_model::{BoundedCurvature, BoundedCurvatureGuess}; use crate::measures::merging::SpikeMerging; -use crate::measures::{DiscreteMeasure, RNDM}; +use crate::measures::DiscreteMeasure; use crate::plot::Plotter; use crate::prox_penalty::StepLengthBoundPair; use crate::regularisation::SlidingRegTerm; @@ -18,13 +18,14 @@ use alg_tools::direct_product::Pair; use alg_tools::error::DynResult; use alg_tools::euclidean::ClosedEuclidean; +use alg_tools::instance::{Instance, Space}; use alg_tools::iterate::AlgIteratorFactory; use alg_tools::linops::{ BoundedLinear, IdOp, SimplyAdjointable, StaticEuclideanOriginGenerator, ZeroOp, AXPY, GEMV, }; -use alg_tools::mapping::{DifferentiableMapping, DifferentiableRealMapping, Instance}; +use alg_tools::mapping::DifferentiableMapping; use alg_tools::nalgebra_support::ToNalgebraRealField; -use alg_tools::norms::L2; +use alg_tools::norms::{Dist, L2}; use anyhow::ensure; use measures::DeltaMeasure; use numeric_literals::replace_float_literals; @@ -64,7 +65,7 @@ } } -type MeasureZ<F, Z, const N: usize> = Pair<RNDM<N, F>, Z>; +type MeasureZ<F, Z, U> = Pair<DiscreteMeasure<U, F>, Z>; /// Iteratively solve the pointsource localisation with an additional variable /// using sliding primal-dual proximal splitting @@ -85,7 +86,7 @@ Plot, /*KOpM, */ KOpZ, H, - const N: usize, + U, >( f: &Dat, reg: &Reg, @@ -93,32 +94,34 @@ config: &SlidingPDPSConfig<F>, iterator: I, mut plotter: Plot, - (μ0, mut z, mut y): (Option<RNDM<N, F>>, Z, Y), + (μ0, mut z, mut y): (Option<DiscreteMeasure<U, F>>, Z, Y), //opKμ : KOpM, opKz: &KOpZ, fnR: &R, fnH: &H, - domain: &mut impl AdmissibleDomain<Loc<N, F>>, -) -> DynResult<MeasureZ<F, Z, N>> + domain: &mut impl AdmissibleDomain<U>, +) -> DynResult<MeasureZ<F, Z, U>> where F: Float + ToNalgebraRealField, + U: Space + Clone + Dist<L2, F>, + for<'a> &'a U: Instance<U>, I: AlgIteratorFactory<IterInfo<F>>, - Dat: DifferentiableMapping<MeasureZ<F, Z, N>, Codomain = F, DerivativeDomain = Pair<S, Z>> + Dat: DifferentiableMapping<MeasureZ<F, Z, U>, Codomain = F, DerivativeDomain = Pair<S, Z>> + BoundedCurvature<F>, - S: DifferentiableRealMapping<N, F> + ClosedMul<F> + MinMaxMapping<Loc<N, F>, F>, + S: DifferentiableMapping<U, Codomain = F> + ClosedMul<F> + MinMaxMapping<U, F>, for<'a> Pair<&'a P, &'a IdOp<Z>>: StepLengthBoundPair<F, Dat>, //Pair<S, Z>: ClosedMul<F>, - RNDM<N, F>: SpikeMerging<F>, - Reg: SlidingRegTerm<Loc<N, F>, F>, - P: TransportProxPenalty<Loc<N, F>, S, Reg, F>, - // KOpM : Linear<RNDM<N, F>, Codomain=Y> - // + GEMV<F, RNDM<N, F>> + DiscreteMeasure<U, F>: SpikeMerging<F>, + Reg: SlidingRegTerm<U, F>, + P: TransportProxPenalty<U, S, Reg, F>, + // KOpM : Linear<DiscreteMeasure<U, F>, Codomain=Y> + // + GEMV<F, DiscreteMeasure<U, F>> // + Preadjointable< - // RNDM<N, F>, Y, + // DiscreteMeasure<U, F>, Y, // PreadjointCodomain = S, // > // + TransportLipschitz<L2Squared, FloatType=F> - // + AdjointProductBoundedBy<RNDM<N, F>, 𝒟, FloatType=F>, + // + AdjointProductBoundedBy<DiscreteMeasure<U, F>, 𝒟, FloatType=F>, // for<'b> KOpM::Preadjoint<'b> : GEMV<F, Y>, // Since Z is Hilbert, we may just as well use adjoints for K_z. KOpZ: BoundedLinear<Z, L2, L2, F, Codomain = Y> @@ -132,7 +135,7 @@ R: Prox<Z, Codomain = F>, H: Conjugable<Y, F, Codomain = F>, for<'b> H::Conjugate<'b>: Prox<Y>, - Plot: Plotter<P::ReturnMapping, S, RNDM<N, F>>, + Plot: Plotter<P::ReturnMapping, S, DiscreteMeasure<U, F>>, { // Check parameters /*ensure!( @@ -206,7 +209,7 @@ let starH = fnH.conjugate(); // Statistics - let full_stats = |μ: &RNDM<N, F>, z: &Z, ε, stats| IterInfo { + let full_stats = |μ: &DiscreteMeasure<U, F>, z: &Z, ε, stats| IterInfo { value: f.apply(Pair(μ, z)) + fnR.apply(z) + reg.apply(μ) @@ -317,7 +320,7 @@ ε, ins, ®, - Some(|μ̃: &RNDM<N, F>| f.apply(Pair(μ̃, &z)) + reg.apply(μ̃)), + Some(|μ̃: &DiscreteMeasure<U, F>| f.apply(Pair(μ̃, &z)) + reg.apply(μ̃)), ); γ.prune_compat(&mut μ, &mut stats); @@ -375,7 +378,7 @@ ε = tolerance.update(ε, iter); } - let fit = |μ̃: &RNDM<N, F>| { + let fit = |μ̃: &DiscreteMeasure<U, F>| { f.apply(Pair(μ̃, &z)) /*+ fnR.apply(z) + reg.apply(μ)*/ + fnH.apply(/* opKμ.apply(&μ̃) + */ opKz.apply(&z)) + reg.apply(μ̃) @@ -392,32 +395,34 @@ /// The implementation uses [`pointsource_sliding_pdps_pair`] with appropriate dummy /// variables, operators, and functions. #[replace_float_literals(F::cast_from(literal))] -pub fn pointsource_sliding_fb_pair_dom<F, I, S, Dat, Reg, P, Z, R, Plot, const N: usize>( +pub fn pointsource_sliding_fb_pair_dom<F, I, S, Dat, Reg, P, Z, R, Plot, U>( f: &Dat, reg: &Reg, prox_penalty: &P, config: &SlidingFBConfig<F>, iterator: I, plotter: Plot, - (μ0, z): (Option<RNDM<N, F>>, Z), + (μ0, z): (Option<DiscreteMeasure<U, F>>, Z), //opKμ : KOpM, fnR: &R, - domain: &mut impl AdmissibleDomain<Loc<N, F>>, -) -> DynResult<MeasureZ<F, Z, N>> + domain: &mut impl AdmissibleDomain<U>, +) -> DynResult<MeasureZ<F, Z, U>> where F: Float + ToNalgebraRealField, + U: Space + Clone + Dist<L2, F>, + for<'a> &'a U: Instance<U>, I: AlgIteratorFactory<IterInfo<F>>, - Dat: DifferentiableMapping<MeasureZ<F, Z, N>, Codomain = F, DerivativeDomain = Pair<S, Z>> + Dat: DifferentiableMapping<MeasureZ<F, Z, U>, Codomain = F, DerivativeDomain = Pair<S, Z>> + BoundedCurvature<F>, - S: DifferentiableRealMapping<N, F> + ClosedMul<F> + MinMaxMapping<Loc<N, F>, F>, - RNDM<N, F>: SpikeMerging<F>, - Reg: SlidingRegTerm<Loc<N, F>, F>, - P: TransportProxPenalty<Loc<N, F>, S, Reg, F>, + S: DifferentiableMapping<U, Codomain = F> + ClosedMul<F> + MinMaxMapping<U, F>, + DiscreteMeasure<U, F>: SpikeMerging<F>, + Reg: SlidingRegTerm<U, F>, + P: TransportProxPenalty<U, S, Reg, F>, for<'a> Pair<&'a P, &'a IdOp<Z>>: StepLengthBoundPair<F, Dat>, Z: ClosedEuclidean<F> + AXPY + Clone, for<'b> &'b Z: Instance<Z>, R: Prox<Z, Codomain = F>, - Plot: Plotter<P::ReturnMapping, S, RNDM<N, F>>, + Plot: Plotter<P::ReturnMapping, S, DiscreteMeasure<U, F>>, // We should not need to explicitly require this: for<'b> &'b Loc<0, F>: Instance<Loc<0, F>>, // Loc<0, F>: StaticEuclidean<Field = F, PrincipalE = Loc<0, F>> @@ -463,7 +468,7 @@ Plot, /*KOpM, */ KOpZ, H, - const N: usize, + U, >( f: &Dat, reg: &Reg, @@ -471,31 +476,33 @@ config: &SlidingPDPSConfig<F>, iterator: I, plotter: Plot, - init: (Option<RNDM<N, F>>, Z, Y), + init: (Option<DiscreteMeasure<U, F>>, Z, Y), //opKμ : KOpM, opKz: &KOpZ, fnR: &R, fnH: &H, -) -> DynResult<MeasureZ<F, Z, N>> +) -> DynResult<MeasureZ<F, Z, U>> where F: Float + ToNalgebraRealField, + U: Space + Clone + Dist<L2, F>, + for<'a> &'a U: Instance<U>, I: AlgIteratorFactory<IterInfo<F>>, - Dat: DifferentiableMapping<MeasureZ<F, Z, N>, Codomain = F, DerivativeDomain = Pair<S, Z>> + Dat: DifferentiableMapping<MeasureZ<F, Z, U>, Codomain = F, DerivativeDomain = Pair<S, Z>> + BoundedCurvature<F>, - S: DifferentiableRealMapping<N, F> + ClosedMul<F> + MinMaxMapping<Loc<N, F>, F>, + S: DifferentiableMapping<U, Codomain = F> + ClosedMul<F> + MinMaxMapping<U, F>, for<'a> Pair<&'a P, &'a IdOp<Z>>: StepLengthBoundPair<F, Dat>, //Pair<S, Z>: ClosedMul<F>, - RNDM<N, F>: SpikeMerging<F>, - Reg: SlidingRegTerm<Loc<N, F>, F>, - P: TransportProxPenalty<Loc<N, F>, S, Reg, F>, - // KOpM : Linear<RNDM<N, F>, Codomain=Y> - // + GEMV<F, RNDM<N, F>> + DiscreteMeasure<U, F>: SpikeMerging<F>, + Reg: SlidingRegTerm<U, F>, + P: TransportProxPenalty<U, S, Reg, F>, + // KOpM : Linear<DiscreteMeasure<U, F>, Codomain=Y> + // + GEMV<F, DiscreteMeasure<U, F>> // + Preadjointable< - // RNDM<N, F>, Y, + // DiscreteMeasure<U, F>, Y, // PreadjointCodomain = S, // > // + TransportLipschitz<L2Squared, FloatType=F> - // + AdjointProductBoundedBy<RNDM<N, F>, 𝒟, FloatType=F>, + // + AdjointProductBoundedBy<DiscreteMeasure<U, F>, 𝒟, FloatType=F>, // for<'b> KOpM::Preadjoint<'b> : GEMV<F, Y>, // Since Z is Hilbert, we may just as well use adjoints for K_z. KOpZ: BoundedLinear<Z, L2, L2, F, Codomain = Y> @@ -509,7 +516,7 @@ R: Prox<Z, Codomain = F>, H: Conjugable<Y, F, Codomain = F>, for<'b> H::Conjugate<'b>: Prox<Y>, - Plot: Plotter<P::ReturnMapping, S, RNDM<N, F>>, + Plot: Plotter<P::ReturnMapping, S, DiscreteMeasure<U, F>>, { pointsource_sliding_pdps_pair_dom( f, @@ -529,31 +536,33 @@ #[deprecated] /// Version of [`pointsource_sliding_fb_pair_dom`] without domain specification #[replace_float_literals(F::cast_from(literal))] -pub fn pointsource_sliding_fb_pair<F, I, S, Dat, Reg, P, Z, R, Plot, const N: usize>( +pub fn pointsource_sliding_fb_pair<F, I, S, Dat, Reg, P, Z, R, Plot, U>( f: &Dat, reg: &Reg, prox_penalty: &P, config: &SlidingFBConfig<F>, iterator: I, plotter: Plot, - (μ0, z): (Option<RNDM<N, F>>, Z), + (μ0, z): (Option<DiscreteMeasure<U, F>>, Z), //opKμ : KOpM, fnR: &R, -) -> DynResult<MeasureZ<F, Z, N>> +) -> DynResult<MeasureZ<F, Z, U>> where F: Float + ToNalgebraRealField, + U: Space + Clone + Dist<L2, F>, + for<'a> &'a U: Instance<U>, I: AlgIteratorFactory<IterInfo<F>>, - Dat: DifferentiableMapping<MeasureZ<F, Z, N>, Codomain = F, DerivativeDomain = Pair<S, Z>> + Dat: DifferentiableMapping<MeasureZ<F, Z, U>, Codomain = F, DerivativeDomain = Pair<S, Z>> + BoundedCurvature<F>, - S: DifferentiableRealMapping<N, F> + ClosedMul<F> + MinMaxMapping<Loc<N, F>, F>, - RNDM<N, F>: SpikeMerging<F>, - Reg: SlidingRegTerm<Loc<N, F>, F>, - P: TransportProxPenalty<Loc<N, F>, S, Reg, F>, + S: DifferentiableMapping<U, Codomain = F> + ClosedMul<F> + MinMaxMapping<U, F>, + DiscreteMeasure<U, F>: SpikeMerging<F>, + Reg: SlidingRegTerm<U, F>, + P: TransportProxPenalty<U, S, Reg, F>, for<'a> Pair<&'a P, &'a IdOp<Z>>: StepLengthBoundPair<F, Dat>, Z: ClosedEuclidean<F> + AXPY + Clone, for<'b> &'b Z: Instance<Z>, R: Prox<Z, Codomain = F>, - Plot: Plotter<P::ReturnMapping, S, RNDM<N, F>>, + Plot: Plotter<P::ReturnMapping, S, DiscreteMeasure<U, F>>, // We should not need to explicitly require this: for<'b> &'b Loc<0, F>: Instance<Loc<0, F>>, // Loc<0, F>: StaticEuclidean<Field = F, PrincipalE = Loc<0, F>>
