Mercurial > repos > pointsource_algs / file comparison

comparison: src/regularisation.rs

src/regularisation.rs

branch
dev
changeset 35
b087e3eab191
parent 34
efa60bc4f743
child 51
0693cc9ba9f0
equal deleted inserted replaced
34:efa60bc4f743 35:b087e3eab191
3 */ 3 */
4 4
5 use numeric_literals::replace_float_literals; 5 use numeric_literals::replace_float_literals;
6 use serde::{Serialize, Deserialize}; 6 use serde::{Serialize, Deserialize};
7 use alg_tools::norms::Norm; 7 use alg_tools::norms::Norm;
8 use alg_tools::linops::Apply; 8 use alg_tools::linops::Mapping;
9 use alg_tools::instance::Instance;
9 use alg_tools::loc::Loc; 10 use alg_tools::loc::Loc;
10 use crate::types::*; 11 use crate::types::*;
11 use crate::measures::{ 12 use crate::measures::{
12 DiscreteMeasure, 13 RNDM,
13 DeltaMeasure, 14 DeltaMeasure,
14 Radon 15 Radon,
15 }; 16 };
16 use crate::fb::FBGenericConfig; 17 use crate::fb::FBGenericConfig;
17 #[allow(unused_imports)] // Used by documentation. 18 #[allow(unused_imports)] // Used by documentation.
18 use crate::fb::pointsource_fb_reg; 19 use crate::fb::pointsource_fb_reg;
19 #[allow(unused_imports)] // Used by documentation. 20 #[allow(unused_imports)] // Used by documentation.
20 use crate::sliding_fb::pointsource_sliding_fb_reg; 21 use crate::sliding_fb::pointsource_sliding_fb_reg;
21 22
22 use nalgebra::{DVector, DMatrix}; 23 use nalgebra::{DVector, DMatrix};
23 use alg_tools::nalgebra_support::ToNalgebraRealField; 24 use alg_tools::nalgebra_support::ToNalgebraRealField;
24 use alg_tools::mapping::Mapping;
25 use alg_tools::bisection_tree::{ 25 use alg_tools::bisection_tree::{
26 BTFN, 26 BTFN,
27 Bounds, 27 Bounds,
28 BTSearch, 28 BTSearch,
29 P2Minimise, 29 P2Minimise,
54 let &NonnegRadonRegTerm(α) = self; 54 let &NonnegRadonRegTerm(α) = self;
55 α 55 α
56 } 56 }
57 } 57 }
58 58
59 impl<'a, F : Float, const N : usize> Apply<&'a DiscreteMeasure<Loc<F, N>, F>> 59 impl<'a, F : Float, const N : usize> Mapping<RNDM<F, N>>
60 for NonnegRadonRegTerm<F> { 60 for NonnegRadonRegTerm<F> {
61 type Output = F; 61 type Codomain = F;
62 62
63 fn apply(&self, μ : &'a DiscreteMeasure<Loc<F, N>, F>) -> F { 63 fn apply<I>(&self, μ : I) -> F
64 self.α() * μ.norm(Radon) 64 where I : Instance<RNDM<F, N>> {
65 self.α() * μ.eval(|x| x.norm(Radon))
65 } 66 }
66 } 67 }
67 68
68 69
69 /// The regularisation term $α\|μ\|_{ℳ(Ω)}$ for [`pointsource_fb_reg`]. 70 /// The regularisation term $α\|μ\|_{ℳ(Ω)}$ for [`pointsource_fb_reg`].
79 let &RadonRegTerm(α) = self; 80 let &RadonRegTerm(α) = self;
80 α 81 α
81 } 82 }
82 } 83 }
83 84
84 impl<'a, F : Float, const N : usize> Apply<&'a DiscreteMeasure<Loc<F, N>, F>> 85 impl<'a, F : Float, const N : usize> Mapping<RNDM<F, N>>
85 for RadonRegTerm<F> { 86 for RadonRegTerm<F> {
86 type Output = F; 87 type Codomain = F;
87 88
88 fn apply(&self, μ : &'a DiscreteMeasure<Loc<F, N>, F>) -> F { 89 fn apply<I>(&self, μ : I) -> F
89 self.α() * μ.norm(Radon) 90 where I : Instance<RNDM<F, N>> {
91 self.α() * μ.eval(|x| x.norm(Radon))
90 } 92 }
91 } 93 }
92 94
93 /// Regularisation term configuration 95 /// Regularisation term configuration
94 #[derive(Clone, Copy, Eq, PartialEq, Serialize, Deserialize, Debug)] 96 #[derive(Clone, Copy, Eq, PartialEq, Serialize, Deserialize, Debug)]
97 Radon(F), 99 Radon(F),
98 /// $α\\|μ\\|\_{ℳ(Ω)} + δ_{≥ 0}(μ)$ 100 /// $α\\|μ\\|\_{ℳ(Ω)} + δ_{≥ 0}(μ)$
99 NonnegRadon(F), 101 NonnegRadon(F),
100 } 102 }
101 103
102 impl<'a, F : Float, const N : usize> Apply<&'a DiscreteMeasure<Loc<F, N>, F>> 104 impl<'a, F : Float, const N : usize> Mapping<RNDM<F, N>>
103 for Regularisation<F> { 105 for Regularisation<F> {
104 type Output = F; 106 type Codomain = F;
105 107
106 fn apply(&self, μ : &'a DiscreteMeasure<Loc<F, N>, F>) -> F { 108 fn apply<I>(&self, μ : I) -> F
109 where I : Instance<RNDM<F, N>> {
107 match *self { 110 match *self {
108 Self::Radon(α) => RadonRegTerm(α).apply(μ), 111 Self::Radon(α) => RadonRegTerm(α).apply(μ),
109 Self::NonnegRadon(α) => NonnegRadonRegTerm(α).apply(μ), 112 Self::NonnegRadon(α) => NonnegRadonRegTerm(α).apply(μ),
110 } 113 }
111 } 114 }
112 } 115 }
113 116
114 /// Abstraction of regularisation terms for [`generic_pointsource_fb_reg`]. 117 /// Abstraction of regularisation terms.
115 pub trait RegTerm<F : Float + ToNalgebraRealField, const N : usize> 118 pub trait RegTerm<F : Float + ToNalgebraRealField, const N : usize>
116 : for<'a> Apply<&'a DiscreteMeasure<Loc<F, N>, F>, Output = F> { 119 : Mapping<RNDM<F, N>, Codomain = F> {
117 /// Approximately solve the problem 120 /// Approximately solve the problem
118 /// <div>$$ 121 /// <div>$$
119 /// \min_{x ∈ ℝ^n} \frac{1}{2} x^⊤Ax - g^⊤ x + τ G(x) 122 /// \min_{x ∈ ℝ^n} \frac{1}{2} x^⊤Ax - g^⊤ x + τ G(x)
120 /// $$</div> 123 /// $$</div>
121 /// for $G$ depending on the trait implementation. 124 /// for $G$ depending on the trait implementation.
169 skip_by_rough_check : bool, 172 skip_by_rough_check : bool,
170 config : &FBGenericConfig<F>, 173 config : &FBGenericConfig<F>,
171 ) -> Option<(Loc<F, N>, F, bool)> 174 ) -> Option<(Loc<F, N>, F, bool)>
172 where BT : BTSearch<F, N, Agg=Bounds<F>>, 175 where BT : BTSearch<F, N, Agg=Bounds<F>>,
173 G : SupportGenerator<F, N, Id=BT::Data>, 176 G : SupportGenerator<F, N, Id=BT::Data>,
174 G::SupportType : Mapping<Loc<F, N>,Codomain=F> 177 G::SupportType : Mapping<Loc<F, N>, Codomain=F> + LocalAnalysis<F, Bounds<F>, N> {
175 + LocalAnalysis<F, Bounds<F>, N> {
176 self.find_tolerance_violation_slack(d, τ, ε, skip_by_rough_check, config, F::ZERO) 178 self.find_tolerance_violation_slack(d, τ, ε, skip_by_rough_check, config, F::ZERO)
177 } 179 }
178 180
179 /// Find a point where `d` may violate the tolerance `ε`. 181 /// Find a point where `d` may violate the tolerance `ε`.
180 /// 182 ///
197 config : &FBGenericConfig<F>, 199 config : &FBGenericConfig<F>,
198 slack : F, 200 slack : F,
199 ) -> Option<(Loc<F, N>, F, bool)> 201 ) -> Option<(Loc<F, N>, F, bool)>
200 where BT : BTSearch<F, N, Agg=Bounds<F>>, 202 where BT : BTSearch<F, N, Agg=Bounds<F>>,
201 G : SupportGenerator<F, N, Id=BT::Data>, 203 G : SupportGenerator<F, N, Id=BT::Data>,
202 G::SupportType : Mapping<Loc<F, N>,Codomain=F> 204 G::SupportType : Mapping<Loc<F, N>, Codomain=F> + LocalAnalysis<F, Bounds<F>, N>;
203 + LocalAnalysis<F, Bounds<F>, N>;
204 205
205 206
206 /// Verify that `d` is in bounds `ε` for a merge candidate `μ` 207 /// Verify that `d` is in bounds `ε` for a merge candidate `μ`
207 /// 208 ///
208 /// `ε` is the current main tolerance and `τ` a scaling factor for the regulariser. 209 /// `ε` is the current main tolerance and `τ` a scaling factor for the regulariser.
209 fn verify_merge_candidate<G, BT>( 210 fn verify_merge_candidate<G, BT>(
210 &self, 211 &self,
211 d : &mut BTFN<F, G, BT, N>, 212 d : &mut BTFN<F, G, BT, N>,
212 μ : &DiscreteMeasure<Loc<F, N>, F>, 213 μ : &RNDM<F, N>,
213 τ : F, 214 τ : F,
214 ε : F, 215 ε : F,
215 config : &FBGenericConfig<F>, 216 config : &FBGenericConfig<F>,
216 ) -> bool 217 ) -> bool
217 where BT : BTSearch<F, N, Agg=Bounds<F>>, 218 where BT : BTSearch<F, N, Agg=Bounds<F>>,
218 G : SupportGenerator<F, N, Id=BT::Data>, 219 G : SupportGenerator<F, N, Id=BT::Data>,
219 G::SupportType : Mapping<Loc<F, N>,Codomain=F> 220 G::SupportType : Mapping<Loc<F, N>, Codomain=F> + LocalAnalysis<F, Bounds<F>, N>;
220 + LocalAnalysis<F, Bounds<F>, N>;
221 221
222 /// Verify that `d` is in bounds `ε` for a merge candidate `μ` 222 /// Verify that `d` is in bounds `ε` for a merge candidate `μ`
223 /// 223 ///
224 /// This version is s used for Radon-norm squared proximal term in [`crate::radon_fb`]. 224 /// This version is s used for Radon-norm squared proximal term in [`crate::radon_fb`].
225 /// The [`DiscreteMeasure`]s `μ` and `radon_μ` are supposed to have same coordinates at 225 /// The [measures][crate::measures::DiscreteMeasure] `μ` and `radon_μ` are supposed to have
226 /// same agreeing indices. 226 /// same coordinates at same agreeing indices.
227 /// 227 ///
228 /// `ε` is the current main tolerance and `τ` a scaling factor for the regulariser. 228 /// `ε` is the current main tolerance and `τ` a scaling factor for the regulariser.
229 fn verify_merge_candidate_radonsq<G, BT>( 229 fn verify_merge_candidate_radonsq<G, BT>(
230 &self, 230 &self,
231 d : &mut BTFN<F, G, BT, N>, 231 d : &mut BTFN<F, G, BT, N>,
232 μ : &DiscreteMeasure<Loc<F, N>, F>, 232 μ : &RNDM<F, N>,
233 τ : F, 233 τ : F,
234 ε : F, 234 ε : F,
235 config : &FBGenericConfig<F>, 235 config : &FBGenericConfig<F>,
236 radon_μ :&DiscreteMeasure<Loc<F, N>, F>, 236 radon_μ :&RNDM<F, N>,
237 ) -> bool 237 ) -> bool
238 where BT : BTSearch<F, N, Agg=Bounds<F>>, 238 where BT : BTSearch<F, N, Agg=Bounds<F>>,
239 G : SupportGenerator<F, N, Id=BT::Data>, 239 G : SupportGenerator<F, N, Id=BT::Data>,
240 G::SupportType : Mapping<Loc<F, N>,Codomain=F> 240 G::SupportType : Mapping<Loc<F, N>, Codomain=F> + LocalAnalysis<F, Bounds<F>, N>;
241 + LocalAnalysis<F, Bounds<F>, N>;
242 241
243 242
244 /// TODO: document this 243 /// TODO: document this
245 fn target_bounds(&self, τ : F, ε : F) -> Option<Bounds<F>>; 244 fn target_bounds(&self, τ : F, ε : F) -> Option<Bounds<F>>;
246 245
256 /// Calculate $τ[w(z) - w(y)]$ for some w in the subdifferential of the regularisation 255 /// Calculate $τ[w(z) - w(y)]$ for some w in the subdifferential of the regularisation
257 /// term, such that $-ε ≤ τw - d ≤ ε$. 256 /// term, such that $-ε ≤ τw - d ≤ ε$.
258 fn goodness<G, BT>( 257 fn goodness<G, BT>(
259 &self, 258 &self,
260 d : &mut BTFN<F, G, BT, N>, 259 d : &mut BTFN<F, G, BT, N>,
261 μ : &DiscreteMeasure<Loc<F, N>, F>, 260 μ : &RNDM<F, N>,
262 y : &Loc<F, N>, 261 y : &Loc<F, N>,
263 z : &Loc<F, N>, 262 z : &Loc<F, N>,
264 τ : F, 263 τ : F,
265 ε : F, 264 ε : F,
266 config : &FBGenericConfig<F>, 265 config : &FBGenericConfig<F>,
267 ) -> F 266 ) -> F
268 where BT : BTSearch<F, N, Agg=Bounds<F>>, 267 where BT : BTSearch<F, N, Agg=Bounds<F>>,
269 G : SupportGenerator<F, N, Id=BT::Data>, 268 G : SupportGenerator<F, N, Id=BT::Data>,
270 G::SupportType : Mapping<Loc<F, N>,Codomain=F> 269 G::SupportType : Mapping<Loc<F, N>, Codomain=F> + LocalAnalysis<F, Bounds<F>, N>;
271 + LocalAnalysis<F, Bounds<F>, N>;
272 270
273 /// Convert bound on the regulariser to a bond on the Radon norm 271 /// Convert bound on the regulariser to a bond on the Radon norm
274 fn radon_norm_bound(&self, b : F) -> F; 272 fn radon_norm_bound(&self, b : F) -> F;
275 } 273 }
276 274
321 config : &FBGenericConfig<F>, 319 config : &FBGenericConfig<F>,
322 slack : F 320 slack : F
323 ) -> Option<(Loc<F, N>, F, bool)> 321 ) -> Option<(Loc<F, N>, F, bool)>
324 where BT : BTSearch<F, N, Agg=Bounds<F>>, 322 where BT : BTSearch<F, N, Agg=Bounds<F>>,
325 G : SupportGenerator<F, N, Id=BT::Data>, 323 G : SupportGenerator<F, N, Id=BT::Data>,
326 G::SupportType : Mapping<Loc<F, N>,Codomain=F> 324 G::SupportType : Mapping<Loc<F, N>, Codomain=F> + LocalAnalysis<F, Bounds<F>, N> {
327 + LocalAnalysis<F, Bounds<F>, N> { 325 let τα = τ * self.α();
328 let τα = τ * self.α(); 326 let keep_above = -τα - slack - ε;
329 let keep_below = τα + slack + ε; 327 let minimise_below = -τα - slack - ε * config.insertion_cutoff_factor;
330 let maximise_above = τα + slack + ε * config.insertion_cutoff_factor;
331 let refinement_tolerance = ε * config.refinement.tolerance_mult; 328 let refinement_tolerance = ε * config.refinement.tolerance_mult;
332 329
333 // If preliminary check indicates that we are in bounds, and if it otherwise matches 330 // If preliminary check indicates that we are in bounds, and if it otherwise matches
334 // the insertion strategy, skip insertion. 331 // the insertion strategy, skip insertion.
335 if skip_by_rough_check && d.bounds().upper() <= keep_below { 332 if skip_by_rough_check && d.bounds().lower() >= keep_above {
336 None 333 None
337 } else { 334 } else {
338 // If the rough check didn't indicate no insertion needed, find maximising point. 335 // If the rough check didn't indicate no insertion needed, find minimising point.
339 d.maximise_above(maximise_above, refinement_tolerance, config.refinement.max_steps) 336 d.minimise_below(minimise_below, refinement_tolerance, config.refinement.max_steps)
340 .map(|(ξ, v_ξ)| (ξ, v_ξ, v_ξ <= keep_below)) 337 .map(|(ξ, v_ξ)| (ξ, v_ξ, v_ξ >= keep_above))
341 } 338 }
342 } 339 }
343 340
344 fn verify_merge_candidate<G, BT>( 341 fn verify_merge_candidate<G, BT>(
345 &self, 342 &self,
346 d : &mut BTFN<F, G, BT, N>, 343 d : &mut BTFN<F, G, BT, N>,
347 μ : &DiscreteMeasure<Loc<F, N>, F>, 344 μ : &RNDM<F, N>,
348 τ : F, 345 τ : F,
349 ε : F, 346 ε : F,
350 config : &FBGenericConfig<F>, 347 config : &FBGenericConfig<F>,
351 ) -> bool 348 ) -> bool
352 where BT : BTSearch<F, N, Agg=Bounds<F>>, 349 where BT : BTSearch<F, N, Agg=Bounds<F>>,
353 G : SupportGenerator<F, N, Id=BT::Data>, 350 G : SupportGenerator<F, N, Id=BT::Data>,
354 G::SupportType : Mapping<Loc<F, N>,Codomain=F> 351 G::SupportType : Mapping<Loc<F, N>, Codomain=F> + LocalAnalysis<F, Bounds<F>, N> {
355 + LocalAnalysis<F, Bounds<F>, N> {
356 let τα = τ * self.α(); 352 let τα = τ * self.α();
357 let refinement_tolerance = ε * config.refinement.tolerance_mult; 353 let refinement_tolerance = ε * config.refinement.tolerance_mult;
358 let merge_tolerance = config.merge_tolerance_mult * ε; 354 let merge_tolerance = config.merge_tolerance_mult * ε;
359 let keep_below = τα + merge_tolerance; 355 let keep_above = -τα - merge_tolerance;
360 let keep_supp_above = τα - merge_tolerance; 356 let keep_supp_below = -τα + merge_tolerance;
361 let bnd = d.bounds(); 357 let bnd = d.bounds();
362 358
363 return ( 359 return (
364 bnd.lower() >= keep_supp_above 360 bnd.upper() <= keep_supp_below
365 || 361 ||
366 μ.iter_spikes().all(|&DeltaMeasure{ α, ref x }| { 362 μ.iter_spikes().all(|&DeltaMeasure{ α, ref x }| {
367 (α == 0.0) || d.apply(x) >= keep_supp_above 363 (α == 0.0) || d.apply(x) <= keep_supp_below
368 }) 364 })
369 ) && ( 365 ) && (
370 bnd.upper() <= keep_below 366 bnd.lower() >= keep_above
371 || 367 ||
372 d.has_upper_bound(keep_below, refinement_tolerance, config.refinement.max_steps) 368 d.has_lower_bound(keep_above, refinement_tolerance, config.refinement.max_steps)
373 ) 369 )
374 } 370 }
375 371
376 fn verify_merge_candidate_radonsq<G, BT>( 372 fn verify_merge_candidate_radonsq<G, BT>(
377 &self, 373 &self,
378 d : &mut BTFN<F, G, BT, N>, 374 d : &mut BTFN<F, G, BT, N>,
379 μ : &DiscreteMeasure<Loc<F, N>, F>, 375 μ : &RNDM<F, N>,
380 τ : F, 376 τ : F,
381 ε : F, 377 ε : F,
382 config : &FBGenericConfig<F>, 378 config : &FBGenericConfig<F>,
383 radon_μ :&DiscreteMeasure<Loc<F, N>, F>, 379 radon_μ :&RNDM<F, N>,
384 ) -> bool 380 ) -> bool
385 where BT : BTSearch<F, N, Agg=Bounds<F>>, 381 where BT : BTSearch<F, N, Agg=Bounds<F>>,
386 G : SupportGenerator<F, N, Id=BT::Data>, 382 G : SupportGenerator<F, N, Id=BT::Data>,
387 G::SupportType : Mapping<Loc<F, N>,Codomain=F> 383 G::SupportType : Mapping<Loc<F, N>, Codomain=F> + LocalAnalysis<F, Bounds<F>, N> {
388 + LocalAnalysis<F, Bounds<F>, N> {
389 let τα = τ * self.α(); 384 let τα = τ * self.α();
390 let refinement_tolerance = ε * config.refinement.tolerance_mult; 385 let refinement_tolerance = ε * config.refinement.tolerance_mult;
391 let merge_tolerance = config.merge_tolerance_mult * ε; 386 let merge_tolerance = config.merge_tolerance_mult * ε;
392 let slack = radon_μ.norm(Radon); 387 let slack = radon_μ.norm(Radon);
393 let bnd = d.bounds(); 388 let bnd = d.bounds();
394 389
395 return { 390 return {
396 μ.both_matching(radon_μ) 391 μ.both_matching(radon_μ)
397 .all(|(α, rα, x)| { 392 .all(|(α, rα, x)| {
398 let v = d.apply(x); 393 let v = -d.apply(x); // TODO: observe ad hoc negation here, after minus_τv
394 // switch to τv.
399 let (l1, u1) = match α.partial_cmp(&0.0).unwrap_or(Equal) { 395 let (l1, u1) = match α.partial_cmp(&0.0).unwrap_or(Equal) {
400 Greater => (τα, τα), 396 Greater => (τα, τα),
401 _ => (F::NEG_INFINITY, τα), 397 _ => (F::NEG_INFINITY, τα),
402 // Less should not happen; treated as Equal 398 // Less should not happen; treated as Equal
403 }; 399 };
408 }; 404 };
409 // TODO: both fail. 405 // TODO: both fail.
410 (l1 + l2 - merge_tolerance <= v) && (v <= u1 + u2 + merge_tolerance) 406 (l1 + l2 - merge_tolerance <= v) && (v <= u1 + u2 + merge_tolerance)
411 }) 407 })
412 } && { 408 } && {
413 let keep_below = τα + slack + merge_tolerance; 409 let keep_above = -τα - slack - merge_tolerance;
414 bnd.upper() <= keep_below 410 bnd.lower() <= keep_above
415 || 411 ||
416 d.has_upper_bound(keep_below, refinement_tolerance, config.refinement.max_steps) 412 d.has_lower_bound(keep_above, refinement_tolerance, config.refinement.max_steps)
417 } 413 }
418 } 414 }
419 415
420 fn target_bounds(&self, τ : F, ε : F) -> Option<Bounds<F>> { 416 fn target_bounds(&self, τ : F, ε : F) -> Option<Bounds<F>> {
421 let τα = τ * self.α(); 417 let τα = τ * self.α();
433 where Cube<F, N> : P2Minimise<Loc<F, N>, F> { 429 where Cube<F, N> : P2Minimise<Loc<F, N>, F> {
434 430
435 fn goodness<G, BT>( 431 fn goodness<G, BT>(
436 &self, 432 &self,
437 d : &mut BTFN<F, G, BT, N>, 433 d : &mut BTFN<F, G, BT, N>,
438 _μ : &DiscreteMeasure<Loc<F, N>, F>, 434 _μ : &RNDM<F, N>,
439 y : &Loc<F, N>, 435 y : &Loc<F, N>,
440 z : &Loc<F, N>, 436 z : &Loc<F, N>,
441 τ : F, 437 τ : F,
442 ε : F, 438 ε : F,
443 _config : &FBGenericConfig<F>, 439 _config : &FBGenericConfig<F>,
444 ) -> F 440 ) -> F
445 where BT : BTSearch<F, N, Agg=Bounds<F>>, 441 where BT : BTSearch<F, N, Agg=Bounds<F>>,
446 G : SupportGenerator<F, N, Id=BT::Data>, 442 G : SupportGenerator<F, N, Id=BT::Data>,
447 G::SupportType : Mapping<Loc<F, N>,Codomain=F> 443 G::SupportType : Mapping<Loc<F, N>, Codomain=F> + LocalAnalysis<F, Bounds<F>, N> {
448 + LocalAnalysis<F, Bounds<F>, N> {
449 let w = |x| 1.0.min((ε + d.apply(x))/(τ * self.α())); 444 let w = |x| 1.0.min((ε + d.apply(x))/(τ * self.α()));
450 w(z) - w(y) 445 w(z) - w(y)
451 } 446 }
452 447
453 fn radon_norm_bound(&self, b : F) -> F { 448 fn radon_norm_bound(&self, b : F) -> F {
498 config : &FBGenericConfig<F>, 493 config : &FBGenericConfig<F>,
499 slack : F, 494 slack : F,
500 ) -> Option<(Loc<F, N>, F, bool)> 495 ) -> Option<(Loc<F, N>, F, bool)>
501 where BT : BTSearch<F, N, Agg=Bounds<F>>, 496 where BT : BTSearch<F, N, Agg=Bounds<F>>,
502 G : SupportGenerator<F, N, Id=BT::Data>, 497 G : SupportGenerator<F, N, Id=BT::Data>,
503 G::SupportType : Mapping<Loc<F, N>,Codomain=F> 498 G::SupportType : Mapping<Loc<F, N>, Codomain=F> + LocalAnalysis<F, Bounds<F>, N> {
504 + LocalAnalysis<F, Bounds<F>, N> {
505 let τα = τ * self.α(); 499 let τα = τ * self.α();
506 let keep_below = τα + slack + ε; 500 let keep_below = τα + slack + ε;
507 let keep_above = -(τα + slack) - ε; 501 let keep_above = -(τα + slack) - ε;
508 let maximise_above = τα + slack + ε * config.insertion_cutoff_factor; 502 let maximise_above = τα + slack + ε * config.insertion_cutoff_factor;
509 let minimise_below = -(τα + slack) - ε * config.insertion_cutoff_factor; 503 let minimise_below = -(τα + slack) - ε * config.insertion_cutoff_factor;
536 } 530 }
537 531
538 fn verify_merge_candidate<G, BT>( 532 fn verify_merge_candidate<G, BT>(
539 &self, 533 &self,
540 d : &mut BTFN<F, G, BT, N>, 534 d : &mut BTFN<F, G, BT, N>,
541 μ : &DiscreteMeasure<Loc<F, N>, F>, 535 μ : &RNDM<F, N>,
542 τ : F, 536 τ : F,
543 ε : F, 537 ε : F,
544 config : &FBGenericConfig<F>, 538 config : &FBGenericConfig<F>,
545 ) -> bool 539 ) -> bool
546 where BT : BTSearch<F, N, Agg=Bounds<F>>, 540 where BT : BTSearch<F, N, Agg=Bounds<F>>,
547 G : SupportGenerator<F, N, Id=BT::Data>, 541 G : SupportGenerator<F, N, Id=BT::Data>,
548 G::SupportType : Mapping<Loc<F, N>,Codomain=F> 542 G::SupportType : Mapping<Loc<F, N>,Codomain=F> + LocalAnalysis<F, Bounds<F>, N> {
549 + LocalAnalysis<F, Bounds<F>, N> {
550 let τα = τ * self.α(); 543 let τα = τ * self.α();
551 let refinement_tolerance = ε * config.refinement.tolerance_mult; 544 let refinement_tolerance = ε * config.refinement.tolerance_mult;
552 let merge_tolerance = config.merge_tolerance_mult * ε; 545 let merge_tolerance = config.merge_tolerance_mult * ε;
553 let keep_below = τα + merge_tolerance; 546 let keep_below = τα + merge_tolerance;
554 let keep_above = -τα - merge_tolerance; 547 let keep_above = -τα - merge_tolerance;
580 } 573 }
581 574
582 fn verify_merge_candidate_radonsq<G, BT>( 575 fn verify_merge_candidate_radonsq<G, BT>(
583 &self, 576 &self,
584 d : &mut BTFN<F, G, BT, N>, 577 d : &mut BTFN<F, G, BT, N>,
585 μ : &DiscreteMeasure<Loc<F, N>, F>, 578 μ : &RNDM<F, N>,
586 τ : F, 579 τ : F,
587 ε : F, 580 ε : F,
588 config : &FBGenericConfig<F>, 581 config : &FBGenericConfig<F>,
589 radon_μ : &DiscreteMeasure<Loc<F, N>, F>, 582 radon_μ : &RNDM<F, N>,
590 ) -> bool 583 ) -> bool
591 where BT : BTSearch<F, N, Agg=Bounds<F>>, 584 where BT : BTSearch<F, N, Agg=Bounds<F>>,
592 G : SupportGenerator<F, N, Id=BT::Data>, 585 G : SupportGenerator<F, N, Id=BT::Data>,
593 G::SupportType : Mapping<Loc<F, N>,Codomain=F> 586 G::SupportType : Mapping<Loc<F, N>,Codomain=F> + LocalAnalysis<F, Bounds<F>, N> {
594 + LocalAnalysis<F, Bounds<F>, N> {
595 let τα = τ * self.α(); 587 let τα = τ * self.α();
596 let refinement_tolerance = ε * config.refinement.tolerance_mult; 588 let refinement_tolerance = ε * config.refinement.tolerance_mult;
597 let merge_tolerance = config.merge_tolerance_mult * ε; 589 let merge_tolerance = config.merge_tolerance_mult * ε;
598 let slack = radon_μ.norm(Radon); 590 let slack = radon_μ.norm(Radon);
599 let bnd = d.bounds(); 591 let bnd = d.bounds();
645 where Cube<F, N> : P2Minimise<Loc<F, N>, F> { 637 where Cube<F, N> : P2Minimise<Loc<F, N>, F> {
646 638
647 fn goodness<G, BT>( 639 fn goodness<G, BT>(
648 &self, 640 &self,
649 d : &mut BTFN<F, G, BT, N>, 641 d : &mut BTFN<F, G, BT, N>,
650 _μ : &DiscreteMeasure<Loc<F, N>, F>, 642 _μ : &RNDM<F, N>,
651 y : &Loc<F, N>, 643 y : &Loc<F, N>,
652 z : &Loc<F, N>, 644 z : &Loc<F, N>,
653 τ : F, 645 τ : F,
654 ε : F, 646 ε : F,
655 _config : &FBGenericConfig<F>, 647 _config : &FBGenericConfig<F>,

Mercurial Repository: pointsource_algs

mercurial