spenso_hep_lib/
lib.rs

1use std::{ops::Neg, sync::LazyLock};
2
3use idenso::{gamma::AGS, representations::initialize};
4
5use spenso::{
6    algebra::complex::Complex,
7    network::{
8        Network,
9        library::{
10            TensorLibraryData,
11            function_lib::{INBUILTS, PanicMissingConcrete, SymbolLib},
12            symbolic::{ExplicitKey, TensorLibrary},
13        },
14        parsing::ShadowedStructure,
15        store::NetworkStore,
16    },
17    structure::{PermutedStructure, TensorStructure, abstract_index::AbstractIndex, slot::AbsInd},
18    tensors::{
19        complex::RealOrComplexTensor,
20        data::{SetTensorData, SparseTensor, StorageTensor},
21        parametric::{MixedTensor, ParamOrConcrete},
22    },
23};
24use symbolica::atom::{Atom, Symbol};
25
26#[allow(clippy::similar_names)]
27pub fn gamma_data_dirac<T, N>(structure: N, one: T, zero: T) -> SparseTensor<Complex<T>, N>
28where
29    T: Clone + Neg<Output = T>,
30    N: TensorStructure,
31{
32    let c1 = Complex::<T>::new(one.clone(), zero.clone());
33    let z = Complex::<T>::new(zero.clone(), zero.clone());
34    let cn1 = Complex::<T>::new(-one.clone(), zero.clone());
35    let ci = Complex::<T>::new(zero.clone(), one.clone());
36    let cni = Complex::<T>::new(zero.clone(), -one.clone());
37    let mut gamma = SparseTensor::empty(structure, z);
38    // ! No check on actual structure, should expext mink,bis,bis
39
40    // dirac gamma matrices
41
42    gamma.set(&[0, 0, 0], c1.clone()).unwrap();
43    gamma.set(&[0, 1, 1], c1.clone()).unwrap();
44    gamma.set(&[0, 2, 2], cn1.clone()).unwrap();
45    gamma.set(&[0, 3, 3], cn1.clone()).unwrap();
46
47    gamma.set(&[1, 0, 3], c1.clone()).unwrap();
48    gamma.set(&[1, 1, 2], c1.clone()).unwrap();
49    gamma.set(&[1, 2, 1], cn1.clone()).unwrap();
50    gamma.set(&[1, 3, 0], cn1.clone()).unwrap();
51
52    gamma.set(&[2, 0, 3], cni.clone()).unwrap();
53    gamma.set(&[2, 1, 2], ci.clone()).unwrap();
54    gamma.set(&[2, 2, 1], ci.clone()).unwrap();
55    gamma.set(&[2, 3, 0], cni.clone()).unwrap();
56
57    gamma.set(&[3, 0, 2], c1.clone()).unwrap();
58    gamma.set(&[3, 1, 3], cn1.clone()).unwrap();
59    gamma.set(&[3, 2, 0], cn1.clone()).unwrap();
60    gamma.set(&[3, 3, 1], c1.clone()).unwrap();
61
62    gamma //.to_dense()
63}
64
65#[allow(clippy::similar_names)]
66pub fn gamma_data_weyl<T, N>(structure: N, one: T, zero: T) -> SparseTensor<Complex<T>, N>
67where
68    T: Neg<Output = T> + Clone,
69    N: TensorStructure,
70{
71    let z = Complex::<T>::new(zero.clone(), zero.clone());
72    let c1 = Complex::<T>::new(one.clone(), zero.clone());
73    let cn1 = Complex::<T>::new(-one.clone(), zero.clone());
74    let ci = Complex::<T>::new(zero.clone(), one.clone());
75    let cni = Complex::<T>::new(zero.clone(), -one.clone());
76    let mut gamma = SparseTensor::empty(structure, z);
77    // ! No check on actual structure, should expext mink,bis,bis
78
79    // dirac gamma matrices
80
81    gamma.set(&[0, 2, 0], c1.clone()).unwrap();
82    gamma.set(&[1, 3, 0], c1.clone()).unwrap();
83    gamma.set(&[2, 0, 0], c1.clone()).unwrap();
84    gamma.set(&[3, 1, 0], c1.clone()).unwrap();
85
86    gamma.set(&[0, 3, 1], c1.clone()).unwrap();
87    gamma.set(&[1, 2, 1], c1.clone()).unwrap();
88    gamma.set(&[2, 1, 1], cn1.clone()).unwrap();
89    gamma.set(&[3, 0, 1], cn1.clone()).unwrap();
90
91    gamma.set(&[0, 3, 2], cni.clone()).unwrap();
92    gamma.set(&[1, 2, 2], ci.clone()).unwrap();
93    gamma.set(&[2, 1, 2], ci.clone()).unwrap();
94    gamma.set(&[3, 0, 2], cni.clone()).unwrap();
95
96    gamma.set(&[0, 2, 3], c1.clone()).unwrap();
97    gamma.set(&[1, 3, 3], cn1.clone()).unwrap();
98    gamma.set(&[2, 0, 3], cn1.clone()).unwrap();
99    gamma.set(&[3, 1, 3], c1.clone()).unwrap();
100
101    gamma //.to_dense()
102}
103
104#[allow(clippy::similar_names)]
105pub fn gamma_transpose_weyl<T, N>(structure: N, one: T, zero: T) -> SparseTensor<Complex<T>, N>
106where
107    T: Neg<Output = T> + Clone,
108    N: TensorStructure,
109{
110    let z = Complex::<T>::new(zero.clone(), zero.clone());
111    let c1 = Complex::<T>::new(one.clone(), zero.clone());
112    let cn1 = Complex::<T>::new(-one.clone(), zero.clone());
113    let ci = Complex::<T>::new(zero.clone(), one.clone());
114    let cni = Complex::<T>::new(zero.clone(), -one.clone());
115    let mut gamma = SparseTensor::empty(structure, z);
116    // ! No check on actual structure, should expext mink,bis,bis
117
118    // dirac gamma matrices
119
120    gamma.set(&[2, 0, 0], c1.clone()).unwrap();
121    gamma.set(&[3, 1, 0], c1.clone()).unwrap();
122    gamma.set(&[0, 2, 0], c1.clone()).unwrap();
123    gamma.set(&[1, 3, 0], c1.clone()).unwrap();
124
125    gamma.set(&[3, 0, 1], c1.clone()).unwrap();
126    gamma.set(&[2, 1, 1], c1.clone()).unwrap();
127    gamma.set(&[1, 2, 1], cn1.clone()).unwrap();
128    gamma.set(&[0, 3, 1], cn1.clone()).unwrap();
129
130    gamma.set(&[3, 0, 2], cni.clone()).unwrap();
131    gamma.set(&[2, 1, 2], ci.clone()).unwrap();
132    gamma.set(&[1, 2, 2], ci.clone()).unwrap();
133    gamma.set(&[0, 3, 2], cni.clone()).unwrap();
134
135    gamma.set(&[2, 0, 3], c1.clone()).unwrap();
136    gamma.set(&[3, 1, 3], cn1.clone()).unwrap();
137    gamma.set(&[0, 2, 3], cn1.clone()).unwrap();
138    gamma.set(&[1, 3, 3], c1.clone()).unwrap();
139
140    gamma //.to_dense()
141}
142
143#[allow(clippy::similar_names)]
144pub fn gamma_conj_data_weyl<T, N>(structure: N, one: T, zero: T) -> SparseTensor<Complex<T>, N>
145where
146    T: Neg<Output = T> + Clone,
147    N: TensorStructure,
148{
149    gamma_data_weyl(structure, one, zero).map_data(|a| {
150        let Complex { re, im } = a;
151        Complex { re, im: -im }
152    })
153}
154
155#[allow(clippy::similar_names)]
156pub fn gamma_adj_data_weyl<T, N>(structure: N, one: T, zero: T) -> SparseTensor<Complex<T>, N>
157where
158    T: Neg<Output = T> + Clone,
159    N: TensorStructure,
160{
161    gamma_transpose_weyl(structure, one, zero).map_data(|a| {
162        let Complex { re, im } = a;
163        Complex { re, im: -im }
164    })
165}
166
167pub fn gamma0_weyl<T, N>(structure: N, one: T, zero: T) -> SparseTensor<Complex<T>, N>
168where
169    T: Clone,
170    N: TensorStructure,
171{
172    let c1 = Complex::<T>::new(one, zero.clone());
173    let z = Complex::<T>::new(zero.clone(), zero.clone());
174    let mut gamma0 = SparseTensor::empty(structure, z);
175    // ! No check on actual structure, should expext bis,bis,lor
176
177    // dirac gamma0 matrices
178
179    gamma0.set(&[0, 2], c1.clone()).unwrap();
180    gamma0.set(&[1, 3], c1.clone()).unwrap();
181    gamma0.set(&[2, 0], c1.clone()).unwrap();
182    gamma0.set(&[3, 1], c1.clone()).unwrap();
183
184    gamma0
185}
186
187pub fn gamma5_dirac_data<T, N>(structure: N, one: T, zero: T) -> SparseTensor<Complex<T>, N>
188where
189    T: Clone,
190    N: TensorStructure,
191{
192    let c1 = Complex::<T>::new(one, zero.clone());
193
194    let z = Complex::<T>::new(zero.clone(), zero.clone());
195    let mut gamma5 = SparseTensor::empty(structure, z);
196
197    gamma5.set(&[0, 2], c1.clone()).unwrap();
198    gamma5.set(&[1, 3], c1.clone()).unwrap();
199    gamma5.set(&[2, 0], c1.clone()).unwrap();
200    gamma5.set(&[3, 1], c1.clone()).unwrap();
201
202    gamma5
203}
204
205pub fn gamma5_weyl_data<T, N>(structure: N, one: T, zero: T) -> SparseTensor<Complex<T>, N>
206where
207    T: Clone + Neg<Output = T>,
208    N: TensorStructure,
209{
210    let z = Complex::<T>::new(zero.clone(), zero.clone());
211    let c1 = Complex::<T>::new(one, zero);
212
213    let mut gamma5 = SparseTensor::empty(structure, z);
214
215    gamma5.set(&[0, 0], -c1.clone()).unwrap();
216    gamma5.set(&[1, 1], -c1.clone()).unwrap();
217    gamma5.set(&[2, 2], c1.clone()).unwrap();
218    gamma5.set(&[3, 3], c1.clone()).unwrap();
219
220    gamma5
221}
222
223#[allow(clippy::similar_names)]
224pub fn proj_m_data_dirac<T, N>(structure: N, half: T, zero: T) -> SparseTensor<Complex<T>, N>
225where
226    T: Clone + Neg<Output = T>,
227    N: TensorStructure,
228{
229    let z = Complex::<T>::new(zero.clone(), zero.clone());
230    // ProjM(1,2) Left chirality projector (( 1−γ5)/ 2 )_s1_s2
231
232    let chalf = Complex::<T>::new(half.clone(), zero.clone());
233    let cnhalf = Complex::<T>::new(-half, zero);
234
235    let mut proj_m = SparseTensor::empty(structure, z);
236
237    proj_m.set(&[0, 0], chalf.clone()).unwrap();
238    proj_m.set(&[1, 1], chalf.clone()).unwrap();
239    proj_m.set(&[2, 2], chalf.clone()).unwrap();
240    proj_m.set(&[3, 3], chalf.clone()).unwrap();
241
242    proj_m.set(&[0, 2], cnhalf.clone()).unwrap();
243    proj_m.set(&[1, 3], cnhalf.clone()).unwrap();
244    proj_m.set(&[2, 0], cnhalf.clone()).unwrap();
245    proj_m.set(&[3, 1], cnhalf.clone()).unwrap();
246
247    proj_m
248}
249
250#[allow(clippy::similar_names)]
251pub fn proj_m_data_weyl<T, N>(structure: N, one: T, zero: T) -> SparseTensor<Complex<T>, N>
252where
253    T: Clone,
254    N: TensorStructure,
255{
256    let z = Complex::<T>::new(zero.clone(), zero.clone());
257    // ProjM(1,2) Left chirality projector (( 1−γ5)/ 2 )_s1_s2
258    let c1 = Complex::<T>::new(one, zero);
259    let mut proj_m = SparseTensor::empty(structure, z);
260
261    proj_m.set(&[0, 0], c1.clone()).unwrap();
262    proj_m.set(&[1, 1], c1.clone()).unwrap();
263
264    proj_m
265}
266
267pub fn proj_p_data_dirac<T, N>(structure: N, half: T, zero: T) -> SparseTensor<Complex<T>, N>
268where
269    T: Clone,
270    N: TensorStructure,
271{
272    let z = Complex::<T>::new(zero.clone(), zero.clone());
273    // ProjP(1,2) Right chirality projector (( 1+γ5)/ 2 )_s1_s2
274    let chalf = Complex::<T>::new(half, zero);
275
276    let mut proj_p = SparseTensor::empty(structure, z);
277
278    proj_p
279        .set(&[0, 0], chalf.clone())
280        .unwrap_or_else(|_| unreachable!());
281    proj_p
282        .set(&[1, 1], chalf.clone())
283        .unwrap_or_else(|_| unreachable!());
284    proj_p
285        .set(&[2, 2], chalf.clone())
286        .unwrap_or_else(|_| unreachable!());
287    proj_p
288        .set(&[3, 3], chalf.clone())
289        .unwrap_or_else(|_| unreachable!());
290
291    proj_p
292        .set(&[0, 2], chalf.clone())
293        .unwrap_or_else(|_| unreachable!());
294    proj_p
295        .set(&[1, 3], chalf.clone())
296        .unwrap_or_else(|_| unreachable!());
297    proj_p
298        .set(&[2, 0], chalf.clone())
299        .unwrap_or_else(|_| unreachable!());
300    proj_p
301        .set(&[3, 1], chalf.clone())
302        .unwrap_or_else(|_| unreachable!());
303
304    proj_p
305}
306
307#[allow(clippy::similar_names)]
308pub fn proj_p_data_weyl<T, N>(structure: N, one: T, zero: T) -> SparseTensor<Complex<T>, N>
309where
310    T: Clone,
311    N: TensorStructure,
312{
313    let z = Complex::<T>::new(zero.clone(), zero.clone());
314    // ProjM(1,2) Left chirality projector (( 1−γ5)/ 2 )_s1_s2
315    let c1 = Complex::<T>::new(one, zero);
316    let mut proj_p = SparseTensor::empty(structure, z);
317
318    proj_p.set(&[2, 2], c1.clone()).unwrap();
319    proj_p.set(&[3, 3], c1.clone()).unwrap();
320
321    proj_p
322}
323
324#[allow(clippy::similar_names)]
325pub fn sigma_data<T, N>(structure: N, one: T, zero: T) -> SparseTensor<Complex<T>, N>
326where
327    T: Clone + Neg<Output = T>,
328    N: TensorStructure,
329{
330    let z = Complex::<T>::new(zero.clone(), zero.clone());
331    let c1 = Complex::<T>::new(one.clone(), zero.clone());
332    let cn1 = Complex::<T>::new(-one.clone(), zero.clone());
333    let ci = Complex::<T>::new(zero.clone(), one.clone());
334    let cni = Complex::<T>::new(zero.clone(), -one.clone());
335
336    let mut sigma = SparseTensor::empty(structure, z);
337    sigma.set(&[0, 2, 0, 1], c1.clone()).unwrap();
338    sigma.set(&[0, 2, 3, 0], c1.clone()).unwrap();
339    sigma.set(&[0, 3, 1, 2], c1.clone()).unwrap();
340    sigma.set(&[1, 0, 2, 2], c1.clone()).unwrap();
341    sigma.set(&[1, 1, 1, 2], c1.clone()).unwrap();
342    sigma.set(&[1, 3, 0, 2], c1.clone()).unwrap();
343    sigma.set(&[2, 2, 1, 0], c1.clone()).unwrap();
344    sigma.set(&[2, 2, 2, 1], c1.clone()).unwrap();
345    sigma.set(&[2, 3, 3, 2], c1.clone()).unwrap();
346    sigma.set(&[3, 0, 0, 2], c1.clone()).unwrap();
347    sigma.set(&[3, 3, 2, 2], c1.clone()).unwrap();
348    sigma.set(&[3, 1, 3, 2], c1.clone()).unwrap();
349    sigma.set(&[0, 1, 3, 0], ci.clone()).unwrap();
350    sigma.set(&[0, 3, 1, 1], ci.clone()).unwrap();
351    sigma.set(&[0, 3, 2, 0], ci.clone()).unwrap();
352    sigma.set(&[1, 0, 3, 3], ci.clone()).unwrap();
353    sigma.set(&[1, 1, 0, 3], ci.clone()).unwrap();
354    sigma.set(&[1, 1, 2, 0], ci.clone()).unwrap();
355    sigma.set(&[2, 1, 1, 0], ci.clone()).unwrap();
356    sigma.set(&[2, 3, 0, 0], ci.clone()).unwrap();
357    sigma.set(&[2, 3, 3, 1], ci.clone()).unwrap();
358    sigma.set(&[3, 0, 1, 3], ci.clone()).unwrap();
359    sigma.set(&[3, 1, 0, 0], ci.clone()).unwrap();
360    sigma.set(&[3, 1, 2, 3], ci.clone()).unwrap();
361    sigma.set(&[0, 0, 3, 2], cn1.clone()).unwrap();
362    sigma.set(&[0, 1, 0, 2], cn1.clone()).unwrap();
363    sigma.set(&[0, 2, 1, 3], cn1.clone()).unwrap();
364    sigma.set(&[1, 2, 0, 3], cn1.clone()).unwrap();
365    sigma.set(&[1, 2, 1, 1], cn1.clone()).unwrap();
366    sigma.set(&[1, 2, 2, 0], cn1.clone()).unwrap();
367    sigma.set(&[2, 0, 1, 2], cn1.clone()).unwrap();
368    sigma.set(&[2, 1, 2, 2], cn1.clone()).unwrap();
369    sigma.set(&[2, 2, 3, 3], cn1.clone()).unwrap();
370    sigma.set(&[3, 2, 0, 0], cn1.clone()).unwrap();
371    sigma.set(&[3, 2, 2, 3], cn1.clone()).unwrap();
372    sigma.set(&[3, 2, 3, 1], cn1.clone()).unwrap();
373    sigma.set(&[0, 0, 2, 3], cni.clone()).unwrap();
374    sigma.set(&[0, 0, 3, 1], cni.clone()).unwrap();
375    sigma.set(&[0, 1, 1, 3], cni.clone()).unwrap();
376    sigma.set(&[1, 0, 2, 1], cni.clone()).unwrap();
377    sigma.set(&[1, 3, 0, 1], cni.clone()).unwrap();
378    sigma.set(&[1, 3, 3, 0], cni.clone()).unwrap();
379    sigma.set(&[2, 0, 0, 3], cni.clone()).unwrap();
380    sigma.set(&[2, 0, 1, 1], cni.clone()).unwrap();
381    sigma.set(&[2, 1, 3, 3], cni.clone()).unwrap();
382    sigma.set(&[3, 0, 0, 1], cni.clone()).unwrap();
383    sigma.set(&[3, 3, 1, 0], cni.clone()).unwrap();
384    sigma.set(&[3, 3, 2, 1], cni.clone()).unwrap();
385
386    sigma
387}
388
389pub fn hep_lib<Aind: AbsInd, T: TensorLibraryData + Clone + Default>(
390    one: T,
391    zero: T,
392) -> TensorLibrary<MixedTensor<T, ExplicitKey<Aind>>, Aind>
393where
394{
395    let mut weyl = TensorLibrary::new();
396    initialize();
397    weyl.update_ids();
398
399    let gamma_key = PermutedStructure::identity(
400        gamma_data_weyl(AGS.gamma_strct::<Aind>(4), one.clone(), zero.clone()).into(),
401    );
402    // println!("permutation{}", gamma_key.rep_permutation);
403    weyl.insert_explicit(gamma_key);
404    let gamma_conj_key = PermutedStructure::identity(
405        gamma_conj_data_weyl(AGS.gamma_conj_strct::<Aind>(4), one.clone(), zero.clone()).into(),
406    );
407    // println!("permutation{}", gamma_key.rep_permutation);
408    weyl.insert_explicit(gamma_conj_key);
409    let gamma_adj_key = PermutedStructure::identity(
410        gamma_adj_data_weyl(AGS.gamma_adj_strct::<Aind>(4), one.clone(), zero.clone()).into(),
411    );
412    // println!("permutation{}", gamma_key.rep_permutation);
413    weyl.insert_explicit(gamma_adj_key);
414    let gamma0_key = PermutedStructure::identity(
415        gamma0_weyl(AGS.gamma0_strct::<Aind>(4), one.clone(), zero.clone()).into(),
416    );
417    // println!("permutation{}", gamma_key.rep_permutation);
418    weyl.insert_explicit(gamma0_key);
419
420    let gamma5_key = PermutedStructure::identity(
421        gamma5_weyl_data(AGS.gamma5_strct::<Aind>(4), one.clone(), zero.clone()).into(),
422    );
423    weyl.insert_explicit(gamma5_key);
424
425    let projm_key = PermutedStructure::identity(
426        proj_m_data_weyl(AGS.projm_strct::<Aind>(4), one.clone(), zero.clone()).into(),
427    );
428    weyl.insert_explicit(projm_key);
429
430    let projp_key = PermutedStructure::identity(
431        proj_p_data_weyl(AGS.projp_strct::<Aind>(4), one.clone(), zero.clone()).into(),
432    );
433    weyl.insert_explicit(projp_key);
434
435    weyl
436}
437
438pub fn hep_lib_atom<Aind: AbsInd, T: TensorLibraryData + Clone + Default>()
439-> TensorLibrary<MixedTensor<T, ExplicitKey<Aind>>, Aind>
440where
441{
442    let mut weyl = TensorLibrary::new();
443    initialize();
444    weyl.update_ids();
445
446    let one = Atom::one();
447    let zero = Atom::Zero;
448
449    let gamma_key = PermutedStructure::identity(ParamOrConcrete::param(
450        gamma_data_weyl(AGS.gamma_strct::<Aind>(4), one.clone(), zero.clone())
451            .map_data(|a| a.re + a.im * Atom::i())
452            .into(),
453    ));
454    // println!("permutation{}", gamma_key.rep_permutation);
455    weyl.insert_explicit(gamma_key);
456    let gamma_conj_key = PermutedStructure::identity(ParamOrConcrete::param(
457        gamma_conj_data_weyl(AGS.gamma_conj_strct::<Aind>(4), one.clone(), zero.clone())
458            .map_data(|a| a.re + a.im * Atom::i())
459            .into(),
460    ));
461    // println!("permutation{}", gamma_key.rep_permutation);
462    weyl.insert_explicit(gamma_conj_key);
463    let gamma_adj_key = PermutedStructure::identity(ParamOrConcrete::param(
464        gamma_adj_data_weyl(AGS.gamma_adj_strct::<Aind>(4), one.clone(), zero.clone())
465            .map_data(|a| a.re + a.im * Atom::i())
466            .into(),
467    ));
468    // println!("permutation{}", gamma_key.rep_permutation);
469    weyl.insert_explicit(gamma_adj_key);
470    let gamma0_key = PermutedStructure::identity(ParamOrConcrete::param(
471        gamma0_weyl(AGS.gamma0_strct::<Aind>(4), one.clone(), zero.clone())
472            .map_data(|a| a.re + a.im * Atom::i())
473            .into(),
474    ));
475    // println!("permutation{}", gamma_key.rep_permutation);
476    weyl.insert_explicit(gamma0_key);
477
478    let gamma5_key = PermutedStructure::identity(ParamOrConcrete::param(
479        gamma5_weyl_data(AGS.gamma5_strct::<Aind>(4), one.clone(), zero.clone())
480            .map_data(|a| a.re + a.im * Atom::i())
481            .into(),
482    ));
483    weyl.insert_explicit(gamma5_key);
484
485    let projm_key = PermutedStructure::identity(ParamOrConcrete::param(
486        proj_m_data_weyl(AGS.projm_strct::<Aind>(4), one.clone(), zero.clone())
487            .map_data(|a| a.re + a.im * Atom::i())
488            .into(),
489    ));
490    weyl.insert_explicit(projm_key);
491
492    let projp_key = PermutedStructure::identity(ParamOrConcrete::param(
493        proj_p_data_weyl(AGS.projp_strct::<Aind>(4), one.clone(), zero.clone())
494            .map_data(|a| a.re + a.im * Atom::i())
495            .into(),
496    ));
497    weyl.insert_explicit(projp_key);
498
499    weyl
500}
501
502pub type HepTensor<Aind> = MixedTensor<f64, ShadowedStructure<Aind>>;
503
504pub type HepNet<Aind> =
505    Network<NetworkStore<HepTensor<Aind>, Atom>, ExplicitKey<Aind>, Symbol, Aind>;
506
507pub static HEP_LIB: LazyLock<
508    TensorLibrary<MixedTensor<f64, ExplicitKey<AbstractIndex>>, AbstractIndex>,
509> = LazyLock::new(|| hep_lib(1., 0.));
510
511pub static FUN_LIB: LazyLock<
512    SymbolLib<RealOrComplexTensor<f64, ShadowedStructure<AbstractIndex>>, PanicMissingConcrete>,
513> = LazyLock::new(|| {
514    let mut lib = PanicMissingConcrete::new_lib();
515    lib.insert(INBUILTS.conj, |a| match a {
516        RealOrComplexTensor::Complex(c) => RealOrComplexTensor::Complex(c.map_data(|x| x.conj())),
517        RealOrComplexTensor::Real(r) => RealOrComplexTensor::Real(r),
518    });
519    lib
520});
521
522#[cfg(test)]
523mod tests {
524
525    use spenso::{
526        network::{
527            Network, SingleSmallestDegree, SmallestDegreeIter, Steps,
528            parsing::{ParseSettings, ShadowedStructure},
529            store::NetworkStore,
530        },
531        structure::{HasStructure, abstract_index::AbstractIndex},
532    };
533    use symbolica::{
534        atom::{Atom, Symbol},
535        parse, parse_lit,
536    };
537
538    use super::*;
539
540    #[test]
541    fn simple_scalar() {
542        initialize();
543        let _a = HEP_LIB.get(&AGS.gamma_strct(4)).unwrap();
544
545        let expr = parse!("gamma(bis(4,l_5),bis(4,l_4),mink(4,l_4))*gamma(bis(4,l_6),bis(4,l_5),mink(4,l_4))*gamma(bis(4,l_4),bis(4,l_6),mink(4,l_5))*p(mink(4,l_5))
546            ",default_namespace="spenso");
547        // let expr = parse!(
548        // "gamma(bis(4,l_4),bis(4,l_6),mink(4,l_5))*p(mink(4,l_5))
549        // ",
550        // "spenso"
551        // );
552        // println!("{}", expr);
553
554        let mut net = Network::<
555            NetworkStore<MixedTensor<f64, ShadowedStructure<AbstractIndex>>, Atom>,
556            _,
557            Symbol,
558        >::try_from_view(expr.as_view(), &*HEP_LIB, &ParseSettings::default())
559        .unwrap();
560
561        println!(
562            "{}",
563            net.dot_display_impl(
564                |a| a.to_string(),
565                |a| Some(format!("{}", a.global_name.unwrap())),
566                |a| a.structure().global_name.unwrap().to_string(),
567                |a| a.to_string()
568            )
569        );
570
571        net.execute::<Steps<1>, SmallestDegreeIter<1>, _, _, _>(&*HEP_LIB, &*FUN_LIB)
572            .unwrap();
573        println!(
574            "{}",
575            net.dot_display_impl(
576                |a| a.to_string(),
577                |a| Some(format!("{}", a.global_name?)),
578                |a| a
579                    .structure()
580                    .global_name
581                    .map(|a| a.to_string())
582                    .unwrap_or("".to_string()),
583                |a| a.to_string()
584            )
585        );
586        net.execute::<Steps<1>, SmallestDegreeIter<2>, _, _, _>(&*HEP_LIB, &*FUN_LIB)
587            .unwrap();
588        println!(
589            "{}",
590            net.dot_display_impl(
591                |a| a.to_string(),
592                |a| Some(format!("{}", a.global_name?)),
593                |a| a
594                    .structure()
595                    .global_name
596                    .map(|a| a.to_string())
597                    .unwrap_or("".to_string()),
598                |a| a.to_string()
599            )
600        );
601
602        println!(
603            "{}",
604            net.dot_display_impl(
605                |a| a.to_string(),
606                |_| None,
607                |a| a.to_string(),
608                |a| a.to_string()
609            )
610        );
611        // if let ExecutionResult::Val(TensorOrScalarOrKey::Tensor { tensor, .. }) =
612        //     net.result().unwrap()
613        // {
614        //     // println!("YaY:{}", (&expr - &tensor.expression).expand());
615        //     // assert_eq!(expr, tensor.expression);
616        // } else {
617        //     panic!("Not tensor")
618        // }
619    }
620
621    #[test]
622    // #[should_panic]
623    fn parse_problem() {
624        initialize();
625        let _a = HEP_LIB.get(&AGS.gamma_strct(4)).unwrap();
626
627        let expr = parse_lit!(
628            (-1 * G
629                ^ 3 * P(0, mink(4, 0))
630                    * P(2, mink(4, 26))
631                    * gamma(bis(4, 3), bis(4, 7), mink(4, 4))
632                    * gamma(bis(4, 7), bis(4, 6), mink(4, 1))
633                    * gamma(bis(4, 6), bis(4, 2), mink(4, 26))
634                    + -1 * G
635                ^ 3 * P(0, mink(4, 26))
636                    * P(1, mink(4, 1))
637                    * gamma(bis(4, 3), bis(4, 7), mink(4, 4))
638                    * gamma(bis(4, 7), bis(4, 6), mink(4, 0))
639                    * gamma(bis(4, 6), bis(4, 2), mink(4, 26))
640                    + -1 * G
641                ^ 3 * P(0, mink(4, 26))
642                    * P(1, mink(4, 5))
643                    * g(mink(4, 0), mink(4, 1))
644                    * gamma(bis(4, 3), bis(4, 7), mink(4, 4))
645                    * gamma(bis(4, 7), bis(4, 6), mink(4, 5))
646                    * gamma(bis(4, 6), bis(4, 2), mink(4, 26))
647                    + -1 * G
648                ^ 3 * P(0, mink(4, 5))
649                    * P(2, mink(4, 26))
650                    * g(mink(4, 0), mink(4, 1))
651                    * gamma(bis(4, 3), bis(4, 7), mink(4, 4))
652                    * gamma(bis(4, 6), bis(4, 2), mink(4, 26))
653                    * gamma(bis(4, 7), bis(4, 6), mink(4, 5))
654                    + -1 * G
655                ^ 3 * P(1, mink(4, 1))
656                    * P(1, mink(4, 26))
657                    * gamma(bis(4, 3), bis(4, 7), mink(4, 4))
658                    * gamma(bis(4, 6), bis(4, 2), mink(4, 26))
659                    * gamma(bis(4, 7), bis(4, 6), mink(4, 0))
660                    + -1 * G
661                ^ 3 * P(1, mink(4, 26))
662                    * P(1, mink(4, 5))
663                    * g(mink(4, 0), mink(4, 1))
664                    * gamma(bis(4, 3), bis(4, 7), mink(4, 4))
665                    * gamma(bis(4, 6), bis(4, 2), mink(4, 26))
666                    * gamma(bis(4, 7), bis(4, 6), mink(4, 5))
667                    + -2 * G
668                ^ 3 * P(0, mink(4, 1))
669                    * P(0, mink(4, 26))
670                    * gamma(bis(4, 3), bis(4, 7), mink(4, 4))
671                    * gamma(bis(4, 6), bis(4, 2), mink(4, 26))
672                    * gamma(bis(4, 7), bis(4, 6), mink(4, 0))
673                    + -2 * G
674                ^ 3 * P(0, mink(4, 1))
675                    * P(1, mink(4, 26))
676                    * gamma(bis(4, 3), bis(4, 7), mink(4, 4))
677                    * gamma(bis(4, 6), bis(4, 2), mink(4, 26))
678                    * gamma(bis(4, 7), bis(4, 6), mink(4, 0))
679                    + -2 * G
680                ^ 3 * P(0, mink(4, 5))
681                    * Q(0, mink(4, 5))
682                    * g(mink(4, 0), mink(4, 1))
683                    * gamma(bis(4, 3), bis(4, 2), mink(4, 4))
684                    + -2 * G
685                ^ 3 * P(1, mink(4, 0))
686                    * P(1, mink(4, 1))
687                    * gamma(bis(4, 3), bis(4, 2), mink(4, 4))
688                    + -2 * G
689                ^ 3 * P(1, mink(4, 0))
690                    * P(2, mink(4, 26))
691                    * gamma(bis(4, 3), bis(4, 7), mink(4, 4))
692                    * gamma(bis(4, 6), bis(4, 2), mink(4, 26))
693                    * gamma(bis(4, 7), bis(4, 6), mink(4, 1))
694                    + -2 * G
695                ^ 3 * P(1, mink(4, 1))
696                    * P(2, mink(4, 0))
697                    * gamma(bis(4, 3), bis(4, 2), mink(4, 4))
698                    + -2 * G
699                ^ 3 * P(1, mink(4, 5))
700                    * P(2, mink(4, 5))
701                    * g(mink(4, 0), mink(4, 1))
702                    * gamma(bis(4, 3), bis(4, 2), mink(4, 4))
703                    + -4 * G
704                ^ 3 * P(0, mink(4, 1))
705                    * P(2, mink(4, 0))
706                    * gamma(bis(4, 3), bis(4, 2), mink(4, 4))
707                    + 2 * G
708                ^ 3 * P(0, mink(4, 0))
709                    * P(0, mink(4, 1))
710                    * gamma(bis(4, 3), bis(4, 2), mink(4, 4))
711                    + 2 * G
712                ^ 3 * P(0, mink(4, 0))
713                    * P(2, mink(4, 1))
714                    * gamma(bis(4, 3), bis(4, 2), mink(4, 4))
715                    + 2 * G
716                ^ 3 * P(0, mink(4, 1))
717                    * P(2, mink(4, 26))
718                    * gamma(bis(4, 3), bis(4, 7), mink(4, 4))
719                    * gamma(bis(4, 6), bis(4, 2), mink(4, 26))
720                    * gamma(bis(4, 7), bis(4, 6), mink(4, 0))
721                    + 2 * G
722                ^ 3 * P(0, mink(4, 26))
723                    * P(1, mink(4, 0))
724                    * gamma(bis(4, 3), bis(4, 7), mink(4, 4))
725                    * gamma(bis(4, 6), bis(4, 2), mink(4, 26))
726                    * gamma(bis(4, 7), bis(4, 6), mink(4, 1))
727                    + 2 * G
728                ^ 3 * P(0, mink(4, 5))
729                    * P(2, mink(4, 5))
730                    * g(mink(4, 0), mink(4, 1))
731                    * gamma(bis(4, 3), bis(4, 2), mink(4, 4))
732                    + 2 * G
733                ^ 3 * P(1, mink(4, 0))
734                    * P(1, mink(4, 26))
735                    * gamma(bis(4, 3), bis(4, 7), mink(4, 4))
736                    * gamma(bis(4, 6), bis(4, 2), mink(4, 26))
737                    * gamma(bis(4, 7), bis(4, 6), mink(4, 1))
738                    + 2 * G
739                ^ 3 * P(1, mink(4, i))
740                ^ 2 * g(mink(4, 0), mink(4, 1)) * gamma(bis(4, 3), bis(4, 2), mink(4, 4)) + G
741                ^ 3 * P(1, mink(4, 1))
742                    * P(2, mink(4, 26))
743                    * gamma(bis(4, 3), bis(4, 7), mink(4, 4))
744                    * gamma(bis(4, 6), bis(4, 2), mink(4, 26))
745                    * gamma(bis(4, 7), bis(4, 6), mink(4, 0))),
746            default_namespace = "spenso"
747        );
748        // println!("{}", expr);
749
750        let mut net = Network::<
751            NetworkStore<MixedTensor<f64, ShadowedStructure<AbstractIndex>>, Atom>,
752            _,
753            Symbol,
754        >::try_from_view(expr.as_view(), &*HEP_LIB, &ParseSettings::default())
755        .unwrap();
756
757        net.merge_ops();
758        println!(
759            "{}",
760            net.dot_display_impl(
761                |a| a.to_string(),
762                |a| Some(format!("{}", a.global_name.unwrap())),
763                |a| a.structure().global_name.unwrap().to_string(),
764                |a| a.to_string()
765            )
766        );
767
768        // net.validate();
769        net.execute::<Steps<1>, SingleSmallestDegree<true>, _, _, _>(&*HEP_LIB, &(*FUN_LIB))
770            .unwrap();
771        // net.validate();
772        // net.execute::<Steps<1>, SmallestDegree, _, _>(&*HEP_LIB);
773        // net.validate();
774        // net.execute::<Steps<1>, SmallestDegree, _, _>(&*HEP_LIB);
775        // net.validate();
776        // net.execute::<Steps<1>, SmallestDegree, _, _>(&*HEP_LIB);
777        // net.validate();
778        // net.execute::<Steps<1>, SmallestDegree, _, _>(&*HEP_LIB);
779        // net.validate();
780        // net.execute::<Steps<1>, SmallestDegree, _, _>(&*HEP_LIB);
781        // net.validate();
782        // net.execute::<Steps<1>, SmallestDegree, _, _>(&*HEP_LIB);
783        // net.validate();
784        // net.execute::<Steps<1>, SmallestDegree, _, _>(&*HEP_LIB);
785        // net.validate();
786        // net.execute::<Steps<1>, SmallestDegree, _, _>(&*HEP_LIB);
787        // net.validate();
788        // net.execute::<Steps<1>, SmallestDegree, _, _>(&*HEP_LIB);
789        // net.validate();
790        // net.execute::<Steps<1>, ContractScalars, _, _>(&*HEP_LIB);
791        // net.execute::<Steps<1>, SmallestDegree, _, _>(&*HEP_LIB);
792        // net.execute::<StepsDebug<1>, SingleSmallestDegree<true>, _, _>(&*HEP_LIB);
793        // net.execute::<Steps<1>, ContractScalars, _, _>(&*HEP_LIB);
794
795        //     .unwrap();
796        // net.execute::<Steps<14>, SingleSmallestDegree<false>, _, _>(&*HEP_LIB)
797        //     .unwrap();
798        // net.execute::<Steps<1>, SingleSmallestDegree<true>, _, _>(&*HEP_LIB)
799        //     .unwrap();
800        // // net.execute::<Sequential, SmallestDegree, _, _>(&*HEP_LIB)
801        //     .unwrap();
802        // println!(
803        //     "{}",
804        //     net.dot_display_impl(|a| a.to_string(), |_| None, |a| a.to_string())
805        // );
806
807        println!(
808            "{}",
809            net.dot_display_impl(
810                |a| a.to_string(),
811                |_| None,
812                |a| a.structure().to_string().replace('\n', "\\n"),
813                |a| a.to_string()
814            )
815        );
816        // if let ExecutionResult::Val(TensorOrScalarOrKey::Tensor { tensor, .. }) =
817        //     net.result().unwrap()
818        // {
819        //     // println!("YaY:{}", (&expr - &tensor.expression).expand());
820        //     // assert_eq!(expr, tensor.expression);
821        // } else {
822        //     panic!("Not tensor")
823        // }
824    }
825
826    #[test]
827    fn transpose_test() {}
828}
spenso_hep_lib/lib.rs

spenso_hep_lib/
lib.rs
