gammalooprs/integrands/
evaluation.rs

1use std::fmt::Display;
2use std::time::Duration;
3
4use bincode::{Decode, Encode};
5use colored::Colorize;
6use serde::{Deserialize, Serialize};
7use spenso::algebra::complex::Complex;
8use tabled::{
9    Table, Tabled,
10    builder::Builder,
11    settings::{
12        Alignment, Modify, Panel, Style,
13        object::{Columns, Object, Rows},
14        style::{HorizontalLine, On, VerticalLine},
15        themes::BorderCorrection,
16    },
17};
18
19use crate::observables::{
20    EventGroupList, GenericEventGroupList, ObservableSnapshotBundle,
21    events::{format_complex_generic, format_optional_real_generic, format_real_generic},
22};
23use crate::{
24    settings::runtime::{IntegrationStatisticsSnapshot, Precision},
25    utils::{
26        ArbPrec, F, FloatLike, duration_from_secs_f64_saturating, f128, format_evaluation_time,
27        newton_solver::RadialRootDiagnostics, normalize_tabled_separator_rows,
28    },
29};
30
31#[derive(Clone, Debug)]
32pub struct GraphEvaluationResult<T: FloatLike> {
33    pub integrand_result: Complex<F<T>>,
34    pub event_groups: GenericEventGroupList<T>,
35    pub event_processing_time: Duration,
36    pub generated_event_count: usize,
37    pub accepted_event_count: usize,
38}
39
40impl<T: FloatLike> GraphEvaluationResult<T> {
41    pub fn zero(zero: F<T>) -> Self {
42        Self {
43            integrand_result: Complex::new_re(zero),
44            event_groups: GenericEventGroupList::default(),
45            event_processing_time: Duration::ZERO,
46            generated_event_count: 0,
47            accepted_event_count: 0,
48        }
49    }
50
51    pub fn merge_in_place(&mut self, mut other: Self) {
52        self.integrand_result += other.integrand_result;
53        self.event_groups.append(&mut other.event_groups);
54        self.event_processing_time += other.event_processing_time;
55        self.generated_event_count += other.generated_event_count;
56        self.accepted_event_count += other.accepted_event_count;
57    }
58
59    pub fn into_f64(self) -> GraphEvaluationResult<f64> {
60        GraphEvaluationResult {
61            integrand_result: Complex::new(
62                self.integrand_result.re.into_ff64(),
63                self.integrand_result.im.into_ff64(),
64            ),
65            event_groups: self.event_groups.to_f64(),
66            event_processing_time: self.event_processing_time,
67            generated_event_count: self.generated_event_count,
68            accepted_event_count: self.accepted_event_count,
69        }
70    }
71}
72
73/// The result of an evaluation of the integrand
74#[derive(Clone, Serialize, Debug)]
75pub struct EvaluationResult {
76    /// Integrand value before any parameterization Jacobian is applied.
77    pub integrand_result: Complex<F<f64>>,
78    pub parameterization_jacobian: Option<F<f64>>,
79    /// Monte Carlo sample weight supplied by the integrator/grid, excluding the parameterization Jacobian.
80    pub integrator_weight: F<f64>,
81    pub event_groups: EventGroupList,
82    pub evaluation_metadata: EvaluationMetaData,
83}
84
85#[derive(Clone, Serialize, Debug)]
86pub struct EvaluationResultOutput {
87    pub integrand_result: Complex<F<f64>>,
88    pub parameterization_jacobian: Option<F<f64>>,
89    pub integrator_weight: F<f64>,
90    pub event_groups: EventGroupList,
91    pub evaluation_metadata: Option<EvaluationMetaData>,
92}
93
94#[derive(Clone, Debug)]
95pub struct GenericEvaluationResult<T: FloatLike> {
96    pub integrand_result: Complex<F<T>>,
97    pub parameterization_jacobian: Option<F<T>>,
98    pub integrator_weight: F<T>,
99    pub event_groups: GenericEventGroupList<T>,
100    pub evaluation_metadata: EvaluationMetaData,
101}
102
103#[derive(Clone, Debug)]
104pub struct GenericEvaluationResultOutput<T: FloatLike> {
105    pub integrand_result: Complex<F<T>>,
106    pub parameterization_jacobian: Option<F<T>>,
107    pub integrator_weight: F<T>,
108    pub event_groups: GenericEventGroupList<T>,
109    pub evaluation_metadata: Option<EvaluationMetaData>,
110}
111
112impl<T: FloatLike> GenericEvaluationResult<T> {
113    pub fn into_output(self, minimal_output: bool) -> GenericEvaluationResultOutput<T> {
114        GenericEvaluationResultOutput {
115            integrand_result: self.integrand_result,
116            parameterization_jacobian: self.parameterization_jacobian,
117            integrator_weight: self.integrator_weight,
118            event_groups: self.event_groups,
119            evaluation_metadata: (!minimal_output).then_some(self.evaluation_metadata),
120        }
121    }
122}
123
124#[derive(Clone, Debug)]
125pub enum PreciseEvaluationResultOutput {
126    Double(GenericEvaluationResultOutput<f64>),
127    Quad(GenericEvaluationResultOutput<f128>),
128    Arb(GenericEvaluationResultOutput<ArbPrec>),
129}
130
131#[derive(Clone, Debug)]
132pub enum PreciseEvaluationResult {
133    Double(GenericEvaluationResult<f64>),
134    Quad(GenericEvaluationResult<f128>),
135    Arb(GenericEvaluationResult<ArbPrec>),
136}
137
138impl PreciseEvaluationResult {
139    pub fn into_output(self, minimal_output: bool) -> PreciseEvaluationResultOutput {
140        match self {
141            PreciseEvaluationResult::Double(result) => {
142                PreciseEvaluationResultOutput::Double(result.into_output(minimal_output))
143            }
144            PreciseEvaluationResult::Quad(result) => {
145                PreciseEvaluationResultOutput::Quad(result.into_output(minimal_output))
146            }
147            PreciseEvaluationResult::Arb(result) => {
148                PreciseEvaluationResultOutput::Arb(result.into_output(minimal_output))
149            }
150        }
151    }
152}
153
154impl PreciseEvaluationResultOutput {
155    pub fn precision(&self) -> Precision {
156        match self {
157            PreciseEvaluationResultOutput::Double(_) => Precision::Double,
158            PreciseEvaluationResultOutput::Quad(_) => Precision::Quad,
159            PreciseEvaluationResultOutput::Arb(_) => Precision::Arb,
160        }
161    }
162
163    pub fn evaluation_metadata(&self) -> Option<&EvaluationMetaData> {
164        match self {
165            PreciseEvaluationResultOutput::Double(result) => result.evaluation_metadata.as_ref(),
166            PreciseEvaluationResultOutput::Quad(result) => result.evaluation_metadata.as_ref(),
167            PreciseEvaluationResultOutput::Arb(result) => result.evaluation_metadata.as_ref(),
168        }
169    }
170}
171
172impl EvaluationResult {
173    pub fn zero() -> Self {
174        Self {
175            integrand_result: Complex::new_zero(),
176            parameterization_jacobian: None,
177            integrator_weight: F(0.0),
178            event_groups: EventGroupList::default(),
179            evaluation_metadata: EvaluationMetaData::new_empty(),
180        }
181    }
182
183    pub fn into_output(self, minimal_output: bool) -> EvaluationResultOutput {
184        EvaluationResultOutput {
185            integrand_result: self.integrand_result,
186            parameterization_jacobian: self.parameterization_jacobian,
187            integrator_weight: self.integrator_weight,
188            event_groups: self.event_groups,
189            evaluation_metadata: (!minimal_output).then_some(self.evaluation_metadata),
190        }
191    }
192}
193
194fn fmt_evaluation_result_output<T: FloatLike>(
195    f: &mut std::fmt::Formatter<'_>,
196    precision: Option<Precision>,
197    integrand_result: &Complex<F<T>>,
198    parameterization_jacobian: Option<&F<T>>,
199    integrator_weight: &F<T>,
200    event_groups: &GenericEventGroupList<T>,
201    evaluation_metadata: Option<&EvaluationMetaData>,
202) -> std::fmt::Result {
203    let mut summary_rows = Vec::with_capacity(5);
204    if let Some(precision) = precision {
205        summary_rows.push(EvaluationSummaryRow {
206            field: "precision".to_string(),
207            value: precision.to_string(),
208        });
209    }
210    summary_rows.extend([
211        EvaluationSummaryRow {
212            field: "integrand result".to_string(),
213            value: format_complex_generic(integrand_result),
214        },
215        EvaluationSummaryRow {
216            field: "parameterization jacobian".to_string(),
217            value: format_optional_real_generic(parameterization_jacobian),
218        },
219        EvaluationSummaryRow {
220            field: "integrator weight".to_string(),
221            value: format_real_generic(integrator_weight),
222        },
223        EvaluationSummaryRow {
224            field: "event groups".to_string(),
225            value: format_count(event_groups.len()),
226        },
227    ]);
228
229    writeln!(f, "{}", "Evaluation result".bold().bright_green())?;
230    writeln!(f, "{}", Table::new(summary_rows).with(Style::rounded()))?;
231    if let Some(metadata) = evaluation_metadata {
232        writeln!(f)?;
233        write!(f, "{metadata}")?;
234    }
235
236    if !event_groups.is_empty() {
237        writeln!(f)?;
238        writeln!(
239            f,
240            "{}",
241            format!(
242                "Generated {} event group(s)",
243                format_count(event_groups.len())
244            )
245            .bold()
246            .bright_cyan()
247        )?;
248        write!(f, "{event_groups}")?;
249    }
250
251    Ok(())
252}
253
254fn fmt_sample_result<D: Display>(
255    f: &mut std::fmt::Formatter<'_>,
256    evaluation: &D,
257    observables: &ObservableSnapshotBundle,
258) -> std::fmt::Result {
259    write!(f, "{evaluation}")?;
260
261    if let Some(observables_table) = summarize_observables(observables) {
262        writeln!(f)?;
263        writeln!(f)?;
264        writeln!(f, "{}", "Observable snapshots".bold().bright_magenta())?;
265        write!(f, "{observables_table}")?;
266    }
267
268    Ok(())
269}
270
271fn fmt_batch_result<D: Display>(
272    f: &mut std::fmt::Formatter<'_>,
273    samples: &[D],
274    observables: &ObservableSnapshotBundle,
275) -> std::fmt::Result {
276    let summary_rows = [EvaluationSummaryRow {
277        field: "samples".to_string(),
278        value: format_count(samples.len()),
279    }];
280
281    writeln!(f, "{}", "Batch evaluation result".bold().bright_green())?;
282    writeln!(f, "{}", Table::new(summary_rows).with(Style::rounded()))?;
283
284    if let Some(observables_table) = summarize_observables(observables) {
285        writeln!(f)?;
286        writeln!(f, "{}", "Observable snapshots".bold().bright_magenta())?;
287        writeln!(f, "{observables_table}")?;
288    }
289
290    for (sample_index, sample) in samples.iter().enumerate() {
291        writeln!(f)?;
292        if sample_index > 0 {
293            writeln!(f)?;
294        }
295        writeln!(
296            f,
297            "{}",
298            format!("Sample {sample_index}").bold().bright_cyan()
299        )?;
300        write!(f, "{sample}")?;
301    }
302
303    Ok(())
304}
305
306impl Display for EvaluationResultOutput {
307    fn fmt(&self, f: &mut std::fmt::Formatter<'_>) -> std::fmt::Result {
308        fmt_evaluation_result_output(
309            f,
310            None,
311            &self.integrand_result,
312            self.parameterization_jacobian.as_ref(),
313            &self.integrator_weight,
314            &self.event_groups,
315            self.evaluation_metadata.as_ref(),
316        )
317    }
318}
319
320impl<T: FloatLike> Display for GenericEvaluationResultOutput<T> {
321    fn fmt(&self, f: &mut std::fmt::Formatter<'_>) -> std::fmt::Result {
322        fmt_evaluation_result_output(
323            f,
324            None,
325            &self.integrand_result,
326            self.parameterization_jacobian.as_ref(),
327            &self.integrator_weight,
328            &self.event_groups,
329            self.evaluation_metadata.as_ref(),
330        )
331    }
332}
333
334impl Display for PreciseEvaluationResultOutput {
335    fn fmt(&self, f: &mut std::fmt::Formatter<'_>) -> std::fmt::Result {
336        match self {
337            PreciseEvaluationResultOutput::Double(result) => fmt_evaluation_result_output(
338                f,
339                Some(Precision::Double),
340                &result.integrand_result,
341                result.parameterization_jacobian.as_ref(),
342                &result.integrator_weight,
343                &result.event_groups,
344                result.evaluation_metadata.as_ref(),
345            ),
346            PreciseEvaluationResultOutput::Quad(result) => fmt_evaluation_result_output(
347                f,
348                Some(Precision::Quad),
349                &result.integrand_result,
350                result.parameterization_jacobian.as_ref(),
351                &result.integrator_weight,
352                &result.event_groups,
353                result.evaluation_metadata.as_ref(),
354            ),
355            PreciseEvaluationResultOutput::Arb(result) => fmt_evaluation_result_output(
356                f,
357                Some(Precision::Arb),
358                &result.integrand_result,
359                result.parameterization_jacobian.as_ref(),
360                &result.integrator_weight,
361                &result.event_groups,
362                result.evaluation_metadata.as_ref(),
363            ),
364        }
365    }
366}
367
368impl Display for SampleEvaluationResult {
369    fn fmt(&self, f: &mut std::fmt::Formatter<'_>) -> std::fmt::Result {
370        write!(f, "{}", self.evaluation)
371    }
372}
373
374#[derive(Clone, Serialize, Debug)]
375pub struct SampleEvaluationResult {
376    pub evaluation: EvaluationResultOutput,
377}
378
379#[derive(Clone, Debug)]
380pub struct PreciseSampleEvaluationResult {
381    pub evaluation: PreciseEvaluationResultOutput,
382}
383
384#[derive(Clone, Serialize, Debug)]
385pub struct SingleSampleEvaluationResult {
386    /// Evaluation record for the requested f64 sample.
387    pub sample: SampleEvaluationResult,
388    /// Persistent runtime-observable snapshot after processing the requested sample.
389    ///
390    /// The snapshot can include samples accumulated by earlier API calls.
391    pub observables: ObservableSnapshotBundle,
392}
393
394#[derive(Clone, Serialize, Debug)]
395pub struct BatchSampleEvaluationResult {
396    /// Evaluation records in the same order as the requested f64 samples.
397    pub samples: Vec<SampleEvaluationResult>,
398    /// Persistent runtime-observable snapshot after processing the requested batch.
399    ///
400    /// The snapshot can include samples accumulated by earlier API calls.
401    pub observables: ObservableSnapshotBundle,
402}
403
404#[derive(Clone, Debug)]
405pub struct RawBatchEvaluationResult {
406    pub samples: Vec<EvaluationResult>,
407    pub statistics: StatisticsCounter,
408}
409
410#[derive(Clone, Debug)]
411pub struct PreciseSingleSampleEvaluationResult {
412    /// Evaluation record retaining the numerical precision selected for the sample.
413    pub sample: PreciseSampleEvaluationResult,
414    /// Observable snapshot produced from the same Monte Carlo sample.
415    pub observables: ObservableSnapshotBundle,
416}
417
418#[derive(Clone, Debug)]
419pub struct PreciseBatchSampleEvaluationResult {
420    /// Precision-preserving evaluation records in requested sample order.
421    pub samples: Vec<PreciseSampleEvaluationResult>,
422    /// Observable snapshot accumulated across the complete sample batch.
423    pub observables: ObservableSnapshotBundle,
424}
425
426#[derive(Clone, Debug)]
427pub struct RawPreciseBatchEvaluationResult {
428    pub samples: Vec<PreciseEvaluationResult>,
429}
430
431impl Display for PreciseSampleEvaluationResult {
432    fn fmt(&self, f: &mut std::fmt::Formatter<'_>) -> std::fmt::Result {
433        write!(f, "{}", self.evaluation)
434    }
435}
436
437impl Display for SingleSampleEvaluationResult {
438    fn fmt(&self, f: &mut std::fmt::Formatter<'_>) -> std::fmt::Result {
439        fmt_sample_result(f, &self.sample.evaluation, &self.observables)
440    }
441}
442
443impl Display for BatchSampleEvaluationResult {
444    fn fmt(&self, f: &mut std::fmt::Formatter<'_>) -> std::fmt::Result {
445        fmt_batch_result(f, &self.samples, &self.observables)
446    }
447}
448
449impl Display for PreciseSingleSampleEvaluationResult {
450    fn fmt(&self, f: &mut std::fmt::Formatter<'_>) -> std::fmt::Result {
451        fmt_sample_result(f, &self.sample.evaluation, &self.observables)
452    }
453}
454
455impl Display for PreciseBatchSampleEvaluationResult {
456    fn fmt(&self, f: &mut std::fmt::Formatter<'_>) -> std::fmt::Result {
457        fmt_batch_result(f, &self.samples, &self.observables)
458    }
459}
460
461#[derive(Tabled)]
462struct EvaluationSummaryRow {
463    field: String,
464    value: String,
465}
466
467#[derive(Tabled)]
468struct HistogramSummaryRow {
469    name: String,
470    bins: usize,
471    #[tabled(rename = "in-range")]
472    in_range_entries: String,
473    phase: String,
474    range: String,
475    underflow: String,
476    overflow: String,
477}
478
479#[derive(Tabled)]
480struct StabilitySummaryRow {
481    level: String,
482    relative_accuracy: String,
483    time: String,
484    status: String,
485}
486
487fn format_duration(duration: Duration) -> String {
488    format_evaluation_time(duration)
489}
490
491fn format_count(value: usize) -> String {
492    if value < 1_000 {
493        return value.to_string();
494    }
495
496    let value = value as f64;
497    for (scale, suffix) in [
498        (1_000_000_000_f64, "B"),
499        (1_000_000_f64, "M"),
500        (1_000_f64, "K"),
501    ] {
502        if value >= scale {
503            let scaled = value / scale;
504            let precision = if scaled >= 100.0 {
505                0
506            } else if scaled >= 10.0 {
507                1
508            } else {
509                2
510            };
511            return format!("{scaled:.precision$}{suffix}");
512        }
513    }
514
515    value.round().to_string()
516}
517
518fn format_percentage(value: f64, significant_digits: usize) -> String {
519    if !value.is_finite() {
520        return "None".red().to_string();
521    }
522
523    if value == 0.0 {
524        let decimals = significant_digits.saturating_sub(1);
525        return format!("{:.*}%", decimals, 0.0);
526    }
527
528    let abs_value = value.abs();
529    let exponent = abs_value.log10().floor() as i32;
530    if exponent < -2 || exponent >= significant_digits as i32 {
531        return format!("{:.*e}%", significant_digits.saturating_sub(1), value);
532    }
533
534    let decimals = (significant_digits as i32 - exponent - 1).max(0) as usize;
535    format!("{value:.decimals$}%")
536}
537
538fn format_status_header(label: &str) -> String {
539    format!("{label:<7}").blue().bold().to_string()
540}
541
542fn format_status_key(label: &str) -> String {
543    format!("{label} :  ")
544}
545
546fn pad_status_value(value: impl Display, width: usize) -> String {
547    format!("{value:<width$}")
548}
549
550fn status_group_separator() -> VerticalLine<On, On, ()> {
551    VerticalLine::new('│').top('┬').bottom('┴')
552}
553
554fn summarize_observables(observables: &ObservableSnapshotBundle) -> Option<String> {
555    if observables.histograms.is_empty() {
556        return None;
557    }
558
559    let rows = observables
560        .histograms
561        .iter()
562        .map(|(name, histogram)| HistogramSummaryRow {
563            name: name.clone(),
564            bins: histogram.bins.len(),
565            in_range_entries: format_count(histogram.statistics.in_range_entry_count),
566            phase: format!("{:?}", histogram.phase).to_lowercase(),
567            range: match histogram.kind {
568                crate::observables::HistogramSnapshotKind::Continuous => format!(
569                    "[{:+.16e}, {:+.16e}]",
570                    histogram.x_min.unwrap_or_default(),
571                    histogram.x_max.unwrap_or_default()
572                ),
573                crate::observables::HistogramSnapshotKind::Discrete => {
574                    let min = histogram.discrete_min_bin_id.unwrap_or_default();
575                    let max = min + histogram.bins.len() as isize - 1;
576                    format!("[{}, {}]", min, max)
577                }
578            },
579            underflow: format_count(histogram.underflow_bin.entry_count),
580            overflow: format_count(histogram.overflow_bin.entry_count),
581        })
582        .collect::<Vec<_>>();
583
584    Some(Table::new(rows).with(Style::rounded()).to_string())
585}
586
587/// Per-precision evaluation details produced during stability checks.
588#[derive(Clone, Copy, Debug, Serialize, PartialEq, Eq)]
589pub enum StabilityFailureReason {
590    ErrorThreshold,
591    ZeroError,
592    WeightThreshold,
593}
594
595#[derive(Clone, Debug, Serialize)]
596pub struct RotatedEvaluation {
597    pub rotation: String,
598    pub result: Complex<F<f64>>,
599}
600
601#[derive(Clone, Debug, Serialize)]
602pub struct StabilityResult {
603    pub precision: Precision,
604    pub estimated_relative_accuracy: Option<F<f64>>,
605    pub status: StabilityStatus,
606    pub total_time: Duration,
607}
608
609#[derive(Clone, Debug, Serialize, PartialEq, Eq)]
610pub enum StabilityStatus {
611    Unknown,
612    Stable(usize),
613    Unstable(usize),
614}
615
616impl StabilityStatus {
617    pub fn from_sample_count(sample_count: usize, accepted_as_stable: bool) -> Self {
618        if sample_count <= 1 {
619            Self::Unknown
620        } else if accepted_as_stable {
621            Self::Stable(sample_count)
622        } else {
623            Self::Unstable(sample_count)
624        }
625    }
626
627    pub fn sample_count(&self) -> usize {
628        match self {
629            Self::Unknown => 1,
630            Self::Stable(sample_count) | Self::Unstable(sample_count) => *sample_count,
631        }
632    }
633
634    fn styled_label(&self) -> String {
635        let label = self.to_string();
636        match self {
637            Self::Unknown => label.yellow().to_string(),
638            Self::Stable(_) => label.green().to_string(),
639            Self::Unstable(_) => label.red().to_string(),
640        }
641    }
642}
643
644impl Display for StabilityStatus {
645    fn fmt(&self, f: &mut std::fmt::Formatter<'_>) -> std::fmt::Result {
646        match self {
647            Self::Unknown => write!(f, "Unknown(1 sample)"),
648            Self::Stable(sample_count) => write!(f, "Stable({sample_count} samples)"),
649            Self::Unstable(sample_count) => write!(f, "Unstable({sample_count} samples)"),
650        }
651    }
652}
653
654#[derive(Clone, Debug, Serialize)]
655pub struct LoopMomentaEscalationMetrics {
656    pub sum_norm: f64,
657    pub threshold: f64,
658}
659
660/// Useful metadata generated during the evaluation, this may be expanded in the future to include more information
661#[derive(Clone, Serialize, Debug)]
662pub struct EvaluationMetaData {
663    pub total_timing: Duration,
664    pub integrand_evaluation_time: Duration,
665    pub evaluator_evaluation_time: Duration,
666    pub parameterization_time: Duration,
667    pub event_processing_time: Duration,
668    pub generated_event_count: usize,
669    pub accepted_event_count: usize,
670    pub relative_instability_error: Complex<F<f64>>,
671    pub is_nan: bool,
672    pub loop_momenta_escalation: Option<LoopMomentaEscalationMetrics>,
673    pub stability_results: Vec<StabilityResult>,
674    #[serde(skip_serializing_if = "Option::is_none")]
675    pub(crate) threshold_counterterm_error: Option<String>,
676    #[serde(skip)]
677    pub(crate) radial_root_diagnostics: RadialRootDiagnostics,
678}
679
680impl Display for EvaluationMetaData {
681    fn fmt(&self, f: &mut std::fmt::Formatter<'_>) -> std::fmt::Result {
682        let rows = vec![
683            EvaluationSummaryRow {
684                field: "nan".to_string(),
685                value: self.is_nan.to_string(),
686            },
687            EvaluationSummaryRow {
688                field: "generated events".to_string(),
689                value: format_count(self.generated_event_count),
690            },
691            EvaluationSummaryRow {
692                field: "accepted events".to_string(),
693                value: format_count(self.accepted_event_count),
694            },
695            EvaluationSummaryRow {
696                field: "parameterization time".to_string(),
697                value: format_duration(self.parameterization_time),
698            },
699            EvaluationSummaryRow {
700                field: "integrand evaluation time".to_string(),
701                value: format_duration(self.integrand_evaluation_time),
702            },
703            EvaluationSummaryRow {
704                field: "evaluator evaluation time".to_string(),
705                value: format_duration(self.evaluator_evaluation_time),
706            },
707            EvaluationSummaryRow {
708                field: "event processing time".to_string(),
709                value: format_duration(self.event_processing_time),
710            },
711            EvaluationSummaryRow {
712                field: "total evaluation time".to_string(),
713                value: format_duration(self.total_timing),
714            },
715        ];
716        writeln!(f, "{}", "Evaluation metadata".bold().bright_yellow())?;
717        writeln!(f, "{}", Table::new(rows).with(Style::rounded()))?;
718
719        if !self.stability_results.is_empty() {
720            let stability_rows = self
721                .stability_results
722                .iter()
723                .map(|result| StabilitySummaryRow {
724                    level: result.precision.to_string(),
725                    relative_accuracy: format_optional_real_generic(
726                        result.estimated_relative_accuracy.as_ref(),
727                    ),
728                    time: format_duration(result.total_time),
729                    status: result.status.styled_label(),
730                })
731                .collect::<Vec<_>>();
732            writeln!(f)?;
733            writeln!(f, "{}", "Stability results".bold().bright_magenta())?;
734            write!(f, "{}", Table::new(stability_rows).with(Style::rounded()))?;
735        }
736
737        Ok(())
738    }
739}
740
741impl EvaluationMetaData {
742    pub(crate) fn new_empty() -> Self {
743        Self {
744            total_timing: Duration::ZERO,
745            integrand_evaluation_time: Duration::ZERO,
746            evaluator_evaluation_time: Duration::ZERO,
747            parameterization_time: Duration::ZERO,
748            event_processing_time: Duration::ZERO,
749            generated_event_count: 0,
750            accepted_event_count: 0,
751            relative_instability_error: Complex::new_zero(),
752            is_nan: false,
753            loop_momenta_escalation: None,
754            stability_results: Vec::new(),
755            threshold_counterterm_error: None,
756            radial_root_diagnostics: RadialRootDiagnostics::default(),
757        }
758    }
759
760    pub(crate) fn clear_threshold_counterterm_error(&mut self) {
761        self.threshold_counterterm_error = None;
762    }
763
764    pub(crate) fn record_threshold_counterterm_error(&mut self, error: impl Into<String>) {
765        if self.threshold_counterterm_error.is_none() {
766            self.threshold_counterterm_error = Some(error.into());
767        }
768    }
769
770    pub(crate) fn final_precision(&self) -> Option<Precision> {
771        self.stability_results.last().map(|result| result.precision)
772    }
773}
774
775/// This struct merges the evaluation metadata of many evaluations into a single struct
776#[derive(Debug, Copy, Clone, Serialize, Deserialize, Encode, Decode)]
777pub struct StatisticsCounter {
778    pub num_evals: usize,
779    num_sample_points: usize,
780    sum_integrand_evaluation_time: Duration,
781    sum_evaluator_evaluation_time: Duration,
782    sum_parameterization_time: Duration,
783    sum_event_time: Duration,
784    sum_integrator_overhead_time: Duration,
785    sum_total_evaluation_time: Duration,
786    sum_relative_instability_error: (F<f64>, F<f64>),
787    num_double_precision_evals: usize,
788    num_quadruple_precision_evals: usize,
789    num_arb_precision_evals: usize,
790    num_nan_evals: usize,
791    num_nan_or_unstable_evals: usize,
792    sum_generated_event_count: usize,
793    sum_accepted_event_count: usize,
794}
795
796impl StatisticsCounter {
797    /// Turn a slice of evaluation results into a statistics counter
798    pub(crate) fn from_evaluation_results(data: &[EvaluationResult]) -> Self {
799        data.iter().fold(
800            StatisticsCounter::new_empty(),
801            |mut accumulator, data_entry| {
802                accumulator.sum_integrand_evaluation_time +=
803                    data_entry.evaluation_metadata.integrand_evaluation_time;
804                accumulator.sum_evaluator_evaluation_time +=
805                    data_entry.evaluation_metadata.evaluator_evaluation_time;
806                accumulator.sum_parameterization_time +=
807                    data_entry.evaluation_metadata.parameterization_time;
808                accumulator.sum_event_time += data_entry.evaluation_metadata.event_processing_time;
809                accumulator.sum_relative_instability_error.0 +=
810                    data_entry.evaluation_metadata.relative_instability_error.re;
811                accumulator.sum_relative_instability_error.1 +=
812                    data_entry.evaluation_metadata.relative_instability_error.im;
813                accumulator.sum_total_evaluation_time +=
814                    data_entry.evaluation_metadata.total_timing;
815                accumulator.sum_generated_event_count +=
816                    data_entry.evaluation_metadata.generated_event_count;
817                accumulator.sum_accepted_event_count +=
818                    data_entry.evaluation_metadata.accepted_event_count;
819
820                accumulator.num_evals += 1;
821                match data_entry
822                    .evaluation_metadata
823                    .final_precision()
824                    .unwrap_or(Precision::Double)
825                {
826                    Precision::Double => accumulator.num_double_precision_evals += 1,
827                    Precision::Quad => accumulator.num_quadruple_precision_evals += 1,
828                    Precision::Arb => accumulator.num_arb_precision_evals += 1,
829                    // _ => (),
830                }
831
832                if data_entry.evaluation_metadata.is_nan {
833                    accumulator.num_nan_evals += 1;
834                }
835                if data_entry.evaluation_metadata.is_nan
836                    || data_entry
837                        .evaluation_metadata
838                        .stability_results
839                        .last()
840                        .map(|result| matches!(result.status, StabilityStatus::Unstable(_)))
841                        .unwrap_or(false)
842                {
843                    accumulator.num_nan_or_unstable_evals += 1;
844                }
845
846                accumulator
847            },
848        )
849    }
850
851    /// Merge two statistics counters into a single one, but keeping the original ones unchanged
852    pub(crate) fn merged(&self, other: &Self) -> Self {
853        Self {
854            sum_integrand_evaluation_time: self.sum_integrand_evaluation_time
855                + other.sum_integrand_evaluation_time,
856            sum_evaluator_evaluation_time: self.sum_evaluator_evaluation_time
857                + other.sum_evaluator_evaluation_time,
858            sum_parameterization_time: self.sum_parameterization_time
859                + other.sum_parameterization_time,
860            sum_integrator_overhead_time: self.sum_integrator_overhead_time
861                + other.sum_integrator_overhead_time,
862            sum_relative_instability_error: (
863                self.sum_relative_instability_error.0 + other.sum_relative_instability_error.0,
864                self.sum_relative_instability_error.1 + other.sum_relative_instability_error.1,
865            ),
866            num_evals: self.num_evals + other.num_evals,
867            num_sample_points: self.num_sample_points + other.num_sample_points,
868            num_double_precision_evals: (self.num_double_precision_evals
869                + other.num_double_precision_evals),
870            num_quadruple_precision_evals: (self.num_quadruple_precision_evals
871                + other.num_quadruple_precision_evals),
872            num_arb_precision_evals: self.num_arb_precision_evals + other.num_arb_precision_evals,
873            sum_total_evaluation_time: self.sum_total_evaluation_time
874                + other.sum_total_evaluation_time,
875            num_nan_evals: self.num_nan_evals + other.num_nan_evals,
876            num_nan_or_unstable_evals: self.num_nan_or_unstable_evals
877                + other.num_nan_or_unstable_evals,
878            sum_event_time: self.sum_event_time + other.sum_event_time,
879            sum_generated_event_count: self.sum_generated_event_count
880                + other.sum_generated_event_count,
881            sum_accepted_event_count: self.sum_accepted_event_count
882                + other.sum_accepted_event_count,
883        }
884    }
885
886    pub(crate) fn new_empty() -> Self {
887        Self {
888            sum_integrand_evaluation_time: Duration::ZERO,
889            sum_evaluator_evaluation_time: Duration::ZERO,
890            sum_parameterization_time: Duration::ZERO,
891            sum_event_time: Duration::ZERO,
892            sum_integrator_overhead_time: Duration::ZERO,
893            sum_relative_instability_error: (F(0.0), F(0.0)),
894            sum_total_evaluation_time: Duration::ZERO,
895            num_evals: 0,
896            num_sample_points: 0,
897            num_double_precision_evals: 0,
898            num_quadruple_precision_evals: 0,
899            num_arb_precision_evals: 0,
900            num_nan_evals: 0,
901            num_nan_or_unstable_evals: 0,
902            sum_generated_event_count: 0,
903            sum_accepted_event_count: 0,
904        }
905    }
906
907    /// Merge a list of statistics counters into a single one.
908    #[allow(dead_code)]
909    pub(crate) fn merge_list(list: Vec<Self>) -> Self {
910        if let Some(merged) = list.into_iter().reduce(|acc, x| acc.merged(&x)) {
911            merged
912        } else {
913            Self::new_empty()
914        }
915    }
916
917    fn avg_duration(sum: Duration, count: usize) -> Duration {
918        if count == 0 {
919            return Duration::ZERO;
920        }
921
922        duration_from_secs_f64_saturating(sum.as_secs_f64() / count as f64)
923    }
924
925    fn eval_percentage(&self, count: usize) -> f64 {
926        if self.num_evals == 0 {
927            0.0
928        } else {
929            count as f64 / self.num_evals as f64 * 100.0
930        }
931    }
932
933    pub(crate) fn add_integrator_overhead(&mut self, duration: Duration, sample_points: usize) {
934        self.sum_integrator_overhead_time += duration;
935        self.num_sample_points += sample_points;
936    }
937
938    /// Compute the average time spent in the evaluate_sample function.
939    pub(crate) fn get_avg_total_timing(&self) -> Duration {
940        Self::avg_duration(self.sum_total_evaluation_time, self.num_evals)
941    }
942
943    /// Compute the average time spent in the original integrand evaluation call.
944    pub(crate) fn get_avg_integrand_timing(&self) -> Duration {
945        Self::avg_duration(self.sum_integrand_evaluation_time, self.num_evals)
946    }
947
948    /// Compute the average time spent inside Symbolica evaluator function calls.
949    pub(crate) fn get_avg_evaluator_timing(&self) -> Duration {
950        Self::avg_duration(self.sum_evaluator_evaluation_time, self.num_evals)
951    }
952
953    /// Compute the average time spent in the parameterization of the integrand. Especaially useful for monitoring the performance of tropical sampling.
954    pub(crate) fn get_avg_param_timing(&self) -> Duration {
955        Self::avg_duration(self.sum_parameterization_time, self.num_evals)
956    }
957
958    pub(crate) fn get_avg_event_timing(&self) -> Duration {
959        Self::avg_duration(self.sum_event_time, self.num_evals)
960    }
961
962    pub(crate) fn get_avg_integrator_overhead_timing(&self) -> Duration {
963        Self::avg_duration(self.sum_integrator_overhead_time, self.num_sample_points)
964    }
965
966    /// Get the average relative error computed during instability checks
967    #[allow(dead_code)]
968    pub(crate) fn get_avg_instabillity_error(&self) -> (F<f64>, F<f64>) {
969        if self.num_evals == 0 {
970            return (
971                self.sum_relative_instability_error.0.zero(),
972                self.sum_relative_instability_error.1.zero(),
973            );
974        }
975
976        (
977            self.sum_relative_instability_error.0
978                / self
979                    .sum_relative_instability_error
980                    .0
981                    .from_usize(self.num_evals),
982            self.sum_relative_instability_error.1
983                / self
984                    .sum_relative_instability_error
985                    .1
986                    .from_usize(self.num_evals),
987        )
988    }
989
990    /// Get the percentage of evaluations that were done in double precision and were stable.
991    pub(crate) fn get_percentage_f64(&self) -> f64 {
992        self.eval_percentage(self.num_double_precision_evals)
993    }
994
995    /// Get the percentage of evaluations that went to quadruple precision and were stable.
996    pub(crate) fn get_percentage_f128(&self) -> f64 {
997        self.eval_percentage(self.num_quadruple_precision_evals)
998    }
999
1000    pub(crate) fn get_percentage_arb(&self) -> f64 {
1001        self.eval_percentage(self.num_arb_precision_evals)
1002    }
1003
1004    pub(crate) fn get_percentage_nan(&self) -> f64 {
1005        self.eval_percentage(self.num_nan_evals)
1006    }
1007
1008    pub(crate) fn get_percentage_nan_or_unstable(&self) -> f64 {
1009        self.eval_percentage(self.num_nan_or_unstable_evals)
1010    }
1011
1012    pub(crate) fn selection_efficiency_percentage(&self) -> Option<f64> {
1013        if self.sum_generated_event_count == 0 {
1014            None
1015        } else {
1016            Some(
1017                self.sum_accepted_event_count as f64 / self.sum_generated_event_count as f64
1018                    * 100.0,
1019            )
1020        }
1021    }
1022
1023    pub(crate) fn snapshot(&self) -> IntegrationStatisticsSnapshot {
1024        IntegrationStatisticsSnapshot {
1025            num_evals: self.num_evals,
1026            average_total_time_seconds: self.get_avg_total_timing().as_secs_f64(),
1027            average_parameterization_time_seconds: self.get_avg_param_timing().as_secs_f64(),
1028            average_integrand_time_seconds: self.get_avg_integrand_timing().as_secs_f64(),
1029            average_evaluator_time_seconds: self.get_avg_evaluator_timing().as_secs_f64(),
1030            average_observable_time_seconds: self.get_avg_event_timing().as_secs_f64(),
1031            average_integrator_time_seconds: self
1032                .get_avg_integrator_overhead_timing()
1033                .as_secs_f64(),
1034            f64_percentage: self.get_percentage_f64(),
1035            f128_percentage: self.get_percentage_f128(),
1036            arb_percentage: self.get_percentage_arb(),
1037            nan_percentage: self.get_percentage_nan(),
1038            nan_or_unstable_percentage: self.get_percentage_nan_or_unstable(),
1039            generated_event_count: self.sum_generated_event_count,
1040            accepted_event_count: self.sum_accepted_event_count,
1041            selection_efficiency_percentage: self.selection_efficiency_percentage(),
1042        }
1043    }
1044
1045    pub(crate) fn build_status_table(&self) -> Table {
1046        let time_integrand_formatted = format_evaluation_time(self.get_avg_integrand_timing());
1047        let time_evaluators_formatted = format_evaluation_time(self.get_avg_evaluator_timing());
1048        let param_time_formatted = format_evaluation_time(self.get_avg_param_timing());
1049        let event_time_formatted = format_evaluation_time(self.get_avg_event_timing());
1050        let integrator_time_formatted =
1051            format_evaluation_time(self.get_avg_integrator_overhead_timing());
1052        let total_time = format_evaluation_time(self.get_avg_total_timing());
1053        let selection_efficiency = self.selection_efficiency_percentage();
1054        let selection_efficiency_display = selection_efficiency
1055            .map(|value| {
1056                pad_status_value(format_percentage(value, 3), 9)
1057                    .green()
1058                    .to_string()
1059            })
1060            .unwrap_or_else(|| pad_status_value("N/A", 9));
1061        let nan_or_unstable = self.get_percentage_nan_or_unstable();
1062        let nan_or_unstable_display = if nan_or_unstable > 0.0 {
1063            pad_status_value(format_percentage(nan_or_unstable, 2), 9)
1064                .red()
1065                .to_string()
1066        } else {
1067            pad_status_value(format_percentage(nan_or_unstable, 2), 9)
1068                .green()
1069                .to_string()
1070        };
1071        let mut table = Builder::new();
1072        table.push_record([
1073            format_status_header("timing"),
1074            format_status_key("total"),
1075            pad_status_value(total_time, 9).green().to_string(),
1076            format_status_key("param"),
1077            pad_status_value(param_time_formatted, 9)
1078                .green()
1079                .to_string(),
1080            format_status_key("itg"),
1081            pad_status_value(time_integrand_formatted, 9)
1082                .green()
1083                .to_string(),
1084            format_status_key("evaluators"),
1085            pad_status_value(time_evaluators_formatted, 9)
1086                .green()
1087                .to_string(),
1088        ]);
1089        table.push_record([
1090            format_status_header("evals"),
1091            format_status_key("f64"),
1092            pad_status_value(format!("{:.2}%", self.get_percentage_f64()), 9)
1093                .green()
1094                .to_string(),
1095            format_status_key("f128"),
1096            pad_status_value(format!("{:.2}%", self.get_percentage_f128()), 9)
1097                .green()
1098                .to_string(),
1099            format_status_key("arb"),
1100            pad_status_value(format!("{:.2}%", self.get_percentage_arb()), 9)
1101                .green()
1102                .to_string(),
1103            format_status_key("nans+unstable"),
1104            nan_or_unstable_display.clone(),
1105        ]);
1106        table.push_record([
1107            format_status_header("events"),
1108            format_status_key("evts #"),
1109            pad_status_value(format_count(self.sum_generated_event_count), 9)
1110                .green()
1111                .to_string(),
1112            format_status_key("sel. %"),
1113            selection_efficiency_display,
1114            format_status_key("obs"),
1115            pad_status_value(event_time_formatted, 9)
1116                .green()
1117                .to_string(),
1118            format_status_key("integrator"),
1119            pad_status_value(integrator_time_formatted, 9)
1120                .green()
1121                .to_string(),
1122        ]);
1123
1124        let mut table = table.build();
1125        table.with(Panel::header(
1126            "Integration statistics".bold().green().to_string(),
1127        ));
1128        table.with(
1129            Style::rounded()
1130                .remove_horizontals()
1131                .verticals([
1132                    (1, status_group_separator()),
1133                    (3, status_group_separator()),
1134                    (5, status_group_separator()),
1135                    (7, status_group_separator()),
1136                ])
1137                .horizontals([(
1138                    1,
1139                    HorizontalLine::new('─')
1140                        .intersection('┬')
1141                        .left('├')
1142                        .right('┤'),
1143                )])
1144                .remove_vertical(),
1145        );
1146        table.with(BorderCorrection::span());
1147        table.with(Modify::new(Rows::new(0..1)).with(Alignment::center()));
1148        table.with(Modify::new(Rows::new(1..)).with(Alignment::left()));
1149        for column in [1usize, 3, 5, 7] {
1150            table.with(
1151                Modify::new(Rows::new(1..).intersect(Columns::one(column)))
1152                    .with(Alignment::right()),
1153            );
1154        }
1155        table
1156    }
1157
1158    pub(crate) fn render_status_table(&self) -> String {
1159        normalize_tabled_separator_rows(&self.build_status_table().to_string())
1160    }
1161}
1162
1163#[cfg(test)]
1164mod tests {
1165    use std::time::Duration;
1166
1167    use crate::settings::runtime::Precision;
1168
1169    use super::{EvaluationMetaData, StabilityResult, StabilityStatus, StatisticsCounter};
1170
1171    #[test]
1172    fn status_table_renders_three_rows_without_header() {
1173        let stats = StatisticsCounter {
1174            num_evals: 10,
1175            num_sample_points: 10,
1176            sum_integrand_evaluation_time: Duration::from_micros(414),
1177            sum_evaluator_evaluation_time: Duration::from_micros(279),
1178            sum_parameterization_time: Duration::from_nanos(5_800),
1179            sum_event_time: Duration::ZERO,
1180            sum_integrator_overhead_time: Duration::from_micros(110),
1181            sum_total_evaluation_time: Duration::from_micros(462),
1182            sum_relative_instability_error: (0.0.into(), 0.0.into()),
1183            num_double_precision_evals: 10,
1184            num_quadruple_precision_evals: 0,
1185            num_arb_precision_evals: 0,
1186            num_nan_evals: 0,
1187            num_nan_or_unstable_evals: 0,
1188            sum_generated_event_count: 0,
1189            sum_accepted_event_count: 0,
1190        };
1191
1192        let rendered = stats.render_status_table();
1193        let lines = rendered.lines().collect::<Vec<_>>();
1194
1195        assert!(
1196            rendered.starts_with('╭') && rendered.ends_with('╯'),
1197            "{rendered}"
1198        );
1199        assert_eq!(lines.len(), 7, "{rendered}");
1200        assert!(rendered.contains("Integration statistics"), "{rendered}");
1201        assert!(lines[2].contains('┬'), "{rendered}");
1202        assert!(lines[6].contains('┴'), "{rendered}");
1203        assert_eq!(lines[3].matches('│').count(), 6, "{rendered}");
1204        assert_eq!(lines[4].matches('│').count(), 6, "{rendered}");
1205        assert_eq!(lines[5].matches('│').count(), 6, "{rendered}");
1206        assert!(rendered.contains("timing"), "{rendered}");
1207        assert!(rendered.contains("evals"), "{rendered}");
1208        assert!(rendered.contains("events"), "{rendered}");
1209        assert!(rendered.contains("nans+unstable :"), "{rendered}");
1210        assert!(rendered.contains("integrator :"), "{rendered}");
1211    }
1212
1213    #[test]
1214    fn statistics_snapshot_reports_selection_efficiency_from_event_counts() {
1215        let stats = StatisticsCounter {
1216            num_evals: 2,
1217            num_sample_points: 2,
1218            sum_integrand_evaluation_time: Duration::ZERO,
1219            sum_evaluator_evaluation_time: Duration::ZERO,
1220            sum_parameterization_time: Duration::ZERO,
1221            sum_event_time: Duration::ZERO,
1222            sum_integrator_overhead_time: Duration::ZERO,
1223            sum_total_evaluation_time: Duration::ZERO,
1224            sum_relative_instability_error: (0.0.into(), 0.0.into()),
1225            num_double_precision_evals: 2,
1226            num_quadruple_precision_evals: 0,
1227            num_arb_precision_evals: 0,
1228            num_nan_evals: 0,
1229            num_nan_or_unstable_evals: 0,
1230            sum_generated_event_count: 10,
1231            sum_accepted_event_count: 4,
1232        };
1233
1234        let snapshot = stats.snapshot();
1235
1236        assert_eq!(snapshot.generated_event_count, 10);
1237        assert_eq!(snapshot.accepted_event_count, 4);
1238        assert_eq!(snapshot.selection_efficiency_percentage, Some(40.0));
1239
1240        let rendered = stats.render_status_table();
1241        assert!(rendered.contains("sel. % :     40.0%"), "{rendered}");
1242    }
1243
1244    #[test]
1245    fn evaluation_metadata_renders_stability_status_with_sample_counts() {
1246        let metadata = EvaluationMetaData {
1247            stability_results: vec![
1248                StabilityResult {
1249                    precision: Precision::Double,
1250                    estimated_relative_accuracy: Some(1.0e-5.into()),
1251                    status: StabilityStatus::Stable(3),
1252                    total_time: Duration::from_micros(120),
1253                },
1254                StabilityResult {
1255                    precision: Precision::Quad,
1256                    estimated_relative_accuracy: None,
1257                    status: StabilityStatus::Unknown,
1258                    total_time: Duration::from_micros(240),
1259                },
1260            ],
1261            ..EvaluationMetaData::new_empty()
1262        };
1263
1264        let rendered = metadata.to_string();
1265        assert!(rendered.contains("Stable(3 samples)"), "{rendered}");
1266        assert!(rendered.contains("Unknown(1 sample)"), "{rendered}");
1267        assert!(rendered.contains("None"), "{rendered}");
1268    }
1269
1270    mod failing {
1271        use super::*;
1272
1273        #[test]
1274        fn zero_eval_statistics_snapshot_and_render_are_sanitized() {
1275            let stats = StatisticsCounter::new_empty();
1276            let snapshot = stats.snapshot();
1277
1278            assert_eq!(snapshot.num_evals, 0);
1279            assert_eq!(snapshot.average_total_time_seconds, 0.0);
1280            assert_eq!(snapshot.average_parameterization_time_seconds, 0.0);
1281            assert_eq!(snapshot.average_integrand_time_seconds, 0.0);
1282            assert_eq!(snapshot.average_evaluator_time_seconds, 0.0);
1283            assert_eq!(snapshot.average_observable_time_seconds, 0.0);
1284            assert_eq!(snapshot.average_integrator_time_seconds, 0.0);
1285            assert_eq!(snapshot.f64_percentage, 0.0);
1286            assert_eq!(snapshot.f128_percentage, 0.0);
1287            assert_eq!(snapshot.arb_percentage, 0.0);
1288            assert_eq!(snapshot.nan_percentage, 0.0);
1289            assert_eq!(snapshot.nan_or_unstable_percentage, 0.0);
1290            assert_eq!(snapshot.selection_efficiency_percentage, None);
1291
1292            let rendered = stats.render_status_table();
1293            assert!(rendered.contains("sel. % :     N/A"), "{rendered}");
1294            assert!(rendered.contains("f64 :     0.00%"), "{rendered}");
1295            assert!(rendered.contains("nans+unstable :     0.0%"), "{rendered}");
1296            assert!(rendered.contains("integrator :   0.00 ns"), "{rendered}");
1297        }
1298    }
1299}
gammalooprs/integrands/evaluation.rs

gammalooprs/integrands/
evaluation.rs
