gammalooprs/momentum/
signature.rs

1#![allow(dead_code)]
2
3use crate::momentum::sample::{
4    ExternalFourMomenta, ExternalIndex, ExternalThreeMomenta, LoopIndex, LoopMomenta,
5};
6use crate::momentum::{FourMomentum, SignOrZero, ThreeMomentum};
7use crate::utils::{F, FloatLike, Length};
8use bincode::{BorrowDecode, Decode, Encode};
9use serde::{Deserialize, Serialize};
10use spenso::algebra::algebraic_traits::RefZero;
11use std::fmt::Display;
12use std::ops::{Add, AddAssign, Index, IndexMut, Neg, SubAssign};
13use symbolica::atom::{Atom, AtomOrView, FunctionBuilder, Symbol};
14use typed_index_collections::TiVec;
15
16#[derive(Debug, Clone, Serialize, Deserialize, PartialEq, PartialOrd, Eq, Ord, Hash)]
17pub struct SignatureLike<T: From<usize>>(TiVec<T, SignOrZero>);
18pub type LoopSignature = SignatureLike<LoopIndex>;
19pub type ExternalSignature = SignatureLike<ExternalIndex>;
20
21// manual implementations because TiVec is not Encode/Devode
22impl<T: Encode + From<usize>> Encode for SignatureLike<T> {
23    fn encode<E: bincode::enc::Encoder>(
24        &self,
25        encoder: &mut E,
26    ) -> Result<(), bincode::error::EncodeError> {
27        self.0.raw.encode(encoder)
28    }
29}
30
31impl<'de, Context, T: Decode<Context> + From<usize>> BorrowDecode<'de, Context>
32    for SignatureLike<T>
33{
34    fn borrow_decode<D: bincode::de::Decoder>(
35        decoder: &mut D,
36    ) -> Result<Self, bincode::error::DecodeError> {
37        Ok(SignatureLike(Vec::decode(decoder)?.into()))
38    }
39}
40
41impl<Context, T: Decode<Context> + From<usize>> Decode<Context> for SignatureLike<T> {
42    fn decode<D: bincode::de::Decoder>(
43        decoder: &mut D,
44    ) -> Result<Self, bincode::error::DecodeError> {
45        Ok(SignatureLike(Vec::decode(decoder)?.into()))
46    }
47}
48
49#[derive(
50    Debug, Clone, Serialize, Deserialize, PartialEq, PartialOrd, Eq, Ord, Hash, Encode, Decode,
51)]
52pub struct LoopExtSignature {
53    pub internal: LoopSignature,
54    pub external: ExternalSignature,
55}
56
57impl LoopExtSignature {
58    pub(crate) fn swap_loops(&mut self, i: LoopIndex, j: LoopIndex) {
59        // println!("i{i},j{j}");÷
60        if !self.internal.is_empty() {
61            self.internal.0.swap(i, j);
62        }
63    }
64
65    pub(crate) fn put_loop_to_ext(&mut self, l: LoopIndex) {
66        let a = self.internal.0.remove(l);
67        self.external.0.push(a);
68    }
69    pub(crate) fn swap_external(&mut self, i: ExternalIndex, j: ExternalIndex) {
70        self.external.0.swap(i, j);
71    }
72
73    pub(crate) fn loop_atom<'a, I>(
74        &self,
75        mom_symbol: Symbol,
76        additional_args: &'a [I],
77        id_map: impl Fn(LoopIndex) -> Atom,
78    ) -> Atom
79    where
80        &'a I: Into<AtomOrView<'a>>,
81    {
82        self.internal.atom(mom_symbol, additional_args, id_map)
83    }
84
85    pub(crate) fn ext_atom<'a, I>(
86        &self,
87        mom_symbol: Symbol,
88        additional_args: &'a [I],
89        id_map: impl Fn(ExternalIndex) -> Atom,
90    ) -> Atom
91    where
92        &'a I: Into<AtomOrView<'a>>,
93    {
94        self.external.atom(mom_symbol, additional_args, id_map)
95    }
96
97    pub(crate) fn equality_up_to_sign(&self, other: &Self) -> bool {
98        self == other
99            || (self.internal.len() == other.internal.len()
100                && self.external.len() == other.external.len()
101                && self
102                    .internal
103                    .iter()
104                    .zip(other.internal.iter())
105                    .all(|(a, b)| *a == -*b)
106                && self
107                    .external
108                    .iter()
109                    .zip(other.external.iter())
110                    .all(|(a, b)| *a == -*b))
111    }
112
113    /// Sampling may identify energies in an exact represented spatial frame.
114    /// This does not change structural routing equality used by raised edges.
115    pub(crate) fn spatial_equality_up_to_sign<T: FloatLike>(
116        &self,
117        other: &Self,
118        externals: &ExternalThreeMomenta<F<T>>,
119    ) -> eyre::Result<bool> {
120        use crate::integrands::process::sampling::maps::SamplingEvaluationError;
121        use rug::{Float, float::Round};
122        const PRECISION: u32 = 2048;
123        if self.internal.len() != other.internal.len()
124            || self.external.len() != other.external.len()
125            || self.external.len() != externals.len()
126        {
127            return Err(eyre::eyre!(
128                "spatial route comparison requires equal routing dimensions and the complete external frame"
129            ));
130        }
131        if externals
132            .iter()
133            .flat_map(|p| [&p.px, &p.py, &p.pz])
134            .any(|value| !value.0.is_finite())
135        {
136            return Err(SamplingEvaluationError::Unrepresentable {
137                operation: "spatial route comparison",
138                detail: "nonfinite original external component".into(),
139            }
140            .into());
141        }
142        if self.equality_up_to_sign(other) {
143            return Ok(true);
144        }
145        let mut uncertain = false;
146        for sign in [1, -1] {
147            if !self
148                .internal
149                .iter()
150                .zip(other.internal.iter())
151                .all(|(a, b)| *a == if sign == 1 { *b } else { -*b })
152            {
153                continue;
154            }
155            // Differences can be +/-2, which SignatureLike cannot represent.
156            let coefficients = self
157                .external
158                .to_momtrop_format()
159                .into_iter()
160                .zip(other.external.to_momtrop_format())
161                .map(|(a, b)| a - sign * b)
162                .collect::<Vec<_>>();
163            let mut sums: [[Float; 2]; 3] = std::array::from_fn(|_| {
164                [Float::with_val(PRECISION, 0), Float::with_val(PRECISION, 0)]
165            });
166            for (&coefficient, momentum) in coefficients.iter().zip(externals) {
167                if coefficient == 0 {
168                    continue;
169                }
170                for (sum, component) in
171                    sums.iter_mut()
172                        .zip([&momentum.px, &momentum.py, &momentum.pz])
173                {
174                    let (mut lo, mut hi) = component.0.mpfr_enclosure(PRECISION);
175                    if coefficient < 0 {
176                        std::mem::swap(&mut lo, &mut hi);
177                    }
178                    let lo = Float::with_val_round(PRECISION, &lo * coefficient, Round::Down).0;
179                    let hi = Float::with_val_round(PRECISION, &hi * coefficient, Round::Up).0;
180                    sum[0] = Float::with_val_round(PRECISION, &sum[0] + &lo, Round::Down).0;
181                    sum[1] = Float::with_val_round(PRECISION, &sum[1] + &hi, Round::Up).0;
182                }
183            }
184            if sums.iter().flatten().any(|value| !value.is_finite()) {
185                return Err(SamplingEvaluationError::Unrepresentable {
186                    operation: "spatial route comparison",
187                    detail: "nonfinite directed external sum".into(),
188                }
189                .into());
190            }
191            if sums.iter().any(|[lo, hi]| lo > &0 || hi < &0) {
192                continue;
193            }
194            if sums.iter().all(|[lo, hi]| lo == &0 && hi == &0) {
195                return Ok(true);
196            }
197            uncertain = true;
198        }
199        if uncertain {
200            return Err(SamplingEvaluationError::UncertainGeometry {
201                detail: format!("spatial routing equality is unresolved at {PRECISION}-bit directed precision: {self:?} versus {other:?}"),
202            }
203            .into());
204        }
205        Ok(false)
206    }
207}
208
209impl From<(Vec<isize>, Vec<isize>)> for LoopExtSignature {
210    fn from(value: (Vec<isize>, Vec<isize>)) -> Self {
211        Self {
212            internal: LoopSignature::from_iter(value.0),
213            external: ExternalSignature::from_iter(value.1),
214        }
215    }
216}
217
218impl<T> Index<T> for SignatureLike<T>
219where
220    usize: From<T>,
221    T: From<usize>,
222{
223    type Output = SignOrZero;
224
225    fn index(&self, index: T) -> &Self::Output {
226        &self.0[index]
227    }
228}
229
230impl<T> IndexMut<T> for SignatureLike<T>
231where
232    usize: From<T>,
233    T: From<usize>,
234{
235    fn index_mut(&mut self, index: T) -> &mut Self::Output {
236        &mut self.0[index]
237    }
238}
239
240impl<T> Display for SignatureLike<T>
241where
242    T: From<usize>,
243{
244    fn fmt(&self, f: &mut std::fmt::Formatter<'_>) -> std::fmt::Result {
245        for sign in &self.0 {
246            write!(f, "{}", sign)?;
247        }
248        Ok(())
249    }
250}
251
252impl<T: From<usize>> Default for SignatureLike<T> {
253    fn default() -> Self {
254        SignatureLike(TiVec::new())
255    }
256}
257
258impl<T> FromIterator<SignOrZero> for SignatureLike<T>
259where
260    T: From<usize>,
261{
262    fn from_iter<I: IntoIterator<Item = SignOrZero>>(iter: I) -> Self {
263        SignatureLike(iter.into_iter().collect())
264    }
265}
266
267impl<T> FromIterator<i8> for SignatureLike<T>
268where
269    T: From<usize>,
270{
271    fn from_iter<I: IntoIterator<Item = i8>>(iter: I) -> Self {
272        SignatureLike(
273            iter.into_iter()
274                .map(|x| match x {
275                    0 => SignOrZero::Zero,
276                    1 => SignOrZero::Plus,
277                    -1 => SignOrZero::Minus,
278                    _ => panic!("Invalid value for Signature"),
279                })
280                .collect(),
281        )
282    }
283}
284
285impl<T> FromIterator<isize> for SignatureLike<T>
286where
287    T: From<usize>,
288{
289    fn from_iter<I: IntoIterator<Item = isize>>(iter: I) -> Self {
290        SignatureLike(
291            iter.into_iter()
292                .map(|x| match x {
293                    0 => SignOrZero::Zero,
294                    1 => SignOrZero::Plus,
295                    -1 => SignOrZero::Minus,
296                    _ => panic!("Invalid value for Signature"),
297                })
298                .collect(),
299        )
300    }
301}
302
303impl<T> From<Vec<i8>> for SignatureLike<T>
304where
305    T: From<usize>,
306{
307    fn from(value: Vec<i8>) -> Self {
308        SignatureLike::from_iter(value)
309    }
310}
311
312impl<T> IntoIterator for SignatureLike<T>
313where
314    T: From<usize>,
315{
316    type Item = SignOrZero;
317    type IntoIter = std::vec::IntoIter<Self::Item>;
318
319    fn into_iter(self) -> Self::IntoIter {
320        self.0.into_iter()
321    }
322}
323
324impl<'a, T> IntoIterator for &'a SignatureLike<T>
325where
326    T: From<usize>,
327{
328    type Item = SignOrZero;
329    type IntoIter = std::iter::Copied<std::slice::Iter<'a, Self::Item>>;
330
331    fn into_iter(self) -> Self::IntoIter {
332        self.0.iter().copied()
333    }
334}
335
336impl<T> AddAssign<()> for SignatureLike<T>
337where
338    T: From<usize> + Copy,
339    usize: From<T>,
340{
341    fn add_assign(&mut self, _rhs: ()) {
342        self.0.push(SignOrZero::Plus);
343    }
344}
345
346impl<T> SignatureLike<T>
347where
348    T: From<usize> + Copy,
349    usize: From<T>,
350{
351    pub(crate) fn atom<'a, I>(
352        &self,
353        mom_symbol: Symbol,
354        additional_args: &'a [I],
355        id_map: impl Fn(T) -> Atom,
356    ) -> Atom
357    where
358        &'a I: Into<AtomOrView<'a>>,
359    {
360        let mut rep = Atom::Zero;
361        for (l, s) in self.iter_enumerated() {
362            let mom = FunctionBuilder::new(mom_symbol)
363                .add_arg(id_map(l))
364                .add_args(additional_args)
365                .finish();
366            // println!("mom: {mom} {s}");
367
368            rep += *s * mom;
369        }
370
371        // println!("rep{rep}");
372
373        rep
374    }
375
376    pub(crate) fn iter_enumerated(&self) -> impl Iterator<Item = (T, &SignOrZero)> {
377        self.0.iter_enumerated()
378    }
379    pub(crate) fn validate_basis<B>(&self, basis: &[B]) -> bool {
380        self.len() == basis.len()
381    }
382
383    pub(crate) fn sum(&mut self, other: &Self) {
384        for (i, sign) in other.iter_enumerated() {
385            match (self[i], sign) {
386                (SignOrZero::Zero, SignOrZero::Zero) => self.0[i] = SignOrZero::Zero,
387                (SignOrZero::Zero, SignOrZero::Plus) => self.0[i] = SignOrZero::Plus,
388                (SignOrZero::Zero, SignOrZero::Minus) => self.0[i] = SignOrZero::Minus,
389                (SignOrZero::Plus, SignOrZero::Zero) => self.0[i] = SignOrZero::Plus,
390                (SignOrZero::Plus, SignOrZero::Plus) => panic!("cannot add two positive signs"),
391                (SignOrZero::Plus, SignOrZero::Minus) => self.0[i] = SignOrZero::Zero,
392                (SignOrZero::Minus, SignOrZero::Zero) => self.0[i] = SignOrZero::Minus,
393                (SignOrZero::Minus, SignOrZero::Plus) => self.0[i] = SignOrZero::Zero,
394                (SignOrZero::Minus, SignOrZero::Minus) => panic!("cannot add two negative signs"),
395            }
396        }
397    }
398
399    pub(crate) fn panic_validate_basis<B>(&self, basis: &[B]) {
400        if !self.validate_basis(basis) {
401            panic!(
402                "Invalid basis for Signature, expected length {}, got length {}",
403                self.len(),
404                basis.len()
405            );
406        }
407    }
408
409    pub(crate) fn to_momtrop_format(&self) -> Vec<isize> {
410        self.0
411            .iter()
412            .map(|x| match x {
413                SignOrZero::Zero => 0,
414                SignOrZero::Plus => 1,
415                SignOrZero::Minus => -1,
416            })
417            .collect()
418    }
419    pub(crate) fn len(&self) -> usize {
420        self.0.len()
421    }
422
423    pub fn iter(&'_ self) -> std::slice::Iter<'_, SignOrZero> {
424        self.0.iter()
425    }
426
427    pub(crate) fn is_empty(&self) -> bool {
428        self.0.is_empty()
429    }
430
431    pub(crate) fn apply<B>(&self, basis: &[B]) -> B
432    where
433        B: RefZero + Clone + Neg<Output = B> + AddAssign<B>,
434    {
435        // self.panic_validate_basis(basis);
436        let mut result = basis[0].ref_zero();
437        for (&sign, t) in self.0.iter().zip(basis.iter().cloned()) {
438            result += sign * t;
439        }
440        result
441    }
442
443    pub(crate) fn try_apply<B>(&self, basis: &[B]) -> Option<B>
444    where
445        B: Clone + Neg<Output = B> + Add<B, Output = B>,
446    {
447        self.0
448            .iter()
449            .zip(basis.iter().cloned())
450            .filter_map(|(sign, t)| match sign {
451                SignOrZero::Zero => None,
452                SignOrZero::Plus => Some(t),
453                SignOrZero::Minus => Some(-t),
454            })
455            .reduce(|sum, t| sum + t)
456    }
457
458    pub(crate) fn apply_typed<O, I, V>(&self, basis: &V) -> O
459    where
460        V: Index<I, Output = O>,
461        O: RefZero + Neg<Output = O> + AddAssign<O> + Clone,
462        I: From<usize>,
463        usize: From<I>,
464    {
465        let mut result = basis[I::from(0)].ref_zero();
466        for (&sign, i) in self.0.iter().zip(0..) {
467            result += sign * basis[I::from(i)].clone();
468        }
469
470        result
471    }
472
473    pub(crate) fn apply_iter<I, O>(&self, basis: I) -> Option<O>
474    where
475        I: IntoIterator,
476        I::Item: RefZero<O>,
477        O: Clone + SubAssign<I::Item> + AddAssign<I::Item>,
478    {
479        let mut basis_iter = basis.into_iter();
480        let mut signature_iter = self.into_iter();
481
482        while let (Some(sign), Some(item)) = (signature_iter.next(), basis_iter.next()) {
483            if sign.is_sign() {
484                // Initialize the result based on the first non-zero sign
485                let mut result = item.ref_zero();
486                match sign {
487                    SignOrZero::Zero => {
488                        panic!("unreachable");
489                        // return None;
490                    }
491                    SignOrZero::Plus => {
492                        result += item;
493                    }
494                    SignOrZero::Minus => {
495                        result -= item;
496                    }
497                }
498
499                // Continue processing the rest of the iterator
500                while let (Some(sign), Some(item)) = (signature_iter.next(), basis_iter.next()) {
501                    match sign {
502                        SignOrZero::Zero => {}
503                        SignOrZero::Plus => {
504                            result += item;
505                        }
506                        SignOrZero::Minus => {
507                            result -= item;
508                        }
509                    }
510                }
511
512                return Some(result);
513            }
514        }
515
516        // Return None if no non-zero sign was found
517        None
518    }
519
520    pub(crate) fn label_with(&self, label: &str) -> String {
521        let mut result = String::new();
522        let mut first = true;
523        for (i, sign) in self.0.iter().enumerate() {
524            if !first {
525                result.push_str(&sign.to_string());
526            } else {
527                first = false;
528            }
529            if sign.is_sign() {
530                result.push_str(&format!("{}_{}", label, i));
531            }
532        }
533        result
534    }
535
536    /// Canonization function to compare two signatures up to an overall sign,
537    /// If the first nonzero entry is positive, it will return itself,
538    /// otherwise it will return the negative of itself.
539    pub(crate) fn first_abs(&self) -> Self {
540        let sign = self.iter().find(|x| x.is_sign());
541
542        if let Some(sign) = sign {
543            if sign.is_positive() {
544                self.clone()
545            } else {
546                self.iter().map(|x| -*x).collect()
547            }
548        } else {
549            self.clone()
550        }
551    }
552
553    pub(crate) fn pop(&mut self) -> Option<SignOrZero> {
554        self.0.pop()
555    }
556}
557
558#[test]
559fn test_signature() {
560    use crate::momentum::FourMomentum;
561    let sig = LoopSignature::from_iter(vec![SignOrZero::Plus, SignOrZero::Minus]);
562    let basis: Vec<i32> = vec![1, 2];
563    assert_eq!(sig.apply(&basis), 1 - 2);
564    assert_eq!(sig.apply_iter(basis.iter()), Some(-1));
565
566    let basis: [FourMomentum<i32>; 4] = [
567        FourMomentum::from_args(1, 1, 0, 0),
568        FourMomentum::from_args(1, 0, 1, 0),
569        FourMomentum::from_args(1, 0, 0, 1),
570        FourMomentum::from_args(1, 1, 1, 1),
571    ];
572
573    let sig = LoopSignature::from_iter(vec![
574        SignOrZero::Plus,
575        SignOrZero::Minus,
576        SignOrZero::Zero,
577        SignOrZero::Plus,
578    ]);
579
580    assert_eq!(sig.apply(&basis), FourMomentum::from_args(1, 2, 0, 1));
581    let sig = ExternalSignature::from_iter(vec![
582        SignOrZero::Zero,
583        SignOrZero::Zero,
584        SignOrZero::Zero,
585        SignOrZero::Zero,
586    ]);
587    assert_eq!(sig.apply_iter(basis.iter()), None);
588    let sig = ExternalSignature::from_iter(vec![
589        SignOrZero::Zero,
590        SignOrZero::Zero,
591        SignOrZero::Zero,
592        SignOrZero::Minus,
593    ]);
594    assert_eq!(sig.apply_iter(basis.iter()), Some(-basis[3]));
595}
596
597impl LoopExtSignature {
598    pub(crate) fn compute_momentum_untyped<'a, 'b: 'a, T>(
599        &self,
600        loop_moms: &'a [T],
601        external_moms: &'b [T],
602    ) -> T
603    where
604        T: RefZero + Clone + Neg<Output = T> + AddAssign<T>,
605    {
606        if loop_moms.is_empty() {
607            return self.external.apply(external_moms);
608        }
609        if external_moms.is_empty() {
610            return self.internal.apply(loop_moms);
611        }
612        let mut res = self.internal.apply(loop_moms);
613        res += self.external.apply(external_moms);
614        res
615    }
616
617    pub(crate) fn try_compute_momentum<'a, 'b: 'a, T>(
618        &self,
619        loop_moms: &'a [T],
620        external_moms: &'b [T],
621    ) -> Option<T>
622    where
623        T: Clone + Neg<Output = T> + Add<T, Output = T>,
624    {
625        let loop_part = self.internal.try_apply(loop_moms);
626        let external_part = self.external.try_apply(external_moms);
627
628        match (loop_part, external_part) {
629            (Some(l), Some(e)) => Some(l + e),
630            (Some(l), None) => Some(l),
631            (None, Some(e)) => Some(e),
632            (None, None) => None,
633        }
634    }
635
636    pub(crate) fn compute_momentum<L, E, M>(&self, loop_momenta: &L, external_momenta: &E) -> M
637    where
638        M: RefZero + Clone + Neg<Output = M> + AddAssign<M>,
639        L: Index<LoopIndex, Output = M> + Length,
640        E: Index<ExternalIndex, Output = M> + Length,
641    {
642        if loop_momenta.is_empty() {
643            return self.external.apply_typed(external_momenta);
644        }
645
646        if external_momenta.is_empty() {
647            return self.internal.apply_typed(loop_momenta);
648        }
649
650        let mut res = self.internal.apply_typed(loop_momenta);
651        res += self.external.apply_typed(external_momenta);
652        res
653    }
654
655    pub(crate) fn to_momtrop_format(&self) -> (Vec<isize>, Vec<isize>) {
656        (
657            self.internal.to_momtrop_format(),
658            self.external.to_momtrop_format(),
659        )
660    }
661
662    /// Usefull for debugging
663    pub(crate) fn format_momentum(&self) -> String {
664        let mut res = String::new();
665        let mut first = true;
666
667        for (i, sign) in (&self.internal).into_iter().enumerate() {
668            if !first {
669                res.push_str(&sign.to_string());
670            } else {
671                first = false;
672            }
673            if sign.is_sign() {
674                res.push_str(&format!("k_{}", i));
675            }
676        }
677
678        for (i, sign) in (&self.external).into_iter().enumerate() {
679            if !first {
680                res.push_str(&sign.to_string());
681            } else {
682                first = false;
683            }
684            if sign.is_sign() {
685                res.push_str(&format!("l_{}", i));
686            }
687        }
688        res
689    }
690
691    #[allow(unused)]
692    pub(crate) fn compute_four_momentum_from_three<T: FloatLike>(
693        &self,
694        loop_moms: &LoopMomenta<F<T>>,
695        external_moms: &ExternalFourMomenta<F<T>>,
696    ) -> FourMomentum<F<T>> {
697        let loop_moms = loop_moms
698            .iter()
699            .map(|m| m.clone().into_on_shell_four_momentum(None))
700            .collect::<TiVec<LoopIndex, _>>();
701
702        self.compute_momentum(&loop_moms, external_moms)
703    }
704
705    pub(crate) fn compute_three_momentum_from_four<T: FloatLike>(
706        &self,
707        loop_moms: &LoopMomenta<F<T>>,
708        external_moms: &ExternalFourMomenta<F<T>>,
709    ) -> ThreeMomentum<F<T>> {
710        let external_moms: ExternalThreeMomenta<F<T>> =
711            external_moms.iter().map(|m| m.spatial.clone()).collect();
712        self.compute_momentum(loop_moms, &external_moms)
713    }
714}
715
716#[cfg(test)]
717mod tests {
718    use super::LoopExtSignature;
719
720    #[test]
721    fn affine_signature_equality_uses_one_global_sign() {
722        let a = LoopExtSignature::from((vec![1], vec![1]));
723        let independent_reversal = LoopExtSignature::from((vec![-1], vec![1]));
724        assert!(!a.equality_up_to_sign(&independent_reversal));
725        assert!(a.equality_up_to_sign(&LoopExtSignature::from((vec![-1], vec![-1]))));
726        for (internal, external) in [
727            (vec![0, 1], vec![0, 0]),
728            (vec![0, 0], vec![-1, 1]),
729            (vec![0, 0], vec![0, 0]),
730            (vec![], vec![]),
731        ] {
732            let reversed = LoopExtSignature::from((
733                internal.iter().map(|value| -*value).collect(),
734                external.iter().map(|value| -*value).collect(),
735            ));
736            assert!(LoopExtSignature::from((internal, external)).equality_up_to_sign(&reversed));
737        }
738        assert!(
739            !LoopExtSignature::from((vec![1], vec![]))
740                .equality_up_to_sign(&LoopExtSignature::from((vec![], vec![1])))
741        );
742        assert!(!a.equality_up_to_sign(&LoopExtSignature::from((vec![1, 0], vec![1]))));
743    }
744
745    #[test]
746    fn spatial_routing_equality_certifies_native_external_frames() {
747        use crate::{
748            integrands::process::sampling::maps::SamplingEvaluationError,
749            momentum::{ThreeMomentum, sample::ExternalThreeMomenta},
750            utils::{ArbPrec, F, FloatLike, QuadFloat},
751        };
752        fn check<T: FloatLike>() {
753            let one = F::<T>::default().one();
754            let zero = one.zero();
755            let left = LoopExtSignature::from((vec![1], vec![0, 0]));
756            let right = LoopExtSignature::from((vec![1], vec![1, 1]));
757            let beam = ThreeMomentum::new(one.from_i64(2), zero.clone(), one.from_i64(3));
758            let cm = ExternalThreeMomenta::from_iter([beam.clone(), -beam]);
759            assert!(!left.equality_up_to_sign(&right));
760            assert!(left.spatial_equality_up_to_sign(&right, &cm).unwrap());
761            let reversed = LoopExtSignature::from((vec![-1], vec![-1, -1]));
762            assert!(left.spatial_equality_up_to_sign(&reversed, &cm).unwrap());
763
764            // A common boost and even a tiny nonzero represented displacement
765            // break the fixed-frame identity; neither has a norm tolerance.
766            for displacement in [one.clone(), one.epsilon().square()] {
767                let mut shifted = cm.clone();
768                for momentum in &mut shifted {
769                    // Preserve the tiny original input on a zero-valued axis.
770                    momentum.py += &displacement;
771                }
772                assert!(!left.spatial_equality_up_to_sign(&right, &shifted).unwrap());
773            }
774            let a = LoopExtSignature::from((vec![1], vec![1, 0]));
775            let wrong_sign = LoopExtSignature::from((vec![-1], vec![1, 0]));
776            assert!(!a.spatial_equality_up_to_sign(&wrong_sign, &cm).unwrap());
777            // With matching internal signs, the external difference is +2.
778            let doubled = LoopExtSignature::from((vec![1], vec![-1, 0]));
779            assert!(!a.spatial_equality_up_to_sign(&doubled, &cm).unwrap());
780
781            let big = &one / one.epsilon().square();
782            let medium = big.sqrt();
783            let values = [big.clone(), medium.clone(), one.clone(), -medium, -big];
784            assert_eq!(
785                values.iter().fold(zero.clone(), |sum, value| sum + value),
786                zero
787            );
788            let cancelled = ExternalThreeMomenta::from_iter(
789                values.map(|x| ThreeMomentum::new(x, zero.clone(), zero.clone())),
790            );
791            let empty = LoopExtSignature::from((vec![1], vec![0; 5]));
792            let sum = LoopExtSignature::from((vec![1], vec![1; 5]));
793            assert!(!empty.spatial_equality_up_to_sign(&sum, &cancelled).unwrap());
794
795            // These original binary64 values span more than the fixed 2048-bit
796            // certificate budget, also when represented exactly in Quad/Arb.
797            let huge = F(T::from_f64_exact_binary(2.0_f64.powi(1023)));
798            let tiny = F(T::from_f64_exact_binary(f64::from_bits(1)));
799            let unresolved = ExternalThreeMomenta::from_iter(
800                [huge.clone(), tiny, -huge]
801                    .map(|x| ThreeMomentum::new(x, zero.clone(), zero.clone())),
802            );
803            let empty = LoopExtSignature::from((vec![1], vec![0; 3]));
804            let sum = LoopExtSignature::from((vec![1], vec![1; 3]));
805            let error = empty
806                .spatial_equality_up_to_sign(&sum, &unresolved)
807                .unwrap_err();
808            assert!(matches!(
809                error.downcast_ref::<SamplingEvaluationError>(),
810                Some(SamplingEvaluationError::UncertainGeometry { .. })
811            ));
812
813            let mut nonfinite = cm.clone();
814            nonfinite.iter_mut().next().unwrap().px = F(T::from_f64_exact_binary(f64::INFINITY));
815            // Input validation precedes even the structural fast path.
816            let error = left
817                .spatial_equality_up_to_sign(&left, &nonfinite)
818                .unwrap_err();
819            assert!(matches!(
820                error.downcast_ref::<SamplingEvaluationError>(),
821                Some(SamplingEvaluationError::Unrepresentable { .. })
822            ));
823            assert!(
824                left.spatial_equality_up_to_sign(&right, &ExternalThreeMomenta::<F<T>>::default())
825                    .is_err()
826            );
827        }
828        check::<f64>();
829        check::<QuadFloat>();
830        check::<ArbPrec>();
831    }
832
833    #[test]
834    fn raised_edge_comparison_preserves_complete_affine_routing() {
835        use crate::{
836            dot,
837            graph::{Graph, parse::from_dot::IntoGraph},
838            initialisation::test_initialise,
839        };
840        use linnet::half_edge::involution::EdgeIndex;
841        test_initialise().unwrap();
842        let graph: Graph = dot!(digraph raised_affine_signature {
843            ext [style=invis]
844            edge [num=1 mass=0]
845            node [num=1]
846            ext -> a [id=0]
847            a -> b [id=1]
848            b -> c [id=2]
849            c -> a [id=3]
850            ext -> c [id=4]
851        })
852        .unwrap();
853        let groups = graph.get_raised_edge_groups();
854        assert!(groups.iter().any(|group| group.len() >= 2));
855        for group in groups {
856            for pair in group.windows(2) {
857                assert!(graph.loop_momentum_basis.edges_are_raised(pair[0], pair[1]));
858            }
859        }
860        // Deliberately supplied affine rows isolate the caller contract; they
861        // are not a new physical routing assignment for this graph.
862        let mut basis = graph.loop_momentum_basis.clone();
863        basis.edge_signatures[EdgeIndex(1)] = LoopExtSignature::from((vec![1], vec![1, 0]));
864        basis.edge_signatures[EdgeIndex(2)] = LoopExtSignature::from((vec![-1], vec![1, 0]));
865        assert!(!basis.edges_are_raised(EdgeIndex(1), EdgeIndex(2)));
866        basis.edge_signatures[EdgeIndex(2)] = LoopExtSignature::from((vec![-1], vec![-1, 0]));
867        assert!(basis.edges_are_raised(EdgeIndex(1), EdgeIndex(2)));
868    }
869}
gammalooprs/momentum/signature.rs

gammalooprs/momentum/
signature.rs
