Sun, 25 Sep 2022 21:45:56 +0300
Attempt at directly referring to measure masses as SliceStorageMut
src/fb.rs | file | annotate | diff | comparison | revisions | |
src/frank_wolfe.rs | file | annotate | diff | comparison | revisions | |
src/measures/discrete.rs | file | annotate | diff | comparison | revisions | |
src/pdps.rs | file | annotate | diff | comparison | revisions | |
src/subproblem.rs | file | annotate | diff | comparison | revisions |
--- a/src/fb.rs Thu Dec 01 23:07:35 2022 +0200 +++ b/src/fb.rs Sun Sep 25 21:45:56 2022 +0300 @@ -457,7 +457,7 @@ iterator : I, plotter : SeqPlotter<F, N> ) -> DiscreteMeasure<Loc<F, N>, F> -where F : Float + ToNalgebraRealField, +where F : Float + ToNalgebraRealField<MixedType = F>, I : AlgIteratorFactory<IterInfo<F, N>>, for<'b> &'b A::Observable : std::ops::Neg<Output=A::Observable>, //+ std::ops::Mul<F, Output=A::Observable>, <-- FIXME: compiler overflow @@ -521,7 +521,7 @@ mut residual : A::Observable, mut specialisation : Spec, ) -> DiscreteMeasure<Loc<F, N>, F> -where F : Float + ToNalgebraRealField, +where F : Float + ToNalgebraRealField<MixedType=F>, I : AlgIteratorFactory<IterInfo<F, N>>, Spec : FBSpecialisation<F, A::Observable, N>, A::Observable : std::ops::MulAssign<F>, @@ -657,7 +657,7 @@ μ.iter_locations() .map(|ζ| minus_τv.apply(ζ) + ω0.apply(ζ)) .map(F::to_nalgebra_mixed)); - let mut x = μ.masses_dvector(); + let mut x = μ.masses_mut(); // The gradient of the forward component of the inner objective is C^*𝒟Cx - g̃. // We have |C^*𝒟Cx|_2 = sup_{|z|_2 ≤ 1} ⟨z, C^*𝒟Cx⟩ = sup_{|z|_2 ≤ 1} ⟨Cz|𝒟Cx⟩ @@ -672,7 +672,7 @@ inner_τ, inner_it); // Update masses of μ based on solution of finite-dimensional subproblem. - μ.set_masses_dvector(&x); + //μ.set_masses_dvector(&x); } // Form d = ω0 - τv - 𝒟μ = -𝒟(μ - μ^k) - τv for checking the proximate optimality
--- a/src/frank_wolfe.rs Thu Dec 01 23:07:35 2022 +0200 +++ b/src/frank_wolfe.rs Sun Sep 25 21:45:56 2022 +0300 @@ -149,13 +149,13 @@ inner : &InnerSettings<F>, iterator : I ) -> usize -where F : Float + ToNalgebraRealField, +where F : Float + ToNalgebraRealField<MixedType=F>, I : AlgIteratorFactory<F>, A : ForwardModel<Loc<F, N>, F> { // Form and solve finite-dimensional subproblem. let (Ã, g̃) = opA.findim_quadratic_model(&μ, b); - let mut x = μ.masses_dvector(); + let mut x = μ.masses_mut(); // `inner_τ1` is based on an estimate of the operator norm of $A$ from ℳ(Ω) to // ℝ^n. This estimate is a good one for the matrix norm from ℝ^m to ℝ^n when the @@ -168,7 +168,7 @@ let inner_τ = inner.τ0 / (findim_data.opAnorm_squared * F::cast_from(μ.len())); let iters = quadratic_nonneg(inner.method, &Ã, &g̃, α, &mut x, inner_τ, iterator); // Update masses of μ based on solution of finite-dimensional subproblem. - μ.set_masses_dvector(&x); + //μ.set_masses_dvector(&x); iters } @@ -193,7 +193,7 @@ iterator : I, mut plotter : SeqPlotter<F, N>, ) -> DiscreteMeasure<Loc<F, N>, F> -where F : Float + ToNalgebraRealField, +where F : Float + ToNalgebraRealField<MixedType=F>, I : AlgIteratorFactory<IterInfo<F, N>>, for<'b> &'b A::Observable : std::ops::Neg<Output=A::Observable>, //+ std::ops::Mul<F, Output=A::Observable>, <-- FIXME: compiler overflow
--- a/src/measures/discrete.rs Thu Dec 01 23:07:35 2022 +0200 +++ b/src/measures/discrete.rs Sun Sep 25 21:45:56 2022 +0300 @@ -7,7 +7,14 @@ }; use std::iter::Sum; use serde::ser::{Serializer, Serialize, SerializeSeq}; -use nalgebra::DVector; +use nalgebra::{ + DVector, + DVectorSliceMut, + Dynamic, + Const, + SliceStorageMut, + Matrix, +}; use alg_tools::norms::Norm; use alg_tools::tabledump::TableDump; @@ -176,6 +183,25 @@ } } +impl<Domain, F : Float + ToNalgebraRealField<MixedType=F>> DiscreteMeasure<Domain, F> { + pub fn masses_mut(&mut self) -> DVectorSliceMut<'_, F::MixedType, Dynamic> { + let n = self.spikes.len(); + unsafe { + let start = self.spikes.as_mut_ptr(); + let next = start.add(1); + let ptr = &mut (*start).α as *mut F; + let ptrnext = &mut (*next).α as *mut F; + let stride = ptrnext.offset_from(ptr); + assert_eq!(start.add(stride as usize), next); + Matrix::from_data(SliceStorageMut::from_raw_parts( + ptr, + (Dynamic::new(n), Const), + (Dynamic::new(stride as usize), Dynamic::new(0)) + )) + } + } +} + impl<Domain, F :Num> Index<usize> for DiscreteMeasure<Domain, F> { type Output = DeltaMeasure<Domain, F>; #[inline]
--- a/src/pdps.rs Thu Dec 01 23:07:35 2022 +0200 +++ b/src/pdps.rs Sun Sep 25 21:45:56 2022 +0300 @@ -311,7 +311,7 @@ plotter : SeqPlotter<F, N>, dataterm : D, ) -> DiscreteMeasure<Loc<F, N>, F> -where F : Float + ToNalgebraRealField, +where F : Float + ToNalgebraRealField<MixedType=F>, I : AlgIteratorFactory<IterInfo<F, N>>, for<'b> &'b A::Observable : std::ops::Neg<Output=A::Observable> + std::ops::Add<A::Observable, Output=A::Observable>,
--- a/src/subproblem.rs Thu Dec 01 23:07:35 2022 +0200 +++ b/src/subproblem.rs Sun Sep 25 21:45:56 2022 +0300 @@ -1,7 +1,19 @@ //! Iterative algorithms for solving finite-dimensional subproblems. use serde::{Serialize, Deserialize}; -use nalgebra::{DVector, DMatrix}; +use nalgebra::{ + Matrix, + Vector, + Storage, + StorageMut, + Const, + Dim, + DVector, + DMatrix, + DefaultAllocator, + Dynamic, +}; +use nalgebra::base::allocator::Allocator; use numeric_literals::replace_float_literals; use itertools::{izip, Itertools}; use colored::Colorize; @@ -80,25 +92,33 @@ /// The `λ` component of the model is handled in the proximal step instead of the gradient step /// for potential performance improvements. #[replace_float_literals(F::cast_from(literal).to_nalgebra_mixed())] -pub fn quadratic_nonneg_fb<F, I>( - mA : &DMatrix<F::MixedType>, - g : &DVector<F::MixedType>, +pub fn quadratic_nonneg_fb<F, I, S, SA, SG, D>( + mA : &Matrix<F::MixedType, D, D, SA>, + g : &Vector<F::MixedType, D, SG>, //c_ : F, λ_ : F, - x : &mut DVector<F::MixedType>, + x : &mut Vector<F::MixedType, D, S>, τ_ : F, iterator : I ) -> usize where F : Float + ToNalgebraRealField, - I : AlgIteratorFactory<F> + I : AlgIteratorFactory<F>, + D : Dim, + S : StorageMut<F::MixedType, D, Const::<1>>, + SA : Storage<F::MixedType, D, D>, + SG : Storage<F::MixedType, D, Const::<1>>, + DefaultAllocator : Allocator<F::MixedType, D> { - let mut xprev = x.clone(); + let mut xprev = x.clone_owned(); //let c = c_.to_nalgebra_mixed(); let λ = λ_.to_nalgebra_mixed(); let τ = τ_.to_nalgebra_mixed(); let τλ = τ * λ; - let mut v = DVector::zeros(x.len()); let mut iters = 0; + let mut v = { + let (r, c) = x.shape_generic(); + Vector::zeros_generic(r, c) + }; iterator.iterate(|state| { // Replace `x` with $x - τ[Ax-g]= [x + τg]- τAx$ @@ -195,26 +215,33 @@ /// forward-backward step. /// </p> #[replace_float_literals(F::cast_from(literal).to_nalgebra_mixed())] -pub fn quadratic_nonneg_ssn<F, I>( - mA : &DMatrix<F::MixedType>, - g : &DVector<F::MixedType>, +pub fn quadratic_nonneg_ssn<F, I, S, SA, SG>( + mA : &Matrix<F::MixedType, Dynamic, Dynamic, SA>, + g : &Vector<F::MixedType, Dynamic, SG>, //c_ : F, λ_ : F, - x : &mut DVector<F::MixedType>, + x : &mut Vector<F::MixedType, Dynamic, S>, τ_ : F, iterator : I ) -> Result<usize, NumericalError> -where F : Float + ToNalgebraRealField, - I : AlgIteratorFactory<F> +where F : Float + ToNalgebraRealField<MixedType=F>, + I : AlgIteratorFactory<F>, + S : StorageMut<F::MixedType, Dynamic, Const::<1>>, + SA : Storage<F::MixedType, Dynamic, Dynamic>, + SG : Storage<F::MixedType, Dynamic, Const::<1>>, { let n = x.len(); - let mut xprev = x.clone(); + let mut xprev = x.clone_owned(); let mut v = DVector::zeros(n); //let c = c_.to_nalgebra_mixed(); let λ = λ_.to_nalgebra_mixed(); let τ = τ_.to_nalgebra_mixed(); let τλ = τ * λ; let mut inact : Vec<bool> = Vec::from_iter(std::iter::repeat(false).take(n)); + let mut v = { + let (r, c) = x.shape_generic(); + Vector::zeros_generic(r, c) + }; let mut s = DVector::zeros(0); let mut decomp = nalgebra::linalg::LU::new(DMatrix::zeros(0, 0)); let mut iters = 0; @@ -346,18 +373,21 @@ /// cone_, <https://doi.org/10.1080/02331934.2014.929680>. /// /// This function returns the number of iterations taken. -pub fn quadratic_nonneg<F, I>( +pub fn quadratic_nonneg<F, I, S, SA, SG>( method : InnerMethod, - mA : &DMatrix<F::MixedType>, - g : &DVector<F::MixedType>, + mA : &Matrix<F::MixedType, Dynamic, Dynamic, SA>, + g : &Vector<F::MixedType, Dynamic, SG>, //c_ : F, λ : F, - x : &mut DVector<F::MixedType>, + x : &mut Vector<F::MixedType, Dynamic, S>, τ : F, iterator : I ) -> usize -where F : Float + ToNalgebraRealField, - I : AlgIteratorFactory<F> +where F : Float + ToNalgebraRealField<MixedType=F>, + I : AlgIteratorFactory<F>, + S : StorageMut<F::MixedType, Dynamic, Const::<1>>, + SA : Storage<F::MixedType, Dynamic, Dynamic>, + SG : Storage<F::MixedType, Dynamic, Const::<1>>, { match method {