test: spincorrelation

feynamp/form/__init__.py

@@ -77,13 +77,12 @@ def compute_squared(
     )
     debug(f"{rr=}")
 
-    ret = sympy.parse_expr(
-        rr.replace("Mom_", "")
-        .replace(".", "_")
+    ret = form.sympyfy(
+        rr
+        # .replace("Mom_", "")
+        # .replace("ms_s", "s")
+        # .replace("ms_u", "u")
+        # .replace("ms_t", "t")
         .replace("PREFACTOR", "1")
-        .replace("^", "**")
-        .replace("mss", "s")
-        .replace("msu", "u")
-        .replace("mst", "t")
     ) * sympy.parse_expr("*".join([*get_color_average(fds), *get_spin_average(fds)]))
     return ret

feynamp/form/form.py

@@ -6,13 +6,18 @@
 import form
 from pqdm.threads import pqdm
 
+from feynamp.leg import get_leg_momentum
+
 count = 0
 dummy = 0
 # TODO auto generate symbols
 init = """
 Symbols Pi,G,ZERO,Tr,Nc,Cf,CA,MC,ee,realpart,PREFACTOR;
 AutoDeclare Index Mu,Spin,Pol,Propagator;
-AutoDeclare Symbol Mass,fd,mss,mst,msu;
+AutoDeclare Symbol Mass,fd;
+* Mandelstamm
+AutoDeclare Symbol ms;
+* Momentum
 AutoDeclare Vector Mom;
 Tensors colorcorrelation,spincorrelation;
 Index scMuMu,scMuNu;
@@ -153,12 +158,16 @@ def sympyfy(string_expr):
 
     ret = simplify(
         parse_expr(
-            string_expr.replace("Mom_", "")
-            .replace(".", "_")
-            .replace("^", "**")
-            .replace("mss", "s")
-            .replace("msu", "u")
-            .replace("mst", "t")
+            string_expr
+            # .replace("Mom_", "")
+            .replace(".", "_").replace("^", "**")
+            # .replace("ms_s", "s")
+            # .replace("ms_u", "u")
+            # .replace("ms_t", "t")
         )
     )
     return simplify(ret.subs("Nc", "3").subs("Cf", "4/3"))
+
+
+def sympy_to_form_string(sympy_expr):
+    return str(sympy_expr).replace("**", "^")

feynamp/form/momentum.py

@@ -100,11 +100,11 @@ def apply_onshell(string_expr, fd, model):
     return run(init + f"Local TMP = {s};" + get_onshell(fd, model))
 
 
-def get_mandelstamm(fds: List[FeynmanDiagram], model: FeynModel, *args, **kwargs):
+def get_mandelstamm(fds: List[FeynmanDiagram], model: FeynModel, **kwargs):
     if fds[0].get_externals_size() == (2, 2):
-        return get_mandelstamm_2_to_2(fds, model, *args, **kwargs)
+        return get_mandelstamm_2_to_2(fds, model, **kwargs)
     elif fds[0].get_externals_size() == (2, 3):
-        return get_mandelstamm_2_to_3(fds, model, *args, **kwargs)
+        return get_mandelstamm_2_to_3(fds, model, **kwargs)
     else:
         raise ValueError("Only 2 to 2 and 2 to 3 Mandelstamm are supported")
 
@@ -115,6 +115,7 @@ def get_mandelstamm_2_to_2(
     replace_s=False,
     replace_t=False,
     replace_u=False,
+    **dead_kwargs,
 ):
     if isinstance(fds, FeynmanDiagram):
         warning(
@@ -146,25 +147,26 @@ def get_mandelstamm_2_to_2(
         p4 = find_leg_in_model(fd, l4, model)
         mom4 = insert_momentum(l4.momentum.name)
         mass4 = insert_mass(string_to_form(p4.mass.name))
-        r += f"id {mom1}.{mom2} = mss/2-{mass1}^2/2-{mass2}^2/2;\n"
-        r += f"id {mom3}.{mom4} = mss/2-{mass3}^2/2-{mass4}^2/2;\n"
-        r += f"id {mom1}.{mom3} = -mst/2+{mass1}^2/2+{mass3}^2/2;\n"
-        r += f"id {mom4}.{mom2} = -mst/2+{mass4}^2/2+{mass2}^2/2;\n"
-        r += f"id {mom1}.{mom4} = -msu/2+{mass1}^2/2+{mass4}^2/2;\n"
-        r += f"id {mom2}.{mom3} = -msu/2+{mass2}^2/2+{mass3}^2/2;\n"
+        r += f"id {mom1}.{mom2} = ms_s/2-{mass1}^2/2-{mass2}^2/2;\n"
+        r += f"id {mom3}.{mom4} = ms_s/2-{mass3}^2/2-{mass4}^2/2;\n"
+        r += f"id {mom1}.{mom3} = -ms_t/2+{mass1}^2/2+{mass3}^2/2;\n"
+        r += f"id {mom4}.{mom2} = -ms_t/2+{mass4}^2/2+{mass2}^2/2;\n"
+        r += f"id {mom1}.{mom4} = -ms_u/2+{mass1}^2/2+{mass4}^2/2;\n"
+        r += f"id {mom2}.{mom3} = -ms_u/2+{mass2}^2/2+{mass3}^2/2;\n"
         if replace_s:
-            r += f"id mss = -msu-mst+{mass2}^2+{mass3}^2+{mass4}^2+{mass1}^2;\n"
+            r += f"id ms_s = -ms_u-ms_t+{mass2}^2+{mass3}^2+{mass4}^2+{mass1}^2;\n"
         if replace_t:
-            r += f"id mst = -mss-msu+{mass2}^2+{mass3}^2+{mass4}^2+{mass1}^2;\n"
+            r += f"id ms_t = -ms_s-ms_u+{mass2}^2+{mass3}^2+{mass4}^2+{mass1}^2;\n"
         if replace_u:
-            r += f"id msu = -mss-mst+{mass2}^2+{mass3}^2+{mass4}^2+{mass1}^2;\n"
+            r += f"id ms_u = -ms_s-ms_t+{mass2}^2+{mass3}^2+{mass4}^2+{mass1}^2;\n"
     return r
 
 
 def get_mandelstamm_2_to_3(
     fds: List[FeynmanDiagram],
-    model: FeynModel
+    model: FeynModel,
     # , replace_s=False, replace_t=False, replace_u=False
+    **dead_kwargs,
 ):
     if isinstance(fds, FeynmanDiagram):
         warning(
@@ -202,18 +204,18 @@ def get_mandelstamm_2_to_3(
         mass5 = insert_mass(string_to_form(p5.mass.name))
 
         # r += f"id {mom5} = {mom1} + {mom2} - {mom3} - {mom4};\n"
-        r += f"id {mom1}.{mom2} = mss12/2-{mass1}^2/2-{mass2}^2/2;\n"
-        r += f"id {mom3}.{mom4} = mss34/2-{mass3}^2/2-{mass4}^2/2;\n"
-        r += f"id {mom3}.{mom5} = mss35/2-{mass3}^2/2-{mass5}^2/2;\n"
-        r += f"id {mom4}.{mom5} = mss45/2-{mass4}^2/2-{mass5}^2/2;\n"
-
-        r += f"id {mom1}.{mom3} = -mst13/2+{mass1}^2/2+{mass3}^2/2;\n"
-        r += f"id {mom1}.{mom4} = -mst14/2+{mass1}^2/2+{mass4}^2/2;\n"
-        r += f"id {mom1}.{mom5} = -mst15/2+{mass1}^2/2+{mass5}^2/2;\n"
-
-        r += f"id {mom2}.{mom3} = -mst23/2+{mass2}^2/2+{mass3}^2/2;\n"
-        r += f"id {mom2}.{mom4} = -mst24/2+{mass2}^2/2+{mass4}^2/2;\n"
-        r += f"id {mom2}.{mom5} = -mst25/2+{mass2}^2/2+{mass5}^2/2;\n"
+        r += f"id {mom1}.{mom2} = ms_s12/2-{mass1}^2/2-{mass2}^2/2;\n"
+        r += f"id {mom3}.{mom4} = ms_s34/2-{mass3}^2/2-{mass4}^2/2;\n"
+        r += f"id {mom3}.{mom5} = ms_s35/2-{mass3}^2/2-{mass5}^2/2;\n"
+        r += f"id {mom4}.{mom5} = ms_s45/2-{mass4}^2/2-{mass5}^2/2;\n"
+
+        r += f"id {mom1}.{mom3} = -ms_t13/2+{mass1}^2/2+{mass3}^2/2;\n"
+        r += f"id {mom1}.{mom4} = -ms_t14/2+{mass1}^2/2+{mass4}^2/2;\n"
+        r += f"id {mom1}.{mom5} = -ms_t15/2+{mass1}^2/2+{mass5}^2/2;\n"
+
+        r += f"id {mom2}.{mom3} = -ms_t23/2+{mass2}^2/2+{mass3}^2/2;\n"
+        r += f"id {mom2}.{mom4} = -ms_t24/2+{mass2}^2/2+{mass4}^2/2;\n"
+        r += f"id {mom2}.{mom5} = -ms_t25/2+{mass2}^2/2+{mass5}^2/2;\n"
 
     return r
 

feynamp/test/colorcorrelation.py

@@ -1,15 +1,13 @@
-from typing import List
-
-from feynml.leg import Leg
-
 from feynamp import sympy
 from feynamp.leg import color_vector_to_casimir, get_color_vector, get_leg_momentum
 
 
-def assert_colorcorrelation(sympy_expr, fd, legs: List[Leg], model):
+def assert_colorcorrelation(sympy_expr, fds, model):
     """
     sympy_expr: sympy expression of cc/born
     """
+    fd = fds[0]
+    legs = fd.legs
     for j in range(len(legs)):
         vec = get_color_vector(fd, legs[j], model)
         if vec is not None:
@@ -20,15 +18,11 @@ def assert_colorcorrelation(sympy_expr, fd, legs: List[Leg], model):
                 momj = get_leg_momentum(legs[j])
                 # colorcorrelation is symmetric, and we only have sorted vertices
                 sum = sum.replace(
-                    sympy.parse_expr(
-                        f"colorcorrelation({momi},{momj})".replace("Mom_", "")
-                    ),
+                    sympy.parse_expr(f"colorcorrelation({momi},{momj})"),
                     1,
                 )
                 sum = sum.replace(
-                    sympy.parse_expr(
-                        f"colorcorrelation({momj},{momi})".replace("Mom_", "")
-                    ),
+                    sympy.parse_expr(f"colorcorrelation({momj},{momi})"),
                     1,
                 )
             # TODO can be optimized

feynamp/test/spincorrelation.py

@@ -1,10 +1,41 @@
-from typing import List
+from feynamp import sympy
+from feynamp.form.form import sympy_to_form_string, sympyfy
+from feynamp.form.momentum import apply, get_kinematics, get_mandelstamm, get_onshell
+from feynamp.leg import (
+    color_vector_to_casimir,
+    get_color_vector,
+    get_leg_momentum,
+    is_vector,
+)
 
-from feynml.leg import Leg
 
-from feynamp import sympy
-from feynamp.leg import color_vector_to_casimir, get_color_vector, get_leg_momentum
+def assert_spincorrelation(sympy_expr, fds, model):
+    """
+    sympy_expr: sympy expression of sc/born
 
+    Checks https://arxiv.org/pdf/1002.2581 Eq. 2.10
+    """
+    fd = fds[0]
+    legs = fd.legs
+    for j, leg in enumerate(legs):
+        if is_vector(fd, leg, model):
+            mom = get_leg_momentum(leg)
+            sum = sympy_expr
+            sum = sum.replace(
+                sympy.parse_expr(f"spincorrelation({mom},scMuMu,scMuNu)"), 1
+            )
+            for i, leg in enumerate(legs):
+                mom = get_leg_momentum(leg)
+                sum = sum.replace(
+                    sympy.parse_expr(f"spincorrelation({mom},scMuMu,scMuNu)"), 0
+                )
+            sum = sum * sympy.parse_expr("d_(MuMu,MuNu)")
+            st = sympy_to_form_string(sum)
+            fs = ""
+            fs += get_kinematics()
+            fs += get_onshell(fds, model)
+            fs += get_mandelstamm(fds, model)
 
-def assert_spincorrelation(sympy_expr, fd, legs: List[Leg], model):
-    assert 1 == 1
+            rr = apply(st, fs)
+            sr = sympyfy(rr)
+            assert sr == sympy.parse_expr("-1")

tests/form/qcd/2_to_2/colorcorrelation/test_colorcorrelation_ee_aa.py

@@ -37,4 +37,4 @@ def test_colorcorrelation_ee_qq():
     cc = compute_squared(fds, fm, colorcorrelated=True)
     assert (cc / born).simplify().equals(sympy.parse_expr("0"))
 
-    assert_colorcorrelation(cc / born, fds[0], fds[0].legs, fm)
+    assert_colorcorrelation(cc / born, fds, fm)

tests/form/qcd/2_to_2/colorcorrelation/test_colorcorrelation_ee_qq.py

@@ -37,4 +37,4 @@ def test_colorcorrelation_ee_qq():
     cc = compute_squared(fds, fm, colorcorrelated=True)
     assert (cc / born).simplify().equals(sympy.parse_expr("Cf*colorcorrelation(p3,p4)"))
 
-    assert_colorcorrelation(cc / born, fds[0], fds[0].legs, fm)
+    assert_colorcorrelation(cc / born, fds, fm)

tests/form/qcd/2_to_2/colorcorrelation/test_colorcorrelation_qq_qq.py

@@ -41,7 +41,7 @@ def test_colorcorrelation_qq_photon_qq():
         .equals("Cf*(colorcorrelation(p1, p2) + colorcorrelation(p3, p4))")
     )
 
-    assert_colorcorrelation(cc / born, fds[0], fds[0].legs, fm)
+    assert_colorcorrelation(cc / born, fds, fm)
 
 
 def test_colorcorrelation_qq_gluon_qq():
@@ -73,4 +73,4 @@ def test_colorcorrelation_qq_gluon_qq():
         )
     )
 
-    assert_colorcorrelation(cc / born, fds[0], fds[0].legs, fm)
+    assert_colorcorrelation(cc / born, fds, fm)

tests/form/qcd/2_to_2/spincorrelation/test_spincorrelation_ee_aa.py

@@ -0,0 +1,40 @@
+import logging
+
+import equation_database.isbn_9780511628788 as ref
+import sympy
+from feynml.interface.qgraf import style
+from feynmodel.interface.qgraf import feynmodel_to_qgraf
+from feynmodel.interface.ufo import load_ufo_model
+from pyfeyn2.feynmandiagram import FeynML
+from pyqgraf import qgraf
+from xsdata.formats.dataclass.parsers import XmlParser
+
+from feynamp.form import compute_squared
+from feynamp.test.colorcorrelation import assert_colorcorrelation
+from feynamp.test.spincorrelation import assert_spincorrelation
+
+logger = logging.getLogger("feynamp")
+logger.setLevel(logging.DEBUG)
+
+
+def test_colorcorrelation_ee_aa():
+    fm = load_ufo_model("ufo_sm")
+    qfm = feynmodel_to_qgraf(fm, True, False)
+
+    qgraf.install()
+    xml_string = qgraf.run(
+        "e_minus[p1], e_plus[p2]",
+        "gamma[p3], gamma[p4]",
+        loops=0,
+        loop_momentum="l",
+        model=qfm,
+        style=style,
+    )
+
+    parser = XmlParser()
+    fml = parser.from_string(xml_string, FeynML)
+    fds = [fml.diagrams[0]]
+    born = compute_squared(fds, fm)
+    sc = compute_squared(fds, fm, spincorrelated=True)
+
+    assert_spincorrelation(sc / born, fds, fm)

tests/form/qcd/2_to_2/spincorrelation/test_spincorrelation_ee_qq.py

@@ -38,4 +38,4 @@ def test_colorcorrelation_ee_qq():
     cc = compute_squared(fds, fm, spincorrelated=True)
     assert (cc / born).simplify().equals(sympy.parse_expr("0"))
 
-    assert_spincorrelation(cc / born, fds[0], fds[0].legs, fm)
+    assert_spincorrelation(cc / born, fds, fm)

tests/form/qcd/2_to_2/tree/test_tree_qqb_qpqpb.py

@@ -73,9 +73,9 @@ def test_qqb_qpqpb():
             rrr.replace("Mom_", "")
             .replace(".", "_")
             .replace("^", "**")
-            .replace("mss", "s")
-            .replace("msu", "u")
-            .replace("mst", "t")
+            .replace("ms_s", "s")
+            .replace("ms_u", "u")
+            .replace("ms_t", "t")
         )
     )
     ret = (

tests/form/qcd/2_to_2/tree/test_tree_qqb_qqb.py

@@ -64,9 +64,9 @@ def test_form_qqb_qqb():
             rr.replace("Mom_", "")
             .replace(".", "_")
             .replace("^", "**")
-            .replace("mss", "s")
-            .replace("msu", "u")
-            .replace("mst", "t")
+            .replace("ms_s", "s")
+            .replace("ms_u", "u")
+            .replace("ms_t", "t")
         )
     )
     # here we use the tags to set the right relative sign