src/prox_penalty/radon_squared.rs

branch
dev
changeset 63
7a8a55fd41c0
parent 61
4f468d35fa29
--- 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()

mercurial