--- a/src/prox_penalty/radon_squared.rs Thu Feb 26 11:38:43 2026 -0500 +++ b/src/prox_penalty/radon_squared.rs Thu Feb 26 11:36:22 2026 -0500 @@ -8,8 +8,8 @@ use crate::dataterm::QuadraticDataTerm; use crate::forward_model::ForwardModel; use crate::measures::merging::SpikeMerging; -use crate::measures::{DeltaMeasure, DiscreteMeasure, Radon}; -use crate::regularisation::RegTerm; +use crate::measures::{DiscreteMeasure, Radon}; +use crate::regularisation::RadonSquaredRegTerm; use crate::types::*; use alg_tools::bounds::MinMaxMapping; use alg_tools::error::DynResult; @@ -18,8 +18,6 @@ use alg_tools::linops::BoundedLinear; use alg_tools::nalgebra_support::ToNalgebraRealField; use alg_tools::norms::{Norm, L2}; -use anyhow::ensure; -use nalgebra::DVector; use numeric_literals::replace_float_literals; use serde::{Deserialize, Serialize}; @@ -35,7 +33,7 @@ for<'a> &'a Domain: Instance<Domain>, F: Float + ToNalgebraRealField, M: MinMaxMapping<Domain, F>, - Reg: RegTerm<Domain, F>, + Reg: RadonSquaredRegTerm<Domain, F>, DiscreteMeasure<Domain, F>: SpikeMerging<F>, { type ReturnMapping = M; @@ -48,8 +46,6 @@ &self, μ: &mut DiscreteMeasure<Domain, F>, τv: &mut M, - μ_base: &DiscreteMeasure<Domain, F>, - ν_delta: Option<&DiscreteMeasure<Domain, F>>, τ: F, ε: F, config: &InsertionConfig<F>, @@ -60,58 +56,8 @@ where I: AlgIterator, { - let mut y = μ_base.masses_dvector(); - - ensure!(μ_base.len() <= μ.len()); - - 'i_and_w: for i in 0..=1 { - // Optimise weights - if μ.len() > 0 { - // Form finite-dimensional subproblem. The subproblem references to the original μ^k - // from the beginning of the iteration are all contained in the immutable c and g. - // TODO: observe negation of -τv after switch from minus_τv: finite-dimensional - // problems have not yet been updated to sign change. - let g̃ = DVector::from_iterator( - μ.len(), - μ.iter_locations() - .map(|ζ| -F::to_nalgebra_mixed(τv.apply(ζ))), - ); - let mut x = μ.masses_dvector(); - y.extend(std::iter::repeat(0.0.to_nalgebra_mixed()).take(0.max(x.len() - y.len()))); - assert_eq!(y.len(), x.len()); - // Solve finite-dimensional subproblem. - // TODO: This assumes that ν_delta has no common locations with μ-μ_base, to - // ignore it. - stats.inner_iters += reg.solve_findim_l1squared(&y, &g̃, τ, &mut x, ε, config); - - // Update masses of μ based on solution of finite-dimensional subproblem. - μ.set_masses_dvector(&x); - } - - if i > 0 { - // Simple debugging test to see if more inserts would be needed. Doesn't seem so. - //let n = μ.dist_matching(μ_base); - //println!("{:?}", reg.find_tolerance_violation_slack(τv, τ, ε, false, config, n)); - break 'i_and_w; - } - - // Calculate ‖μ - μ_base‖_ℳ - // TODO: This assumes that ν_delta has no common locations with μ-μ_base. - let n = μ.dist_matching(μ_base) + ν_delta.map_or(0.0, |ν| ν.norm(Radon)); - - // Find a spike to insert, if needed. - // This only check the overall tolerances, not tolerances on support of μ-μ_base or μ, - // which are supposed to have been guaranteed by the finite-dimensional weight optimisation. - match reg.find_tolerance_violation_slack(τv, τ, ε, false, config, n) { - None => break 'i_and_w, - Some((ξ, _v_ξ, _in_bounds)) => { - // Weight is found out by running the finite-dimensional optimisation algorithm - // above - *μ += DeltaMeasure { x: ξ, α: 0.0 }; - stats.inserted += 1; - } - }; - } + let violation = reg.find_tolerance_violation(τv, τ, ε, true, config); + reg.solve_oc_radonsq(μ, τv, τ, ε, violation, config, stats); Ok((None, true)) } @@ -121,7 +67,6 @@ μ: &mut DiscreteMeasure<Domain, F>, τv: &mut M, μ_base: &DiscreteMeasure<Domain, F>, - ν_delta: Option<&DiscreteMeasure<Domain, F>>, τ: F, ε: F, config: &InsertionConfig<F>, @@ -141,10 +86,7 @@ // after μ_candidate's extra points. // TODO: doesn't seem to work, maybe need to merge μ_base as well? // Although that doesn't seem to make sense. - let μ_radon = match ν_delta { - None => μ_candidate.sub_matching(μ_base), - Some(ν) => μ_candidate.sub_matching(μ_base) - ν, - }; + let μ_radon = μ_candidate.sub_matching(μ_base); reg.verify_merge_candidate_radonsq(τv, μ_candidate, τ, ε, &config, &μ_radon) //let n = μ_candidate.dist_matching(μ_base); //reg.find_tolerance_violation_slack(τv, τ, ε, false, config, n).is_none()