pointsource_algs: comparison src/fb.rs

-:c32171f7cce5
+:ea3ca78873e8
 let quiet = iterator.is_quiet();
 let op𝒟norm = op𝒟.opnorm_bound();
 // We multiply tolerance by τ for FB since
 // our subproblems depending on tolerances are scaled by τ compared to the conditional
 // gradient approach.
-let mut tolerance = config.tolerance * τ * α;
+let tolerance = config.tolerance * τ * α;
 let mut ε = tolerance.initial();
 // Initialise operators
 let preadjA = opA.preadjoint();
 // Run the algorithm
 iterator.iterate(|state| {
 // Calculate subproblem tolerances, and update main tolerance for next iteration
 let τα = τ * α;
-// if μ.len() == 0 /*state.iteration() == 1*/ {
-//     let t = minus_τv.bounds().upper() * 0.001;
-//     if t > 0.0 {
-//         let (ξ, v_ξ) = minus_τv.maximise(t, config.refinement.max_steps);
-//         if τα + ε > v_ξ && v_ξ > τα {
-//             // The zero measure is already within bounds, so improve them
-//             tolerance = config.tolerance * (v_ξ - τα);
-//             ε = tolerance.initial();
-//         }
-//         μ += DeltaMeasure { x : ξ, α : 0.0 };
-//     } else {
-//         // Zero is the solution.
-//         return Step::Terminated
-//     }
-// }
 let target_bounds = Bounds(τα - ε,  τα + ε);
 let merge_tolerance = config.merge_tolerance_mult * ε;
 let merge_target_bounds = Bounds(τα - merge_tolerance,  τα + merge_tolerance);
 let inner_tolerance = ε * config.inner.tolerance_mult;
 let refinement_tolerance = ε * config.refinement.tolerance_mult;
 let maximise_above = τα + ε * config.insertion_cutoff_factor;
-let mut ε1 = ε;
 let ε_prev = ε;
 ε = tolerance.update(ε, state.iteration());
 // Maximum insertion count and measure difference calculation depend on insertion style.
 let (m, warn_insertions) = match (state.iteration(), config.bootstrap_insertions) {
 // the calculations in the beginning of the loop that v_ξ = (ω0-τv-𝒟μ)(ξ) = d(ξ),
 // where 𝒟μ is now distinct from μ0 after the insertions already performed.
 // We do not need to check lower bounds, as a solution of the finite-dimensional
 // subproblem should always satisfy them.
-// // Find the mimimum over the support of μ.
+// If μ has some spikes, only find a maximum of d if it is above a threshold
-// let d_min_supp = d_max;μ.iter_spikes().filter_map(|&DeltaMeasure{ α, ref x }| {
+// defined by the refinment tolerance.
-//    (α != F::ZERO).then(|| d.value(x))
+let (ξ, v_ξ) =  match d.maximise_above(maximise_above, refinement_tolerance,
-// }).reduce(F::min).unwrap_or(0.0);
+config.refinement.max_steps) {
+None => break 'insertion (true, d),
-let (ξ, v_ξ) = if false /* μ.len() == 0*/ /*count == 0 &&*/ {
+Some(res) => res,
-// If μ has no spikes, just find the maximum of d. Then adjust the tolerance, if
-// necessary, to adapt it to the problem.
-let (ξ, v_ξ) = d.maximise(refinement_tolerance, config.refinement.max_steps);
-//dbg!((τα, v_ξ, target_bounds.upper(), maximise_above));
-if τα < v_ξ  && v_ξ < target_bounds.upper() {
-ε1 = v_ξ - τα;
-ε *= ε1 / ε_prev;
-tolerance *= ε1 / ε_prev;
-}
-(ξ, v_ξ)
-} else {
-// If μ has some spikes, only find a maximum of d if it is above a threshold
-// defined by the refinment tolerance.
-match d.maximise_above(maximise_above, refinement_tolerance,
-config.refinement.max_steps) {
-None => break 'insertion (true, d),
-Some(res) => res,
-}
 };
-// // Do a one final check whether we can stop already without inserting more points
-// // because `d` actually in bounds based on a more refined estimate.
-// if may_break && target_bounds.upper() >= v_ξ {
-//     break (true, d)
-// }
 // Break if maximum insertion count reached
 if count >= max_insertions {
 let in_bounds2 = target_bounds.upper() >= v_ξ;
 break 'insertion (in_bounds2, d)
 // Merge spikes
 if state.iteration() % config.merge_every == 0 {
 let n_before_merge = μ.len();
 μ.merge_spikes(config.merging, |μ_candidate| {
-//println!("Merge attempt!");
 let mut d = &minus_τv + op𝒟.preapply(μ_diff(&μ_candidate, &μ_base));
 if merge_target_bounds.superset(&d.bounds()) {
-//println!("…Early Ok");
 return Some(())
 }
 let d_min_supp = μ_candidate.iter_spikes().filter_map(|&DeltaMeasure{ α, ref x }| {
 (α != 0.0).then(|| d.apply(x))
 }).reduce(F::min);
 if d_min_supp.map_or(true, |b| b >= merge_target_bounds.lower()) &&
 d.has_upper_bound(merge_target_bounds.upper(), refinement_tolerance,
 config.refinement.max_steps) {
-//println!("…Ok");
 Some(())
 } else {
-//println!("…Fail");
 None
 }
 });
 debug_assert!(μ.len() >= n_before_merge);
 merged += μ.len() - n_before_merge;
 plotter.plot_spikes(
 format!("iter {} end; {}", state.iteration(), within_tolerances), &d,
 "start".to_string(), Some(&minus_τv),
 Some(target_bounds), value_μ,
 );
-// Calculate mean inner iterations and reset relevant counters
+// Calculate mean inner iterations and reset relevant counters.
 // Return the statistics
 let res = IterInfo {
 value : specialisation.calculate_fit(&μ, &residual) + α * value_μ.norm(Radon),
 n_spikes : value_μ.len(),
 inner_iters,
 this_iters,
 merged,
 pruned,
 ε : ε_prev,
-maybe_ε1 : Some(ε1),
 postprocessing: config.postprocessing.then(|| value_μ.clone()),
 };
 inner_iters = 0;
 this_iters = 0;
 merged = 0;

Mercurial > repos > pointsource_algs / file comparison

comparison: src/fb.rs

src/fb.rs