Source code for auto_contractor.operators

#    Qlattice (https://github.com/jinluchang/qlattice)
#
#    Copyright (C) 2021
#
#    Author: Luchang Jin (ljin.luchang@gmail.com)
#    Author: Masaaki Tomii
#
#    This program is free software; you can redistribute it and/or modify
#    it under the terms of the GNU General Public License as published by
#    the Free Software Foundation; either version 3 of the License, or
#    (at your option) any later version.
#
#    This program is distributed in the hope that it will be useful,
#    but WITHOUT ANY WARRANTY; without even the implied warranty of
#    MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
#    GNU General Public License for more details.
#
#    You should have received a copy of the GNU General Public License along
#    with this program; if not, write to the Free Software Foundation, Inc.,
#    51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA.

"""\n
Operators for pi and K follows Eq.(103,122) of the following reference\n
@article{Christ:2019sah,
    author = "Christ, Norman H. and Kelly, Christopher and Zhang, Daiqian",
    title = "{Lattice simulations with G-parity Boundary Conditions}",
    eprint = "1908.08640",
    archivePrefix = "arXiv",
    primaryClass = "hep-lat",
    doi = "10.1103/PhysRevD.101.014506",
    journal = "Phys. Rev. D",
    volume = "101",
    number = "1",
    pages = "014506",
    year = "2020"
}\n
"""

try:
    from .wick import *
    from .compile import *
except:
    from wick import *
    from compile import *

import sympy

spin_index_counter = 0

color_index_counter = 0

def new_spin_index():
    global spin_index_counter
    spin_index_counter += 1
    return f"a_s_{spin_index_counter}"

def new_color_index():
    global color_index_counter
    color_index_counter += 1
    return f"a_c_{color_index_counter}"

saved_sc_indices = []
# saved_sc_indices = None

def save_sc_indices():
    if saved_sc_indices is None:
        return
    saved_sc_indices.append([spin_index_counter, color_index_counter])

def restore_sc_indices():
    if saved_sc_indices is None:
        return
    global spin_index_counter
    global color_index_counter
    spin_index_counter, color_index_counter = saved_sc_indices.pop()

def rsc_call(x, *args):
    restore_sc_indices()
    return x(*args)

def jump_sc_indices(step=100):
    if saved_sc_indices is None:
        return
    global spin_index_counter
    global color_index_counter
    spin_index_counter += 100
    color_index_counter += 100

###################

def show_dagger(is_dagger):
    if not is_dagger:
        return ""
    else:
        return "^dag"

def show_parity(parity):
    if parity is None:
        return ""
    elif parity == "even":
        return ",e"
    elif parity == "odd":
        return ",o"
    else:
        return parity

###################

[docs] def mk_scalar(f1: str, f2: str, p: str, is_dagger=False): """ q1bar q2 """ s = new_spin_index() c = new_color_index() if not is_dagger: return Qb(f1, p, s, c) * Qv(f2, p, s, c) + f"({f1}bar {f2})({p})" else: return Qb(f2, p, s, c) * Qv(f1, p, s, c) + f"({f2}bar {f1})({p})"
[docs] def mk_scalar5(f1: str, f2: str, p: str, is_dagger=False): """ q1bar g5 q2 """ s1 = new_spin_index() s2 = new_spin_index() c = new_color_index() if not is_dagger: return ( Qb(f1, p, s1, c) * G(5, s1, s2) * Qv(f2, p, s2, c) + f"({f1}bar g5 {f2})({p})" ) else: return ( -Qb(f2, p, s1, c) * G(5, s1, s2) * Qv(f1, p, s2, c) + f"(-{f2}bar g5 {f1})({p})" )
[docs] def mk_vec_mu(f1: str, f2: str, p: str, mu, is_dagger=False): """ q1bar gmu q2 """ s1 = new_spin_index() s2 = new_spin_index() c = new_color_index() if not is_dagger: return ( Qb(f1, p, s1, c) * G(mu, s1, s2) * Qv(f2, p, s2, c) + f"({f1}bar g{mu} {f2})({p})" ) else: if mu in [0, 1, 2, 5]: return ( -Qb(f2, p, s1, c) * G(mu, s1, s2) * Qv(f1, p, s2, c) + f"(-{f2}bar g{mu} {f1})({p})" ) else: assert mu in [ 3, ] return ( Qb(f2, p, s1, c) * G(mu, s1, s2) * Qv(f1, p, s2, c) + f"({f2}bar g{mu} {f1})({p})" )
[docs] def mk_vec5_mu(f1: str, f2: str, p: str, mu, is_dagger=False): """ q1bar gmu g5 q2 """ s1 = new_spin_index() s2 = new_spin_index() s3 = new_spin_index() c = new_color_index() if not is_dagger: return ( Qb(f1, p, s1, c) * G(mu, s1, s2) * G(5, s2, s3) * Qv(f2, p, s3, c) + f"({f1}bar g{mu} g5 {f2})({p})" ) else: if mu in [ 0, 1, 2, ]: return ( -Qb(f2, p, s1, c) * G(mu, s1, s2) * G(5, s2, s3) * Qv(f1, p, s3, c) + f"(-{f2}bar g{mu} {f1})({p})" ) else: assert mu in [ 3, 5, ] return ( Qb(f2, p, s1, c) * G(mu, s1, s2) * G(5, s2, s3) * Qv(f1, p, s3, c) + f"({f2}bar g{mu} {f1})({p})" )
[docs] def mk_meson(f1: str, f2: str, p: str, is_dagger=False): """ i q1bar g5 q2 #dag: i q2bar g5 q1 """ if not is_dagger: return sympy.I * mk_scalar5(f1, f2, p, is_dagger) + f"(i {f1}bar g5 {f2})({p})" else: return -sympy.I * mk_scalar5(f1, f2, p, is_dagger) + f"(i {f2}bar g5 {f1})({p})"
###################
[docs] def mk_pi_0(p: str, is_dagger=False): """ i/sqrt(2) * (ubar g5 u - dbar g5 d) #dag: same """ return ( 1 / sympy.sqrt(2) * (mk_meson("u", "u", p, is_dagger) - mk_meson("d", "d", p, is_dagger)) + f"pi0({p}){show_dagger(is_dagger)}" )
[docs] def mk_pi_p(p: str, is_dagger=False): """ i ubar g5 d #dag: i dbar g5 u """ return mk_meson("u", "d", p, is_dagger) + f"pi+({p}){show_dagger(is_dagger)}"
[docs] def mk_pi_m(p: str, is_dagger=False): """ -i dbar g5 u #dag: -i ubar g5 d """ return -mk_meson("d", "u", p, is_dagger) + f"pi-({p}){show_dagger(is_dagger)}"
[docs] def mk_a0_0(p: str, is_dagger=False): """ 1/sqrt(2) * (ubar u - dbar d) """ return ( 1 / sympy.sqrt(2) * (mk_scalar("u", "u", p, is_dagger) - mk_scalar("d", "d", p, is_dagger)) + f"a0_0({p}){show_dagger(is_dagger)}" )
[docs] def mk_a0_p(p: str, is_dagger=False): """ ubar d """ return mk_scalar("u", "d", p, is_dagger) + f"a0_+({p}){show_dagger(is_dagger)}"
[docs] def mk_a0_m(p: str, is_dagger=False): """ dbar u """ return mk_scalar("d", "u", p, is_dagger) + f"a0_-({p}){show_dagger(is_dagger)}"
[docs] def mk_k_p(p: str, is_dagger=False): """ i ubar g5 s #dag: i sbar g5 u """ return mk_meson("u", "s", p, is_dagger) + f"K+({p}){show_dagger(is_dagger)}"
[docs] def mk_k_m(p: str, is_dagger=False): """ -i sbar g5 u #dag: -i ubar g5 s """ return -mk_meson("s", "u", p, is_dagger) + f"K-({p}){show_dagger(is_dagger)}"
[docs] def mk_k_0(p: str, is_dagger=False): """ i dbar g5 s #dag: i sbar g5 d """ return mk_meson("d", "s", p, is_dagger) + f"K0({p}){show_dagger(is_dagger)}"
[docs] def mk_k_0_bar(p: str, is_dagger=False): """ -i sbar g5 d #dag: -i dbar g5 s """ return -mk_meson("s", "d", p, is_dagger) + f"K0b({p}){show_dagger(is_dagger)}"
def mk_eta_l(p: str, is_dagger=False): """ i/sqrt(2) * (ubar g5 u + dbar g5 d) #dag: same """ return ( 1 / sympy.sqrt(2) * (mk_meson("u", "u", p, is_dagger) + mk_meson("d", "d", p, is_dagger)) + f"eta_l({p}){show_dagger(is_dagger)}" ) def mk_eta_s(p: str, is_dagger=False): """ i sbar g5 s #dag: i sbar g5 s """ return mk_meson("s", "s", p, is_dagger) + f"eta_s({p}){show_dagger(is_dagger)}"
[docs] def mk_sigma(p: str, is_dagger=False): """ 1/sqrt(2) * (ubar u + dbar d) """ s = new_spin_index() c = new_color_index() return ( 1 / sympy.sqrt(2) * (Qb("u", p, s, c) * Qv("u", p, s, c) + Qb("d", p, s, c) * Qv("d", p, s, c)) + f"sigma({p})" )
[docs] def mk_kappa_p(p: str, is_dagger=False): """ ubar s """ return mk_scalar("u", "s", p, is_dagger) + f"kappa+({p}){show_dagger(is_dagger)}"
[docs] def mk_kappa_m(p: str, is_dagger=False): """ sbar u """ return mk_scalar("s", "u", p, is_dagger) + f"kappa-({p}){show_dagger(is_dagger)}"
[docs] def mk_kappa_0(p: str, is_dagger=False): """ dbar s """ return mk_scalar("d", "s", p, is_dagger) + f"kappa0({p}){show_dagger(is_dagger)}"
[docs] def mk_kappa_0_bar(p: str, is_dagger=False): """ sbar u """ return mk_scalar("s", "u", p, is_dagger) + f"kappa0bar({p}){show_dagger(is_dagger)}"
def mk_k_p_star_mu(p: str, mu, is_dagger=False): """ ubar gmu s """ return mk_vec_mu("u", "s", p, mu, is_dagger) def mk_k_m_star_mu(p: str, mu, is_dagger=False): """ sbar gmu u """ return mk_vec_mu("s", "u", p, mu, is_dagger)
[docs] def mk_k_0_star_mu(p: str, mu, is_dagger=False): """ dbar gmu s """ return mk_vec_mu("d", "s", p, mu, is_dagger)
[docs] def mk_k_0_star_bar_mu(p: str, mu, is_dagger=False): """ sbar gmu d """ return mk_vec_mu("s", "d", p, mu, is_dagger)
###################
[docs] def mk_pipi_i22(p1: str, p2: str, is_dagger=False): return ( mk_pi_p(p1, is_dagger) * mk_pi_p(p2, is_dagger) + f"pipi_I22({p1},{p2}){show_dagger(is_dagger)}" )
[docs] def mk_pipi_i21(p1: str, p2: str, is_dagger=False): return ( 1 / sympy.sqrt(2) * ( mk_pi_p(p1, is_dagger) * mk_pi_0(p2, is_dagger) + mk_pi_0(p1, is_dagger) * mk_pi_p(p2, is_dagger) ) + f"pipi_I21({p1},{p2}){show_dagger(is_dagger)}" )
[docs] def mk_pipi_i11(p1: str, p2: str, is_dagger=False): return ( 1 / sympy.sqrt(2) * ( mk_pi_p(p1, is_dagger) * mk_pi_0(p2, is_dagger) - mk_pi_0(p1, is_dagger) * mk_pi_p(p2, is_dagger) ) + f"pipi_I11({p1},{p2}){show_dagger(is_dagger)}" )
[docs] def mk_pipi_i20(p1: str, p2: str, is_dagger=False): return ( 1 / sympy.sqrt(6) * ( 2 * mk_pi_0(p1, is_dagger) * mk_pi_0(p2, is_dagger) + mk_pi_m(p1, is_dagger) * mk_pi_p(p2, is_dagger) + mk_pi_p(p1, is_dagger) * mk_pi_m(p2, is_dagger) ) + f"pipi_I20({p1},{p2}){show_dagger(is_dagger)}" )
[docs] def mk_pipi_i10(p1: str, p2: str, is_dagger=False): return ( 1 / sympy.sqrt(2) * ( mk_pi_p(p1, is_dagger) * mk_pi_m(p2, is_dagger) - mk_pi_m(p1, is_dagger) * mk_pi_p(p2, is_dagger) ) + f"pipi_I10({p1},{p2}){show_dagger(is_dagger)}" )
[docs] def mk_pipi_i0(p1: str, p2: str, is_dagger=False): return ( 1 / sympy.sqrt(3) * ( -mk_pi_0(p1, is_dagger) * mk_pi_0(p2, is_dagger) + mk_pi_m(p1, is_dagger) * mk_pi_p(p2, is_dagger) + mk_pi_p(p1, is_dagger) * mk_pi_m(p2, is_dagger) ) + f"pipi_I0({p1},{p2}){show_dagger(is_dagger)}" )
[docs] def mk_kk_i11(p1: str, p2: str, is_dagger=False, *, is_sym=False): if is_sym: return ( 1 / sympy.sqrt(2) * (mk_kk_i11(p1, p2, is_dagger) + mk_kk_i11(p2, p1, is_dagger)) + f"KK_I11({p1},{p2},sym){show_dagger(is_dagger)}" ) else: return ( mk_k_p(p1, is_dagger) * mk_k_0_bar(p2, is_dagger) + f"KK_I11({p1},{p2}){show_dagger(is_dagger)}" )
[docs] def mk_kk_i10(p1: str, p2: str, is_dagger=False, *, is_sym=False): if is_sym: return ( 1 / sympy.sqrt(2) * (mk_kk_i10(p1, p2, is_dagger) + mk_kk_i10(p2, p1, is_dagger)) + f"KK_I10({p1},{p2},sym){show_dagger(is_dagger)}" ) else: return ( 1 / sympy.sqrt(2) * ( -mk_k_0(p1, is_dagger) * mk_k_0_bar(p2, is_dagger) + mk_k_p(p1, is_dagger) * mk_k_m(p2, is_dagger) ) + f"KK_I10({p1},{p2}){show_dagger(is_dagger)}" )
[docs] def mk_kk_i0(p1: str, p2: str, is_dagger=False, *, is_sym=False): if is_sym: return ( 1 / sympy.sqrt(2) * (mk_kk_i0(p1, p2, is_dagger) + mk_kk_i0(p2, p1, is_dagger)) + f"KK_I0({p1},{p2},sym){show_dagger(is_dagger)}" ) else: return ( 1 / sympy.sqrt(2) * ( mk_k_0(p1, is_dagger) * mk_k_0_bar(p2, is_dagger) + mk_k_p(p1, is_dagger) * mk_k_m(p2, is_dagger) ) + f"KK_I0({p1},{p2}){show_dagger(is_dagger)}" )
[docs] def mk_k0k0bar(p1: str, p2: str, is_dagger=False, *, is_sym=False): if is_sym: return ( 1 / sympy.sqrt(2) * (mk_k0k0bar(p1, p2, is_dagger) + mk_k0k0bar(p2, p1, is_dagger)) + f"K0K0b({p1},{p2},sym){show_dagger(is_dagger)}" ) else: return ( mk_k_0(p1, is_dagger) * mk_k_0_bar(p2, is_dagger) + f"K0K0b({p1},{p2}){show_dagger(is_dagger)}" )
[docs] def mk_k0pi0(p1: str, p2: str, is_dagger=False, *, is_sym=False): if is_sym: return ( 1 / sympy.sqrt(2) * (mk_k0pi0(p1, p2, is_dagger) + mk_k0pi0(p2, p1, is_dagger)) + f"K0pi0({p1},{p2},sym){show_dagger(is_dagger)}" ) else: return ( mk_k_0(p1, is_dagger) * mk_pi_0(p2, is_dagger) + f"K0pi0({p1},{p2}){show_dagger(is_dagger)}" )
[docs] def mk_k0barpi0(p1: str, p2: str, is_dagger=False, *, is_sym=False): if is_sym: return ( 1 / sympy.sqrt(2) * (mk_k0barpi0(p1, p2, is_dagger) + mk_k0barpi0(p2, p1, is_dagger)) + f"K0barpi0({p1},{p2},sym){show_dagger(is_dagger)}" ) else: return ( mk_k_0_bar(p1, is_dagger) * mk_pi_0(p2, is_dagger) + f"K0barpi0({p1},{p2}){show_dagger(is_dagger)}" )
[docs] def mk_kppim(p1: str, p2: str, is_dagger=False, *, is_sym=False): if is_sym: return ( 1 / sympy.sqrt(2) * (mk_kppim(p1, p2, is_dagger) + mk_kppim(p2, p1, is_dagger)) + f"K+pi-({p1},{p2},sym){show_dagger(is_dagger)}" ) else: return ( mk_k_p(p1, is_dagger) * mk_pi_m(p2, is_dagger) + f"K+pi-({p1},{p2}){show_dagger(is_dagger)}" )
[docs] def mk_kmpip(p1: str, p2: str, is_dagger=False, *, is_sym=False): if is_sym: return ( 1 / sympy.sqrt(2) * (mk_kmpip(p1, p2, is_dagger) + mk_kmpip(p2, p1, is_dagger)) + f"K-pi+({p1},{p2},sym){show_dagger(is_dagger)}" ) else: return ( mk_k_m(p1, is_dagger) * mk_pi_p(p2, is_dagger) + f"K-pi+({p1},{p2}){show_dagger(is_dagger)}" )
[docs] def mk_kpi_0_i1half( p1: str, p2: str, is_dagger=False, *, is_sym=False ): # strangeness = -1 if is_sym: return ( 1 / sympy.sqrt(2) * (mk_kpi_0_i1half(p1, p2, is_dagger) + mk_kpi_0_i1half(p2, p1, is_dagger)) + f"Kpi_0_I1half({p1},{p2},sym){show_dagger(is_dagger)}" ) else: return ( simplified( sympy.simplify(1) / sympy.sqrt(3) * mk_k_0(p1, is_dagger) * mk_pi_0(p2, is_dagger) + sympy.sqrt(2) / sympy.sqrt(3) * mk_k_p(p1, is_dagger) * mk_pi_m(p2, is_dagger) ) + f"Kpi_0_I1half({p1},{p2}){show_dagger(is_dagger)}" )
[docs] def mk_kpi_p_i1half( p1: str, p2: str, is_dagger=False, *, is_sym=False ): # strangeness = -1 if is_sym: return ( 1 / sympy.sqrt(2) * (mk_kpi_m_i1half(p1, p2, is_dagger) + mk_kpi_m_i1half(p2, p1, is_dagger)) + f"Kpi_+_I1half({p1},{p2},sym){show_dagger(is_dagger)}" ) else: return ( simplified( sympy.sqrt(2) / sympy.sqrt(3) * mk_k_0(p1, is_dagger) * mk_pi_p(p2, is_dagger) + sympy.simplify(1) / sympy.sqrt(3) * mk_k_p(p1, is_dagger) * mk_pi_0(p2, is_dagger) ) + f"Kpi_+_I1half({p1},{p2}){show_dagger(is_dagger)}" )
[docs] def mk_kpi_m_i3halves( p1: str, p2: str, is_dagger=False, *, is_sym=False ): # strangeness = -1 if is_sym: return ( 1 / sympy.sqrt(2) * ( mk_kpi_m_i3halves(p1, p2, is_dagger) + mk_kpi_m_i3halves(p2, p1, is_dagger) ) + f"Kpi_-_I3halves({p1},{p2},sym){show_dagger(is_dagger)}" ) else: return ( mk_k_0(p1, is_dagger) * mk_pi_m(p2, is_dagger) + f"Kpi_-_I3halves({p1},{p2}){show_dagger(is_dagger)}" )
[docs] def mk_kpi_0_i3halves( p1: str, p2: str, is_dagger=False, *, is_sym=False ): # strangeness = -1 if is_sym: return ( 1 / sympy.sqrt(2) * ( mk_kpi_0_i3halves(p1, p2, is_dagger) + mk_kpi_0_i3halves(p2, p1, is_dagger) ) + f"Kpi_0_I3halves({p1},{p2},sym){show_dagger(is_dagger)}" ) else: return ( simplified( -sympy.sqrt(2) / sympy.sqrt(3) * mk_k_0(p1, is_dagger) * mk_pi_0(p2, is_dagger) + sympy.simplify(1) / sympy.sqrt(3) * mk_k_p(p1, is_dagger) * mk_pi_m(p2, is_dagger) ) + f"Kpi_0_I3halves({p1},{p2}){show_dagger(is_dagger)}" )
[docs] def mk_kpi_p1_i3halves( p1: str, p2: str, is_dagger=False, *, is_sym=False ): # strangeness = -1 if is_sym: return ( 1 / sympy.sqrt(2) * ( mk_kpi_p1_i3halves(p1, p2, is_dagger) + mk_kpi_p1_i3halves(p2, p1, is_dagger) ) + f"Kpi_+_I3halves({p1},{p2},sym){show_dagger(is_dagger)}" ) else: return ( simplified( -sympy.simplify(1) / sympy.sqrt(3) * mk_k_0(p1, is_dagger) * mk_pi_p(p2, is_dagger) + sympy.sqrt(2) / sympy.sqrt(3) * mk_k_p(p1, is_dagger) * mk_pi_0(p2, is_dagger) ) + f"Kpi_+_I3halves({p1},{p2}){show_dagger(is_dagger)}" )
[docs] def mk_kpi_p2_i3halves( p1: str, p2: str, is_dagger=False, *, is_sym=False ): # strangeness = -1 if is_sym: return ( 1 / sympy.sqrt(2) * ( mk_kpi_p2_i3halves(p1, p2, is_dagger) + mk_kpi_p2_i3halves(p2, p1, is_dagger) ) + f"Kpi_++_I3halves({p1},{p2},sym){show_dagger(is_dagger)}" ) else: return ( mk_k_p(p1, is_dagger) * mk_pi_p(p2, is_dagger) + f"Kpi_++_I3halves({p1},{p2}){show_dagger(is_dagger)}" )
###################
[docs] def mk_m(f: str, p: str, is_dagger=False): return mk_scalar(f, f, p, is_dagger) + f"{f}bar{f}({p}){show_dagger(is_dagger)}"
[docs] def mk_j5pi_mu(p: str, mu, is_dagger=False): return ( mk_vec5_mu("d", "u", p, mu, is_dagger) + f"j5pi_mu({p},{mu}){show_dagger(is_dagger)}" )
[docs] def mk_j5k_mu(p: str, mu, is_dagger=False): return ( mk_vec5_mu("s", "u", p, mu, is_dagger) + f"j5k_mu({p},{mu}){show_dagger(is_dagger)}" )
[docs] def mk_j5km_mu(p: str, mu, is_dagger=False): return ( -mk_vec5_mu("u", "s", p, mu, is_dagger) + f"j5km_mu({p},{mu}){show_dagger(is_dagger)}" )
def mk_j5eta_l_mu(p: str, mu, is_dagger=False): return ( 1 / sympy.sqrt(2) * ( mk_vec5_mu("u", "u", p, mu, is_dagger) + mk_vec5_mu("d", "d", p, mu, is_dagger) ) + f"j5eta_l_mu({p},{mu}){show_dagger(is_dagger)}" ) def mk_j5eta_s_mu(p: str, mu, is_dagger=False): return ( mk_vec5_mu("s", "s", p, mu, is_dagger) + f"j5eta_s_mu({p},{mu}){show_dagger(is_dagger)}" )
[docs] def mk_jpi_mu(p: str, mu, is_dagger=False): return ( mk_vec_mu("d", "u", p, mu, is_dagger) + f"jpi_mu({p},{mu}){show_dagger(is_dagger)}" )
[docs] def mk_jk_mu(p: str, mu, is_dagger=False): return ( mk_vec_mu("s", "u", p, mu, is_dagger) + f"jk_mu({p},{mu}){show_dagger(is_dagger)}" )
[docs] def mk_j_mu(p: str, mu, is_dagger=False): return ( sympy.simplify(1) * 2 / 3 * mk_vec_mu("u", "u", p, mu, is_dagger) - sympy.simplify(1) * 1 / 3 * mk_vec_mu("d", "d", p, mu, is_dagger) - sympy.simplify(1) * 1 / 3 * mk_vec_mu("s", "s", p, mu, is_dagger) + f"j_mu({p},{mu}){show_dagger(is_dagger)}" )
def mk_j_prime_mu(p: str, mu, is_dagger=False): return ( sympy.simplify(1) * 2 / 3 * mk_vec_mu("u'", "u'", p, mu, is_dagger) - sympy.simplify(1) * 1 / 3 * mk_vec_mu("d'", "d'", p, mu, is_dagger) - sympy.simplify(1) * 1 / 3 * mk_vec_mu("s'", "s'", p, mu, is_dagger) + f"j'_mu({p},{mu}){show_dagger(is_dagger)}" ) def mk_j_prime2_mu(p: str, mu, is_dagger=False): return ( sympy.simplify(1) * 2 / 3 * mk_vec_mu("u''", "u''", p, mu, is_dagger) - sympy.simplify(1) * 1 / 3 * mk_vec_mu("d''", "d''", p, mu, is_dagger) - sympy.simplify(1) * 1 / 3 * mk_vec_mu("s''", "s''", p, mu, is_dagger) + f"j''_mu({p},{mu}){show_dagger(is_dagger)}" )
[docs] def mk_jl_mu(p: str, mu, is_dagger=False): """ jl = sqrt(2)/6 * (j0 + 3 * j10) if no s quark """ return ( sympy.simplify(1) * 2 / 3 * mk_vec_mu("u", "u", p, mu, is_dagger) - sympy.simplify(1) * 1 / 3 * mk_vec_mu("d", "d", p, mu, is_dagger) + f"jl_mu({p},{mu}){show_dagger(is_dagger)}" )
def mk_js_mu(p: str, mu, is_dagger=False): return ( -sympy.simplify(1) * 1 / 3 * mk_vec_mu("s", "s", p, mu, is_dagger) + f"js_mu({p},{mu}){show_dagger(is_dagger)}" )
[docs] def mk_j0_mu(p: str, mu, is_dagger=False): """ I=0 Gparity - """ return ( sympy.simplify(1) * 1 / sympy.sqrt(2) * ( mk_vec_mu("u", "u", p, mu, is_dagger) + mk_vec_mu("d", "d", p, mu, is_dagger) ) + f"j0_mu({p},{mu}){show_dagger(is_dagger)}" )
[docs] def mk_j10_mu(p: str, mu, is_dagger=False): """ I=1 Gparity + """ return ( sympy.simplify(1) * 1 / sympy.sqrt(2) * ( mk_vec_mu("u", "u", p, mu, is_dagger) - mk_vec_mu("d", "d", p, mu, is_dagger) ) + f"j10_mu({p},{mu}){show_dagger(is_dagger)}" )
[docs] def mk_j11_mu(p: str, mu, is_dagger=False): """ I=1 Gparity + """ return ( mk_vec_mu("u", "d", p, mu, is_dagger) + f"j11_mu({p},{mu}){show_dagger(is_dagger)}" )
[docs] def mk_j1n1_mu(p: str, mu, is_dagger=False): """ I=1 Gparity + """ return ( -mk_vec_mu("d", "u", p, mu, is_dagger) + f"j1n1_mu({p},{mu}){show_dagger(is_dagger)}" )
################### def mk_jw_v_mu(p, mu): return mk_jpi_mu(p, mu) + mk_jk_mu(p, mu) def mk_jw_a_mu(p, mu): return mk_j5pi_mu(p, mu) + mk_j5k_mu(p, mu) def mk_sw5(p): return mk_pi_p(p, is_dagger=True) + mk_k_p(p, is_dagger=True) ###################
[docs] def mk_4qOp_VV( f1: str, f2: str, f3: str, f4: str, p, is_scalar=False, parity=None, is_dagger=False ): if parity == "odd": return 0 if is_scalar: return mk_4qOp_SS(f1, f2, f3, f4, p, is_dagger) s = 0 for mu in range(4): save_sc_indices() for mu in range(4): restore_sc_indices() s = s + mk_vec_mu(f1, f2, p, mu, is_dagger) * mk_vec_mu( f3, f4, p, mu, is_dagger ) s.simplify() jump_sc_indices() return s
[docs] def mk_4qOp_VA( f1: str, f2: str, f3: str, f4: str, p, is_scalar=False, parity=None, is_dagger=False ): if parity == "even": return 0 if is_scalar: return mk_4qOp_SP(f1, f2, f3, f4, p, is_dagger) s = 0 for mu in range(4): save_sc_indices() for mu in range(4): restore_sc_indices() s = s + mk_vec_mu(f1, f2, p, mu, is_dagger) * mk_vec5_mu( f3, f4, p, mu, is_dagger ) s.simplify() jump_sc_indices() return s
[docs] def mk_4qOp_AV( f1: str, f2: str, f3: str, f4: str, p, is_scalar=False, parity=None, is_dagger=False ): if parity == "even": return 0 if is_scalar: return mk_4qOp_PS(f1, f2, f3, f4, p, is_dagger) s = 0 for mu in range(4): save_sc_indices() for mu in range(4): restore_sc_indices() s = s + mk_vec5_mu(f1, f2, p, mu, is_dagger) * mk_vec_mu( f3, f4, p, mu, is_dagger ) s.simplify() jump_sc_indices() return s
[docs] def mk_4qOp_AA( f1: str, f2: str, f3: str, f4: str, p, is_scalar=False, parity=None, is_dagger=False ): if parity == "odd": return 0 if is_scalar: return mk_4qOp_PP(f1, f2, f3, f4, p, is_dagger) s = 0 for mu in range(4): save_sc_indices() for mu in range(4): restore_sc_indices() s = s + mk_vec5_mu(f1, f2, p, mu, is_dagger) * mk_vec5_mu( f3, f4, p, mu, is_dagger ) s.simplify() jump_sc_indices() return s
[docs] def mk_4qOp_SS(f1: str, f2: str, f3: str, f4: str, p, is_dagger=False): return mk_scalar(f1, f2, p, is_dagger) * mk_scalar(f3, f4, p, is_dagger)
[docs] def mk_4qOp_SP(f1: str, f2: str, f3: str, f4: str, p, is_dagger=False): return mk_scalar(f1, f2, p, is_dagger) * mk_scalar5(f3, f4, p, is_dagger)
[docs] def mk_4qOp_PS(f1: str, f2: str, f3: str, f4: str, p, is_dagger=False): return mk_scalar5(f1, f2, p, is_dagger) * mk_scalar(f3, f4, p, is_dagger)
[docs] def mk_4qOp_PP(f1: str, f2: str, f3: str, f4: str, p, is_dagger=False): return mk_scalar5(f1, f2, p, is_dagger) * mk_scalar5(f3, f4, p, is_dagger)
[docs] def mk_4qOp_LL(*args): for mu in range(4): save_sc_indices() expr = simplified( rsc_call(mk_4qOp_VV, *args) - rsc_call(mk_4qOp_VA, *args) - rsc_call(mk_4qOp_AV, *args) + rsc_call(mk_4qOp_AA, *args) ) jump_sc_indices() return expr
[docs] def mk_4qOp_LR(*args): for mu in range(4): save_sc_indices() expr = simplified( rsc_call(mk_4qOp_VV, *args) + rsc_call(mk_4qOp_VA, *args) - rsc_call(mk_4qOp_AV, *args) - rsc_call(mk_4qOp_AA, *args) ) jump_sc_indices() return expr
[docs] def mk_4qOp_RL(*args): for mu in range(4): save_sc_indices() expr = simplified( rsc_call(mk_4qOp_VV, *args) - rsc_call(mk_4qOp_VA, *args) + rsc_call(mk_4qOp_AV, *args) - rsc_call(mk_4qOp_AA, *args) ) jump_sc_indices() return expr
[docs] def mk_4qOp_RR(*args): for mu in range(4): save_sc_indices() expr = simplified( rsc_call(mk_4qOp_VV, *args) + rsc_call(mk_4qOp_VA, *args) + rsc_call(mk_4qOp_AV, *args) + rsc_call(mk_4qOp_AA, *args) ) jump_sc_indices() return expr
[docs] def mk_4qOp_LL_cmix(f1, f2, f3, f4, p, is_scalar=False, parity=None, is_dagger=False): assert not is_scalar return mk_4qOp_LL(f1, f4, f3, f2, p, is_scalar, parity, is_dagger)
[docs] def mk_4qOp_LR_cmix(f1, f2, f3, f4, p, is_scalar=False, parity=None, is_dagger=False): assert not is_scalar return -2 * mk_4qOp_RL(f1, f4, f3, f2, p, True, parity, is_dagger)
[docs] def mk_Qsub(p, parity=None, is_dagger=False): if parity is None: expr = mk_Qsub(p, "even", is_dagger) + mk_Qsub(p, "odd", is_dagger) elif parity == "even": expr = mk_scalar("s", "d", p, is_dagger) elif parity == "odd": expr = -mk_scalar5("s", "d", p, is_dagger) return expr + f"Qsub({p}{show_parity(parity)}{show_dagger(is_dagger)})"
[docs] def mk_Q1(p, parity=None, is_dagger=False): return ( mk_4qOp_LL("s", "d", "u", "u", p, False, parity, is_dagger) + f"Q1({p}{show_parity(parity)}{show_dagger(is_dagger)})" )
[docs] def mk_Q2(p, parity=None, is_dagger=False): return ( mk_4qOp_LL_cmix("s", "d", "u", "u", p, False, parity, is_dagger) + f"Q2({p}{show_parity(parity)}{show_dagger(is_dagger)})" )
[docs] def mk_Q3(p, parity=None, is_dagger=False): for mu in range(3): save_sc_indices() expr = simplified( rsc_call(mk_4qOp_LL, "s", "d", "u", "u", p, False, parity, is_dagger) + rsc_call(mk_4qOp_LL, "s", "d", "d", "d", p, False, parity, is_dagger) + rsc_call(mk_4qOp_LL, "s", "d", "s", "s", p, False, parity, is_dagger) ) jump_sc_indices() return expr + f"Q3({p}{show_parity(parity)}{show_dagger(is_dagger)})"
[docs] def mk_Q4(p, parity=None, is_dagger=False): for mu in range(3): save_sc_indices() expr = simplified( rsc_call(mk_4qOp_LL_cmix, "s", "d", "u", "u", p, False, parity, is_dagger) + rsc_call(mk_4qOp_LL_cmix, "s", "d", "d", "d", p, False, parity, is_dagger) + rsc_call(mk_4qOp_LL_cmix, "s", "d", "s", "s", p, False, parity, is_dagger) ) jump_sc_indices() return expr + f"Q4({p}{show_parity(parity)}{show_dagger(is_dagger)})"
[docs] def mk_Q5(p, parity=None, is_dagger=False): for mu in range(3): save_sc_indices() expr = simplified( rsc_call(mk_4qOp_LR, "s", "d", "u", "u", p, False, parity, is_dagger) + rsc_call(mk_4qOp_LR, "s", "d", "d", "d", p, False, parity, is_dagger) + rsc_call(mk_4qOp_LR, "s", "d", "s", "s", p, False, parity, is_dagger) ) jump_sc_indices() return expr + f"Q5({p}{show_parity(parity)}{show_dagger(is_dagger)})"
[docs] def mk_Q6(p, parity=None, is_dagger=False): for mu in range(3): save_sc_indices() expr = simplified( rsc_call(mk_4qOp_LR_cmix, "s", "d", "u", "u", p, False, parity, is_dagger) + rsc_call(mk_4qOp_LR_cmix, "s", "d", "d", "d", p, False, parity, is_dagger) + rsc_call(mk_4qOp_LR_cmix, "s", "d", "s", "s", p, False, parity, is_dagger) ) jump_sc_indices() return expr + f"Q6({p}{show_parity(parity)}{show_dagger(is_dagger)})"
[docs] def mk_Q7(p, parity=None, is_dagger=False): for mu in range(3): save_sc_indices() expr = simplified( rsc_call(mk_4qOp_LR, "s", "d", "u", "u", p, False, parity, is_dagger) - sympy.simplify(1) / 2 * rsc_call(mk_4qOp_LR, "s", "d", "d", "d", p, False, parity, is_dagger) - sympy.simplify(1) / 2 * rsc_call(mk_4qOp_LR, "s", "d", "s", "s", p, False, parity, is_dagger) ) jump_sc_indices() return expr + f"Q7({p}{show_parity(parity)}{show_dagger(is_dagger)})"
[docs] def mk_Q8(p, parity=None, is_dagger=False): for mu in range(3): save_sc_indices() expr = simplified( rsc_call(mk_4qOp_LR_cmix, "s", "d", "u", "u", p, False, parity, is_dagger) - sympy.simplify(1) / 2 * rsc_call(mk_4qOp_LR_cmix, "s", "d", "d", "d", p, False, parity, is_dagger) - sympy.simplify(1) / 2 * rsc_call(mk_4qOp_LR_cmix, "s", "d", "s", "s", p, False, parity, is_dagger) ) jump_sc_indices() return expr + f"Q8({p}{show_parity(parity)}{show_dagger(is_dagger)})"
[docs] def mk_Q9(p, parity=None, is_dagger=False): for mu in range(3): save_sc_indices() expr = simplified( rsc_call(mk_4qOp_LL, "s", "d", "u", "u", p, False, parity, is_dagger) - sympy.simplify(1) / 2 * rsc_call(mk_4qOp_LL, "s", "d", "d", "d", p, False, parity, is_dagger) - sympy.simplify(1) / 2 * rsc_call(mk_4qOp_LL, "s", "d", "s", "s", p, False, parity, is_dagger) ) jump_sc_indices() return expr + f"Q9({p}{show_parity(parity)}{show_dagger(is_dagger)})"
[docs] def mk_Q10(p, parity=None, is_dagger=False): for mu in range(3): save_sc_indices() expr = simplified( rsc_call(mk_4qOp_LL_cmix, "s", "d", "u", "u", p, False, parity, is_dagger) - sympy.simplify(1) / 2 * rsc_call(mk_4qOp_LL_cmix, "s", "d", "d", "d", p, False, parity, is_dagger) - sympy.simplify(1) / 2 * rsc_call(mk_4qOp_LL_cmix, "s", "d", "s", "s", p, False, parity, is_dagger) ) jump_sc_indices() return expr + f"Q10({p}{show_parity(parity)}{show_dagger(is_dagger)})"
################### # 3-flavor operators in (8,1) representation # mk_Q1_b81 # mk_Q2_b81 # mk_Q3_b81 # mk_Q4_b81 # # subtraction operators # mk_Q0_b81 ( = mk_Qsub ) # # charm-contained operators in (8,1) representation # mk_Q5_b81 # mk_Q6_b81 # mk_Q7_b81 # mk_Q8_b81 # # Qa^{e/o} = Aa^{e/o} Q0^{e/o} + Mai Qi^{e/o} ( i = 1, ... ,4; a = 5, ... ,8 )
[docs] def mk_Q0_b81(p, parity=None, is_dagger=False): return mk_Qsub(p, parity, is_dagger)
[docs] def mk_Q1_b81(p, parity=None, is_dagger=False): for mu in range(4): save_sc_indices() expr = simplified( sympy.simplify(1) / sympy.sqrt(10) * rsc_call(mk_4qOp_LL, "s", "d", "u", "u", p, False, parity, is_dagger) + sympy.simplify(1) / sympy.sqrt(10) * rsc_call(mk_4qOp_LL_cmix, "s", "d", "u", "u", p, False, parity, is_dagger) + sympy.simplify(2) / sympy.sqrt(10) * rsc_call(mk_4qOp_LL, "s", "d", "d", "d", p, False, parity, is_dagger) + sympy.simplify(2) / sympy.sqrt(10) * rsc_call(mk_4qOp_LL, "s", "d", "s", "s", p, False, parity, is_dagger) ) jump_sc_indices() return expr + f"Q1_b81({p}{show_parity(parity)}{show_dagger(is_dagger)})"
[docs] def mk_Q2_b81(p, parity=None, is_dagger=False): for mu in range(2): save_sc_indices() expr = simplified( sympy.simplify(1) / sympy.sqrt(2) * rsc_call(mk_4qOp_LL, "s", "d", "u", "u", p, False, parity, is_dagger) - sympy.simplify(1) / sympy.sqrt(2) * rsc_call(mk_4qOp_LL_cmix, "s", "d", "u", "u", p, False, parity, is_dagger) ) jump_sc_indices() return expr + f"Q2_b81({p}{show_parity(parity)}{show_dagger(is_dagger)})"
[docs] def mk_Q3_b81(p, parity=None, is_dagger=False): for mu in range(4): save_sc_indices() expr = simplified( sympy.simplify(1) / sympy.sqrt(3) * rsc_call(mk_4qOp_LR, "s", "d", "u", "u", p, False, parity, is_dagger) + sympy.simplify(1) / sympy.sqrt(3) * rsc_call(mk_4qOp_LR, "s", "d", "d", "d", p, False, parity, is_dagger) + sympy.simplify(1) / sympy.sqrt(3) * rsc_call(mk_4qOp_LR, "s", "d", "s", "s", p, False, parity, is_dagger) ) jump_sc_indices() return expr + f"Q3_b81({p}{show_parity(parity)}{show_dagger(is_dagger)})"
[docs] def mk_Q4_b81(p, parity=None, is_dagger=False): for mu in range(4): save_sc_indices() expr = simplified( sympy.simplify(1) / sympy.sqrt(3) * rsc_call(mk_4qOp_LR_cmix, "s", "d", "u", "u", p, False, parity, is_dagger) + sympy.simplify(1) / sympy.sqrt(3) * rsc_call(mk_4qOp_LR_cmix, "s", "d", "d", "d", p, False, parity, is_dagger) + sympy.simplify(1) / sympy.sqrt(3) * rsc_call(mk_4qOp_LR_cmix, "s", "d", "s", "s", p, False, parity, is_dagger) ) jump_sc_indices() return expr + f"Q4_b81({p}{show_parity(parity)}{show_dagger(is_dagger)})"
[docs] def mk_Q5_b81(p, parity=None, is_dagger=False): return ( mk_4qOp_LL("s", "d", "c", "c", p, False, parity, is_dagger) + f"Q5_b81({p}{show_parity(parity)}{show_dagger(is_dagger)})" )
[docs] def mk_Q6_b81(p, parity=None, is_dagger=False): return ( mk_4qOp_LL_cmix("s", "d", "c", "c", p, False, parity, is_dagger) + f"Q6_b81({p}{show_parity(parity)}{show_dagger(is_dagger)})" )
[docs] def mk_Q7_b81(p, parity=None, is_dagger=False): return ( mk_4qOp_LR("s", "d", "c", "c", p, False, parity, is_dagger) + f"Q7_b81({p}{show_parity(parity)}{show_dagger(is_dagger)})" )
[docs] def mk_Q8_b81(p, parity=None, is_dagger=False): return ( mk_4qOp_LR_cmix("s", "d", "c", "c", p, False, parity, is_dagger) + f"Q8_b81({p}{show_parity(parity)}{show_dagger(is_dagger)})" )
################### def mk_vec_uu_mu(f1, f2, p1, p2, mu, is_dagger=False): """ q1bar gmu Umu q2 """ s1 = new_spin_index() s2 = new_spin_index() c1 = new_color_index() c2 = new_color_index() if not is_dagger: return ( Qb(f1, p1, s1, c1) * G(mu, s1, s2) * U("gauge", p1, mu, c1, c2) * Qv(f2, p2, s2, c2) + f"{f1}bar({p1}) g{mu} U_{mu}({p1}) {f2}({p2})" ) else: if mu in [0, 1, 2, 5]: return ( -Qb(f2, p2, s1, c1) * G(mu, s1, s2) * U("gauge-dagger", p1, mu) * Qv(f1, p1, s2, c2) + f"(-){f2}bar({p2}) g{mu} U^dag_{mu}({p1}) {f1}({p1})" ) else: assert mu in [ 3, ] return ( Qb(f2, p2, s1, c1) * G(mu, s1, s2) * U("gauge-dagger", p1, mu) * Qv(f1, p1, s2, c2) + f"{f2}bar({p2}) g{mu} U^dag_{mu}({p1}) {f1}({p1})" ) ################### def mk_baryon( f1: str, f2: str, f3: str, p: str, spin: str, baryon_type: str = "std", is_dagger=False, ): """ spin in [ "u", "d", ] baryon_type in [ "std", "pos", ] """ assert spin in [ "u", "d", ] assert baryon_type in [ "std", "pos", ] s1 = new_spin_index() s2 = new_spin_index() s3 = new_spin_index() c1 = new_color_index() c2 = new_color_index() c3 = new_color_index() tag = f"{baryon_type}-{spin}" if not is_dagger: v = Bfield(tag, s1, s2, s3, c1, c2, c3) v = v * Qb(f1, p, s1, c1) * Qb(f2, p, s2, c2) * Qb(f3, p, s3, c3) v = v + f"Bb_{baryon_type}({f1},{f2},{f3},{spin})({p})" return v else: v = Bfield(tag, s1, s2, s3, c1, c2, c3) v = v * Qv(f3, p, s3, c3) * Qv(f2, p, s2, c2) * Qv(f1, p, s1, c1) v = v + f"Bv_{baryon_type}({f1},{f2},{f3},{spin})({p})" return v def mk_proton(p: str, spin: str, baryon_type: str = "std", is_dagger=False): """ spin in [ "u", "d", ] baryon_type in [ "std", "pos", ] """ return mk_baryon( "u", "u", "d", p, spin, baryon_type=baryon_type, is_dagger=is_dagger ) def mk_neutron(p: str, spin: str, baryon_type: str = "std", is_dagger=False): """ spin in [ "u", "d", ] baryon_type in [ "std", "pos", ] """ return mk_baryon( "d", "d", "u", p, spin, baryon_type=baryon_type, is_dagger=is_dagger ) def mk_baryon3( f1: str, f2: str, f3: str, p: str, spin: str, baryon_type: str = "std3", is_dagger=False, ): """ spin in [ "u3", "u1", "d1", "d3", ] baryon_type in [ "std3", "pos3", ] """ assert spin in [ "u3", "u1", "d1", "d3", ] assert baryon_type in [ "std3", "pos3", ] s1 = new_spin_index() s2 = new_spin_index() s3 = new_spin_index() c1 = new_color_index() c2 = new_color_index() c3 = new_color_index() tag = f"{baryon_type}-{spin}" if not is_dagger: v = Bfield(tag, s1, s2, s3, c1, c2, c3) v = v * Qb(f1, p, s1, c1) * Qb(f2, p, s2, c2) * Qb(f3, p, s3, c3) v = v + f"Bb_{baryon_type}({f1},{f2},{f3},{spin})({p})" return v else: v = Bfield(tag, s1, s2, s3, c1, c2, c3) v = v * Qv(f3, p, s3, c3) * Qv(f2, p, s2, c2) * Qv(f1, p, s1, c1) v = v + f"Bv_{baryon_type}({f1},{f2},{f3},{spin})({p})" return v def mk_omega(p: str, spin: str, baryon_type: str = "std3", is_dagger=False): """ spin in [ "u3", "u1", "d1", "d3", ] baryon_type in [ "std3", "pos3", ] """ return mk_baryon3( "s", "s", "s", p, spin, baryon_type=baryon_type, is_dagger=is_dagger ) ################### def test(): print("test") args = ["x", None] for mu in range(11): save_sc_indices() expr1 = simplified( rsc_call(mk_Q1, *args) - rsc_call(mk_Q2, *args) - rsc_call(mk_Q3, *args) + rsc_call(mk_Q4, *args) ) expr2 = simplified( 3 * rsc_call(mk_Q1, *args) - rsc_call(mk_Q3, *args) - 2 * rsc_call(mk_Q9, *args) ) expr3 = simplified( rsc_call(mk_Q1, *args) + 2 * rsc_call(mk_Q2, *args) - rsc_call(mk_Q3, *args) - 2 * rsc_call(mk_Q10, *args) ) jump_sc_indices() expr1 = mk_pipi_i0("x1_1", "x1_2", True) * expr1 * mk_k_0("x2") expr2 = mk_pipi_i0("x1_1", "x1_2", True) * expr2 * mk_k_0("x2") expr3 = mk_pipi_i0("x1_1", "x1_2", True) * expr3 * mk_k_0("x2") print(display_cexpr(contract_simplify_compile(expr1, expr2, expr3))) # print(display_cexpr(contract_simplify_compile( # mk_pipi_i0("x1_1", "x1_2", True) * mk_Q1("x") * mk_k_0("x2"), # mk_pipi_i0("x1_1", "x1_2", True) * mk_Q2("x") * mk_k_0("x2"), # mk_pipi_i0("x1_1", "x1_2", True) * mk_Q3("x") * mk_k_0("x2"), # mk_pipi_i0("x1_1", "x1_2", True) * mk_Q4("x") * mk_k_0("x2"), # mk_pipi_i0("x1_1", "x1_2", True) * mk_Q5("x") * mk_k_0("x2"), # mk_pipi_i0("x1_1", "x1_2", True) * mk_Q6("x") * mk_k_0("x2"), # mk_pipi_i0("x1_1", "x1_2", True) * mk_Q7("x") * mk_k_0("x2"), # mk_pipi_i0("x1_1", "x1_2", True) * mk_Q8("x") * mk_k_0("x2"), # mk_pipi_i0("x1_1", "x1_2", True) * mk_Q9("x") * mk_k_0("x2"), # mk_pipi_i0("x1_1", "x1_2", True) * mk_Q10("x") * mk_k_0("x2"), # ))) def test1(): def A(j_p, pi_p, is_dagger=False): # I=21 Gparity + return mk_j10_mu(j_p, 0, is_dagger) * mk_pi_p(pi_p, is_dagger) + mk_j11_mu( j_p, 0, is_dagger ) * mk_pi_0(pi_p, is_dagger) # def B(j_p, pi_p, is_dagger=False): # I=20 Gparity + return ( 2 * mk_j10_mu(j_p, 0, is_dagger) * mk_pi_0(pi_p, is_dagger) + mk_j1n1_mu(j_p, 0, is_dagger) * mk_pi_p(pi_p, is_dagger) + mk_j11_mu(j_p, 0, is_dagger) * mk_pi_m(pi_p, is_dagger) ) # def C(j_p, pi_p, is_dagger=False): # I=11 Gparity + return mk_j10_mu(j_p, 0, is_dagger) * mk_pi_p(pi_p, is_dagger) - mk_j11_mu( j_p, 0, is_dagger ) * mk_pi_0(pi_p, is_dagger) # def D(j_p, pi_p, is_dagger=False): # I=10 Gparity + return mk_j1n1_mu(j_p, 0, is_dagger) * mk_pi_p(pi_p, is_dagger) - mk_j11_mu( j_p, 0, is_dagger ) * mk_pi_m(pi_p, is_dagger) # def E(j_p, pi_p, is_dagger=False): # I=0 Gparity + return ( -mk_j10_mu(j_p, 0, is_dagger) * mk_pi_0(pi_p, is_dagger) + mk_j1n1_mu(j_p, 0, is_dagger) * mk_pi_p(pi_p, is_dagger) + mk_j11_mu(j_p, 0, is_dagger) * mk_pi_m(pi_p, is_dagger) ) # def F(j_p, pi_p, is_dagger=False): # I=11 Gparity - return mk_j0_mu(j_p, 0, is_dagger) * mk_pi_p(pi_p, is_dagger) # def Gf(j_p, pi_p, is_dagger=False): # I=10 Gparity - return mk_j0_mu(j_p, 0, is_dagger) * mk_pi_0(pi_p, is_dagger) # def jpi_p(j_p, pi_p, is_dagger=False): # charged return mk_jl_mu(j_p, 0, is_dagger) * mk_pi_p(pi_p, is_dagger) # def jpi_0(j_p, pi_p, is_dagger=False): # neutral return mk_jl_mu(j_p, 0, is_dagger) * mk_pi_0(pi_p, is_dagger) # expr1 = ( sympy.simplify(1) * 1 / 4 * ( A("xj_2", "x_2", True) * A("xj_1", "x_1") + A("xj_2", "x_2") * A("xj_1", "x_1", True) ) ) expr2 = ( sympy.simplify(1) * 1 / 12 * ( B("xj_2", "x_2", True) * B("xj_1", "x_1") + B("xj_2", "x_2") * B("xj_1", "x_1", True) ) ) expr3 = ( sympy.simplify(1) * 1 / 4 * ( C("xj_2", "x_2", True) * C("xj_1", "x_1") + C("xj_2", "x_2") * C("xj_1", "x_1", True) ) ) expr4 = ( sympy.simplify(1) * 1 / 4 * ( D("xj_2", "x_2", True) * D("xj_1", "x_1") + D("xj_2", "x_2") * D("xj_1", "x_1", True) ) ) expr5 = ( sympy.simplify(1) * 1 / 6 * ( E("xj_2", "x_2", True) * E("xj_1", "x_1") + E("xj_2", "x_2") * E("xj_1", "x_1", True) ) ) expr6 = ( sympy.simplify(1) * 1 / 2 * ( F("xj_2", "x_2", True) * F("xj_1", "x_1") + F("xj_2", "x_2") * F("xj_1", "x_1", True) ) ) expr7 = ( sympy.simplify(1) * 1 / 2 * ( Gf("xj_2", "x_2", True) * Gf("xj_1", "x_1") + Gf("xj_2", "x_2") * Gf("xj_1", "x_1", True) ) ) expr8 = ( sympy.simplify(1) * 1 / (4 * sympy.sqrt(2)) * ( C("xj_2", "x_2", True) * F("xj_1", "x_1") + F("xj_2", "x_2", True) * C("xj_1", "x_1") + C("xj_2", "x_2") * F("xj_1", "x_1", True) + F("xj_2", "x_2") * C("xj_1", "x_1", True) ) ) expr9 = ( sympy.simplify(1) * 1 / (4 * sympy.sqrt(2)) * ( D("xj_2", "x_2", True) * Gf("xj_1", "x_1") + Gf("xj_2", "x_2", True) * D("xj_1", "x_1") + D("xj_2", "x_2") * Gf("xj_1", "x_1", True) + Gf("xj_2", "x_2") * D("xj_1", "x_1", True) ) ) expr_p = ( -sympy.simplify(1) * 1 / 2 * ( jpi_p("xj_2", "x_2", True) * jpi_p("xj_1", "x_1") + jpi_p("xj_2", "x_2") * jpi_p("xj_1", "x_1", True) ) ) expr_0 = ( -sympy.simplify(1) * 1 / 2 * ( jpi_0("xj_2", "x_2", True) * jpi_0("xj_1", "x_1") + jpi_0("xj_2", "x_2") * jpi_0("xj_1", "x_1", True) ) ) exprs = [ expr1, expr1 - expr2, expr3, expr3 - expr4, expr5, expr6, expr6 - expr7, expr8, expr8 - expr9, ] diagram_type_dict = dict() diagram_type_dict[ ( (("x_1", "xj_1"), 1), (("x_2", "xj_2"), 1), (("xj_1", "x_1"), 1), (("xj_2", "x_2"), 1), ) ] = "Type1" diagram_type_dict[ ( (("x_1", "xj_1"), 1), (("x_2", "xj_2"), 1), (("xj_1", "x_2"), 1), (("xj_2", "x_1"), 1), ) ] = "Type2" diagram_type_dict[ ( (("x_1", "x_2"), 1), (("x_2", "xj_1"), 1), (("xj_1", "xj_2"), 1), (("xj_2", "x_1"), 1), ) ] = "Type3" diagram_type_dict[ ( (("x_1", "x_2"), 1), (("x_2", "x_1"), 1), (("xj_1", "xj_2"), 1), (("xj_2", "xj_1"), 1), ) ] = "Type4" print( display_cexpr( contract_simplify_compile(*exprs, diagram_type_dict=diagram_type_dict) ) ) expr_i2_gm = expr1 expr_i1_gm = expr3 expr_i0_gm = expr5 expr_i1_gp = expr6 exprs1 = [ # expr_p + sympy.simplify(1)*1/18*(sympy.simplify(1)*9/2*expr1 + sympy.simplify(1)*9/2*expr3 + expr6 + 3*sympy.sqrt(2)*expr8), expr_p + sympy.simplify(1) * 1 / 18 * ( sympy.simplify(1) * 9 / 2 * expr_i2_gm + sympy.simplify(1) * 9 / 2 * expr_i1_gm + expr_i1_gp ), expr_0 + sympy.simplify(1) * 1 / 18 * (6 * expr_i2_gm + 3 * expr_i0_gm + expr_i1_gp), # (expr_p - expr_0), ] print( display_cexpr( contract_simplify_compile(*exprs1, diagram_type_dict=diagram_type_dict) ) ) etype1n = Term( [ Tr([G(0), S("l", "xj_1", "x_1"), G(5), S("l", "x_1", "xj_1")], "sc"), Tr([G(0), S("l", "xj_2", "x_2"), G(5), S("l", "x_2", "xj_2")], "sc"), ], [], 1, ) etype1r = Term( [ Tr([G(0), S("l", "xj_1", "x_2"), G(5), S("l", "x_2", "xj_1")], "sc"), Tr([G(0), S("l", "xj_2", "x_1"), G(5), S("l", "x_1", "xj_2")], "sc"), ], [], 1, ) etype2 = ( sympy.simplify(1) / 2 * ( Term( [ Tr( [ G(0), S("l", "xj_1", "x_1"), G(5), S("l", "x_1", "xj_2"), G(0), S("l", "xj_2", "x_2"), G(5), S("l", "x_2", "xj_1"), ], "sc", ) ], [], 1, ) + Term( [ Tr( [ G(0), S("l", "xj_1", "x_2"), G(5), S("l", "x_2", "xj_2"), G(0), S("l", "xj_2", "x_1"), G(5), S("l", "x_1", "xj_1"), ], "sc", ) ], [], 1, ) ) ) etype3n = ( sympy.simplify(1) / 2 * ( Term( [ Tr( [ G(0), S("l", "xj_1", "x_1"), G(5), S("l", "x_1", "x_2"), G(5), S("l", "x_2", "xj_2"), G(0), S("l", "xj_2", "xj_1"), ], "sc", ) ], [], 1, ) + Term( [ Tr( [ G(0), S("l", "xj_1", "xj_2"), G(0), S("l", "xj_2", "x_2"), G(5), S("l", "x_2", "x_1"), G(5), S("l", "x_1", "xj_1"), ], "sc", ) ], [], 1, ) ) ) etype3r = ( sympy.simplify(1) / 2 * ( Term( [ Tr( [ G(0), S("l", "xj_1", "x_2"), G(5), S("l", "x_2", "x_1"), G(5), S("l", "x_1", "xj_2"), G(0), S("l", "xj_2", "xj_1"), ], "sc", ) ], [], 1, ) + Term( [ Tr( [ G(0), S("l", "xj_1", "xj_2"), G(0), S("l", "xj_2", "x_1"), G(5), S("l", "x_1", "x_2"), G(5), S("l", "x_2", "xj_1"), ], "sc", ) ], [], 1, ) ) ) etype4 = Term( [ Tr([G(0), S("l", "xj_1", "xj_2"), G(0), S("l", "xj_2", "xj_1")], "sc"), Tr([G(5), S("l", "x_1", "x_2"), G(5), S("l", "x_2", "x_1")], "sc"), ], [], 1, ) exprs2 = [ expr_i2_gm - (etype1r - 2 * etype3r + etype4), # I=2 Gparity + expr_i1_gm - (-etype1r + 2 * etype2 - 2 * etype3n + etype4), # I=1 Gparity + expr_i0_gm - ( 3 * etype1n + etype1r - 3 * etype2 - 3 * etype3n + etype3r + etype4 ), # I=0 Gparity + expr_i1_gp - (-etype2 - etype3n - etype3r + etype4), # I=1 Gparity - ] print( display_cexpr( contract_simplify_compile(*exprs2, diagram_type_dict=diagram_type_dict) ) ) def test_kk(): expr1 = mk_kk_i0("x2_1", "x2_2", True) * mk_kk_i0("x1_1", "x1_2") expr2 = mk_kk_i0("x2_1", "x2_2", True) * mk_pipi_i0("x1_1", "x1_2") expr3 = mk_kk_i0("x2_1", "x2_2", True) * mk_sigma("x1", "x1") expr4 = mk_k0k0bar("x2_1", "x2_2", True) * mk_k0k0bar("x1_1", "x1_2") expr5 = mk_k0k0bar("x2_1", "x2_2", True) * mk_pipi_i0("x1_1", "x1_2") expr6 = mk_k0k0bar("x2_1", "x2_2", True) * mk_sigma("x1", "x1") all_exprs = [ [ expr1, expr4, ], [ expr2, expr5, ], [ expr3, expr6, ], ] names = [ "kk kk", "kk pipi", "kk sigma", ] for name, exprs in zip(names, all_exprs): print(f"\n-- {name} --") cexpr = contract_simplify_compile(*exprs) print(display_cexpr(cexpr)) cexpr.collect_op() print(display_cexpr(cexpr)) print(f"-- {name} --\n") def test_kk_sym(): expr1 = mk_kk_i0("x2_1", "x2_2", True, is_sym=True) * mk_kk_i0( "x1_1", "x1_2", is_sym=True ) expr2 = mk_kk_i0("x2_1", "x2_2", True, is_sym=True) * mk_pipi_i0("x1_1", "x1_2") expr3 = mk_kk_i0("x2_1", "x2_2", True, is_sym=True) * mk_sigma("x1", "x1") expr4 = mk_k0k0bar("x2_1", "x2_2", True, is_sym=True) * mk_k0k0bar( "x1_1", "x1_2", is_sym=True ) expr5 = mk_k0k0bar("x2_1", "x2_2", True, is_sym=True) * mk_pipi_i0("x1_1", "x1_2") expr6 = mk_k0k0bar("x2_1", "x2_2", True, is_sym=True) * mk_sigma("x1", "x1") all_exprs = [ [ expr1, expr4, ], [ expr2, expr5, ], [ expr3, expr6, ], ] names = [ "kk kk", "kk pipi", "kk sigma", ] for name, exprs in zip(names, all_exprs): print(f"\n-- {name} --") cexpr = contract_simplify_compile(*exprs) print(display_cexpr(cexpr)) cexpr.collect_op() print(display_cexpr(cexpr)) print(f"-- {name} --\n") def test_kpipi(): s1 = new_spin_index() s2 = new_spin_index() c = new_color_index() expr1 = Qb("s", "x", s1, c) * G("G1", s1, s2) * Qv("d", "x", s2, c) + "sb_G1_d" s1 = new_spin_index() s2 = new_spin_index() c = new_color_index() expr2 = Qb("u", "x", s1, c) * G("G2", s1, s2) * Qv("u", "x", s2, c) + "ub_G2_u" s1 = new_spin_index() s2 = new_spin_index() c = new_color_index() expr3 = Qb("d", "y", s1, c) * G("5", s1, s2) * Qv("s", "y", s2, c) + "db_g5_s" print(expr1) print(expr2) print(expr3) print( display_cexpr( contract_simplify_compile( expr1 * expr2 * expr3, is_isospin_symmetric_limit=False ) ) ) def test_prop(): s1 = new_spin_index() s2 = new_spin_index() c1 = new_color_index() c2 = new_color_index() expr = Qv("q", "x", s1, c1) * Qb("q", "y", s2, c2) + "q(x) qb(y)" print(expr) print( display_cexpr(contract_simplify_compile(expr, is_isospin_symmetric_limit=False)) ) def test_pipi(): expr1 = ( mk_pi_0("x2_1", True) * mk_pi_0("x2_2", True) * mk_pi_0("x1_1") * mk_pi_0("x1_2") ) expr2 = ( mk_pi_p("x2_1", True) * mk_pi_m("x2_2", True) * mk_pi_p("x1_1") * mk_pi_m("x1_2") ) expr2p = ( mk_pi_m("x2_1", True) * mk_pi_p("x2_2", True) * mk_pi_p("x1_1") * mk_pi_m("x1_2") ) expr = (expr2 + expr2p) - expr1 print(expr) print( display_cexpr(contract_simplify_compile(expr, is_isospin_symmetric_limit=True)) ) def test_meson_corr(): diagram_type_dict = dict() diagram_type_dict[((("t_1", "t_2"), 1), (("t_2", "t_1"), 1))] = "Type1" exprs = [ mk_pi_0("t_1", True) * mk_pi_0("t_2"), mk_k_0("t_1", True) * mk_k_0("t_2"), ] cexpr = contract_simplify_compile( *exprs, is_isospin_symmetric_limit=True, diagram_type_dict=diagram_type_dict ) print() print("meson_corr") print(display_cexpr(cexpr)) def test_meson_f_corr(): diagram_type_dict = dict() diagram_type_dict[((("t_1", "x_2"), 1), (("x_2", "t_1"), 1))] = "Type1" exprs = [ mk_j5pi_mu("x_2", 3) * mk_pi_p("t_1") + "(a_pi * pi)", mk_j5k_mu("x_2", 3) * mk_k_p("t_1") + "(a_k * k )", ] cexpr = contract_simplify_compile( *exprs, is_isospin_symmetric_limit=True, diagram_type_dict=diagram_type_dict ) print() print("meson_f_corr") print(display_cexpr(cexpr)) def test_meson_m(): diagram_type_dict = dict() diagram_type_dict[ ((("t_1", "t_2"), 1), (("t_2", "x_1"), 1), (("x_1", "t_1"), 1)) ] = "Type1" diagram_type_dict[ ((("t_1", "t_2"), 1), (("t_2", "t_1"), 1), (("x_1", "x_1"), 1)) ] = None exprs = [ mk_pi_0("t_1", True) * mk_m("u", "x_1") * mk_pi_0("t_2"), mk_pi_0("t_1", True) * mk_m("d", "x_1") * mk_pi_0("t_2"), mk_pi_p("t_1", True) * mk_m("u", "x_1") * mk_pi_p("t_2"), mk_pi_p("t_1", True) * mk_m("d", "x_1") * mk_pi_p("t_2"), mk_k_0("t_1", True) * mk_m("d", "x_1") * mk_k_0("t_2"), mk_k_0("t_1", True) * mk_m("s", "x_1") * mk_k_0("t_2"), mk_k_p("t_1", True) * mk_m("u", "x_1") * mk_k_p("t_2"), mk_k_p("t_1", True) * mk_m("s", "x_1") * mk_k_p("t_2"), ] cexpr = contract_simplify_compile( *exprs, is_isospin_symmetric_limit=True, diagram_type_dict=diagram_type_dict ) print() print("meson_m") print(display_cexpr(cexpr)) def test_meson_jt(): diagram_type_dict = dict() diagram_type_dict[ ((("t_1", "t_2"), 1), (("t_2", "x_1"), 1), (("x_1", "t_1"), 1)) ] = "Type1" diagram_type_dict[ ((("t_1", "t_2"), 1), (("t_2", "t_1"), 1), (("x_1", "x_1"), 1)) ] = None diagram_type_dict[ ((("t_1", "t_2"), 1), (("t_2", "x_2"), 1), (("x_2", "t_1"), 1)) ] = "Type2" diagram_type_dict[ ((("t_1", "t_2"), 1), (("t_2", "t_1"), 1), (("x_2", "x_2"), 1)) ] = None exprs = [ mk_pi_p("t_1", True) * mk_vec_mu("u", "u", "x_1", 3) * mk_pi_p("t_2"), mk_k_p("t_1", True) * mk_vec_mu("u", "u", "x_1", 3) * mk_k_p("t_2"), mk_k_m("t_1", True) * mk_vec_mu("s", "s", "x_1", 3) * mk_k_m("t_2"), mk_pi_p("t_1p", True) * mk_vec_mu("u", "u", "x_1", 3) * mk_pi_p("t_2p"), mk_k_p("t_1p", True) * mk_vec_mu("u", "u", "x_1", 3) * mk_k_p("t_2p"), mk_k_m("t_1p", True) * mk_vec_mu("s", "s", "x_1", 3) * mk_k_m("t_2p"), ] cexpr = contract_simplify_compile( *exprs, is_isospin_symmetric_limit=True, diagram_type_dict=diagram_type_dict ) print() print("meson_jt") print(display_cexpr(cexpr)) def test_meson_jj(): diagram_type_dict = dict() diagram_type_dict[((("x_1", "x_2"), 1), (("x_2", "x_1"), 1))] = "Type0" diagram_type_dict[((("x_1", "x_1"), 1), (("x_2", "x_2"), 1))] = None diagram_type_dict[ ( (("t_1", "x_1"), 1), (("t_2", "x_2"), 1), (("x_1", "t_1"), 1), (("x_2", "t_2"), 1), ) ] = "Type1" diagram_type_dict[ ( (("t_1", "x_1"), 1), (("t_2", "x_2"), 1), (("x_1", "t_2"), 1), (("x_2", "t_1"), 1), ) ] = "Type2" diagram_type_dict[ ( (("t_1", "t_2"), 1), (("t_2", "x_1"), 1), (("x_1", "x_2"), 1), (("x_2", "t_1"), 1), ) ] = "Type3" diagram_type_dict[ ( (("t_1", "t_2"), 1), (("t_2", "t_1"), 1), (("x_1", "x_2"), 1), (("x_2", "x_1"), 1), ) ] = "Type4" diagram_type_dict[ ( (("t_1", "t_2"), 1), (("t_2", "x_1"), 1), (("x_1", "t_1"), 1), (("x_2", "x_2"), 1), ) ] = None diagram_type_dict[ ( (("t_1", "t_2"), 1), (("t_2", "t_1"), 1), (("x_1", "x_1"), 1), (("x_2", "x_2"), 1), ) ] = None exprs = [ mk_pi_0("t_1", True) * mk_j_mu("x_1", "mu") * mk_j_mu("x_2", "nu") * mk_pi_0("t_2"), mk_pi_p("t_1", True) * mk_j_mu("x_1", "mu") * mk_j_mu("x_2", "nu") * mk_pi_p("t_2"), mk_k_0("t_1", True) * mk_j_mu("x_1", "mu") * mk_j_mu("x_2", "nu") * mk_k_0("t_2"), mk_k_p("t_1", True) * mk_j_mu("x_1", "mu") * mk_j_mu("x_2", "nu") * mk_k_p("t_2"), mk_j_mu("x_1", "mu") * mk_j_mu("x_2", "nu"), ] cexpr = contract_simplify_compile( *exprs, is_isospin_symmetric_limit=True, diagram_type_dict=diagram_type_dict ) print() print("meson_jj") print(display_cexpr(cexpr)) def test_meson_fj(): diagram_type_dict = dict() diagram_type_dict[((("t", "x_1"), 1), (("x_1", "x_2"), 1), (("x_2", "t"), 1))] = ( "Type1" ) diagram_type_dict[((("t", "x_1"), 1), (("x_1", "t"), 1), (("x_2", "x_2"), 1))] = ( None ) exprs = [ mk_j5pi_mu("x_1", "mu") * mk_j_mu("x_2", "nu") * mk_pi_p("t"), mk_j5k_mu("x_1", "mu") * mk_j_mu("x_2", "nu") * mk_k_p("t"), mk_jpi_mu("x_1", "mu") * mk_j_mu("x_2", "nu") * mk_pi_p("t"), mk_jk_mu("x_1", "mu") * mk_j_mu("x_2", "nu") * mk_k_p("t"), ] cexpr = contract_simplify_compile( *exprs, is_isospin_symmetric_limit=True, diagram_type_dict=diagram_type_dict ) print() print("meson_fj") print(display_cexpr(cexpr)) def test_meson_jj_all(): print("< pi+(x_2)^dag j_mu(xj_1) j_nu(xj_2) pi+(x_1) >") print("< pi-(x_2)^dag j_mu(xj_1) j_nu(xj_2) pi-(x_1) >") print("< pi0(x_2)^dag j_mu(xj_1) j_nu(xj_2) pi0(x_1) >") pi_p = ( mk_pi_p("x_2", True) * mk_jl_mu("xj_1", "mu") * mk_jl_mu("xj_2", "nu") * mk_pi_p("x_1") + "< pi+(x_2)^dag j_mu(xj_1) j_nu(xj_2) pi+(x_1) >" ) pi_m = ( mk_pi_m("x_2", True) * mk_jl_mu("xj_1", "mu") * mk_jl_mu("xj_2", "nu") * mk_pi_m("x_1") + "< pi-(x_2)^dag j_mu(xj_1) j_nu(xj_2) pi-(x_1) >" ) pi_0 = ( mk_pi_0("x_1", True) * mk_jl_mu("xj_1", "mu") * mk_jl_mu("xj_2", "nu") * mk_pi_0("x_2") + "< pi0(x_2)^dag j_mu(xj_1) j_nu(xj_2) pi0(x_1) >" ) exprs = [ pi_p, pi_m, mk_sym(1) / 2 * (pi_p + pi_m), pi_0, ] print(display_cexpr(contract_simplify_compile(*exprs))) if __name__ == "__main__": test_meson_jj_all() exit() # test() # test1() # test_pipi() # test_prop() # test_kpipi() # test_kk() # test_kk_sym() test_meson_corr() test_meson_f_corr() test_meson_m() test_meson_jt() test_meson_jj() test_meson_fj() exit() # print("pi+(x1):\n", mk_pi_p("x1")) print("pi+(x2)^dag pi+(x1):\n", mk_pi_p("x2", True) * mk_pi_p("x1")) print("< pi+(x2)^dag pi+(x1) >:") expr = mk_pi_p("x2", True) * mk_pi_p("x1") print(display_cexpr(contract_simplify_compile(expr))) print() print("pi0(x1):\n", mk_pi_0("x1")) print("pi0(x2)^dag pi0(x1):\n", mk_pi_0("x2", True) * mk_pi_0("x1")) print("< pi0(x2)^dag pi0(x1) >:") expr = mk_pi_0("x2", True) * mk_pi_0("x1") print( display_cexpr(contract_simplify_compile(expr, is_isospin_symmetric_limit=False)) ) print(display_cexpr(contract_simplify_compile(expr))) print() print( "< pi0(x2)^dag pi0(x1) > - < pi+(x2)^dag pi+(x1) >: ", simplified( contract_expr( mk_pi_0("x2", True) * mk_pi_0("x1") - mk_pi_p("x2", True) * mk_pi_p("x1") ), is_isospin_symmetric_limit=True, ), ) print() print("< pipiI22(x2_1,x2_2)^dag pipiI22(x1_1,x1_2) >:") expr = mk_pipi_i22("x2_1", "x2_2", True) * mk_pipi_i22("x1_1", "x1_2") print(display_cexpr(contract_simplify_compile(expr))) print( "< pipiI21(x2_1,x2_2)^dag pipiI21(x1_1,x1_2) - pipiI22(x2_1,x2_2)^dag pipiI22(x1_1,x1_2) >:" ) expr = mk_pipi_i21("x2_1", "x2_2", True) * mk_pipi_i21( "x1_1", "x1_2" ) - mk_pipi_i22("x2_1", "x2_2", True) * mk_pipi_i22("x1_1", "x1_2") print(display_cexpr(contract_simplify_compile(expr))) print( "< pipiI20(x2_1,x2_2)^dag pipiI20(x1_1,x1_2) - pipiI22(x2_1,x2_2)^dag pipiI22(x1_1,x1_2) >:" ) expr = mk_pipi_i20("x2_1", "x2_2", True) * mk_pipi_i20( "x1_1", "x1_2" ) - mk_pipi_i22("x2_1", "x2_2", True) * mk_pipi_i22("x1_1", "x1_2") print(display_cexpr(contract_simplify_compile(expr))) print() print("< pipiI11(x2_1,x2_2)^dag pipiI11(x1_1,x1_2) >:") expr = mk_pipi_i11("x2_1", "x2_2", True) * mk_pipi_i11("x1_1", "x1_2") print(display_cexpr(contract_simplify_compile(expr))) print( "< pipiI10(x2_1,x2_2)^dag pipiI10(x1_1,x1_2) - pipiI11(x2_1,x2_2)^dag pipiI11(x1_1,x1_2) >:" ) expr = mk_pipi_i10("x2_1", "x2_2", True) * mk_pipi_i10( "x1_1", "x1_2" ) - mk_pipi_i11("x2_1", "x2_2", True) * mk_pipi_i11("x1_1", "x1_2") print(display_cexpr(contract_simplify_compile(expr))) print() print("< pipiI0(x1_1,x1_2) >:") expr = mk_pipi_i0("x1_1", "x1_2") print(display_cexpr(contract_simplify_compile(expr))) print("< pipiI0(x2_1,x2_2)^dag pipiI0(x1_1,x1_2) >:") expr = mk_pipi_i0("x2_1", "x2_2", True) * mk_pipi_i0("x1_1", "x1_2") print(display_cexpr(contract_simplify_compile(expr))) print() print("< sigma(x1) >:") expr = mk_sigma("x1") print(display_cexpr(contract_simplify_compile(expr))) print("< pipiI0(x2_1,x2_2)^dag sigma(x1) >:") expr = mk_pipi_i0("x2_1", "x2_2", True) * mk_sigma("x1") print(display_cexpr(contract_simplify_compile(expr))) print() print("< K+(x2)^dag K+(x1)>:") expr = mk_k_p("x2", True) * mk_k_p("x1") print(display_cexpr(contract_simplify_compile(expr))) print() print("< pi0(x2)^dag j_mu(xj_1) j_nu(xj_2) pi0(x1) >:") expr = ( mk_pi_0("x1", True) * mk_j_mu("xj_1", "mu") * mk_j_mu("xj_2", "nu") * mk_pi_0("x2") ) print(display_cexpr(contract_simplify_compile(expr))) print() print( "< pi+(x2)^dag j_mu(xj_1) j_nu(xj_2) pi+(x1) / 2 + pi-(x2)^dag j_mu(xj_1) j_nu(xj_2) pi-(x1) / 2 >:" ) expr = sympy.simplify(1) / 2 * mk_pi_p("x2", True) * mk_j_mu( "xj_1", "mu" ) * mk_j_mu("xj_2", "nu") * mk_pi_p("x1") + sympy.simplify(1) / 2 * mk_pi_m( "x2", True ) * mk_j_mu("xj_1", "mu") * mk_j_mu("xj_2", "nu") * mk_pi_m("x1") print(display_cexpr(contract_simplify_compile(expr))) print() print( "< pi+(x2)^dag j_mu(xj_1) j_nu(xj_2) pi+(x1) / 2 + pi-(x2)^dag j_mu(xj_1) j_nu(xj_2) pi-(x1) / 2 - pi0(x2)^dag j_mu(xj_1) j_nu(xj_2) pi0(x1) >:" ) expr = ( sympy.simplify(1) / 2 * mk_pi_p("x2", True) * mk_j_mu("xj_1", "mu") * mk_j_mu("xj_2", "nu") * mk_pi_p("x1") + sympy.simplify(1) / 2 * mk_pi_m("x2", True) * mk_j_mu("xj_1", "mu") * mk_j_mu("xj_2", "nu") * mk_pi_m("x1") - mk_pi_0("x1", True) * mk_j_mu("xj_1", "mu") * mk_j_mu("xj_2", "nu") * mk_pi_0("x2") ) print(display_cexpr(contract_simplify_compile(expr))) print() print( "< pi+(x2)^dag dm(xj_1) dm(xj_2) pi+(x1) / 2 + pi-(x2)^dag dm(xj_1) dm(xj_2) pi-(x1) / 2 - pi0(x2)^dag dm(xj_1) dm(xj_2) pi0(x1) >:" ) expr_dm1 = mk_scalar("d", "d", "xj_1") - mk_scalar("u", "u", "xj_1") expr_dm2 = mk_scalar("d", "d", "xj_2") - mk_scalar("u", "u", "xj_2") expr_dm = expr_dm1 * expr_dm2 expr = ( sympy.simplify(1) / 2 * mk_pi_p("x2", True) * expr_dm * mk_pi_p("x1") + sympy.simplify(1) / 2 * mk_pi_m("x2", True) * expr_dm * mk_pi_m("x1") - mk_pi_0("x1", True) * expr_dm * mk_pi_0("x2") ) print(display_cexpr(contract_simplify_compile(expr))) print() print( "< pi+(x_2)^dag j_mu(xj_1) j_nu(xj_2) pi+(x_1) + pi-(x_2)^dag j_mu(xj_1) j_nu(xj_2) pi-(x_1) + pi0(x_2)^dag j_mu(xj_1) j_nu(xj_2) pi0(x_1) >:" ) expr = ( sympy.simplify(1) * mk_pi_p("x_2", True) * mk_jl_mu("xj_1", "mu") * mk_jl_mu("xj_2", "nu") * mk_pi_p("x_1") + sympy.simplify(1) * mk_pi_m("x_2", True) * mk_jl_mu("xj_1", "mu") * mk_jl_mu("xj_2", "nu") * mk_pi_m("x_1") + mk_pi_0("x_1", True) * mk_jl_mu("xj_1", "mu") * mk_jl_mu("xj_2", "nu") * mk_pi_0("x_2") ) diagram_type_dict = dict() diagram_type_dict[ ( (("x_1", "xj_1"), 1), (("x_2", "xj_2"), 1), (("xj_1", "x_1"), 1), (("xj_2", "x_2"), 1), ) ] = "Type1" diagram_type_dict[ ( (("x_1", "xj_1"), 1), (("x_2", "xj_2"), 1), (("xj_1", "x_2"), 1), (("xj_2", "x_1"), 1), ) ] = "Type2" diagram_type_dict[ ( (("x_1", "x_2"), 1), (("x_2", "xj_1"), 1), (("xj_1", "xj_2"), 1), (("xj_2", "x_1"), 1), ) ] = "Type3" diagram_type_dict[ ( (("x_1", "x_2"), 1), (("x_2", "x_1"), 1), (("xj_1", "xj_2"), 1), (("xj_2", "xj_1"), 1), ) ] = "Type4" print( display_cexpr( contract_simplify_compile(expr, diagram_type_dict=diagram_type_dict) ) ) # print() print("< KpiI3/2_Iz1/2(x2_1,x2_2)^dag KpiI1/2_Iz1/2(x1_1,x1_2) >:") expr = mk_kpi_0_i3halves("x2_1", "x2_2", True) * mk_kpi_0_i1half("x1_1", "x1_2") print(display_cexpr(contract_simplify_compile(expr))) print() print( "< KpiI3/2_Iz3/2(x2_1,x2_2)^dag KpiI3/2_Iz3/2(x1_1,x1_2) > - KpiI3/2_Iz1/2(x2_1,x2_2)^dag KpiI3/2_Iz1/2(x1_1,x1_2) >:" ) expr = mk_kpi_p2_i3halves("x2_1", "x2_2", True) * mk_kpi_p2_i3halves( "x1_1", "x1_2" ) - mk_kpi_p1_i3halves("x2_1", "x2_2", True) * mk_kpi_p1_i3halves("x1_1", "x1_2") print(display_cexpr(contract_simplify_compile(expr))) print() print( "< KpiI3/2_Iz3/2(x2_1,x2_2)^dag KpiI3/2_Iz3/2(x1_1,x1_2) > - KpiI3/2_Iz-1/2(x2_1,x2_2)^dag KpiI3/2_Iz-1/2(x1_1,x1_2) >:" ) expr = mk_kpi_p2_i3halves("x2_1", "x2_2", True) * mk_kpi_p2_i3halves( "x1_1", "x1_2" ) - mk_kpi_0_i3halves("x2_1", "x2_2", True) * mk_kpi_0_i3halves("x1_1", "x1_2") print(display_cexpr(contract_simplify_compile(expr))) print() print( "< KpiI3/2_Iz3/2(x2_1,x2_2)^dag KpiI3/2_Iz3/2(x1_1,x1_2) > - KpiI3/2_Iz-3/2(x2_1,x2_2)^dag KpiI3/2_Iz-3/2(x1_1,x1_2) >:" ) expr = mk_kpi_p2_i3halves("x2_1", "x2_2", True) * mk_kpi_p2_i3halves( "x1_1", "x1_2" ) - mk_kpi_m_i3halves("x2_1", "x2_2", True) * mk_kpi_m_i3halves("x1_1", "x1_2") print(display_cexpr(contract_simplify_compile(expr))) print() print( "< KpiI1/2_Iz1/2(x2_1,x2_2)^dag KpiI1/2_Iz1/2(x1_1,x1_2) > - < KpiI1/2_Iz-1/2(x2_1,x2_2)^dag KpiI1/2_Iz-1/2(x1_1,x1_2) >:" ) expr = mk_kpi_0_i1half("x2_1", "x2_2", True) * mk_kpi_0_i1half( "x1_1", "x1_2" ) - mk_kpi_p_i1half("x2_1", "x2_2", True) * mk_kpi_p_i1half("x1_1", "x1_2") print(display_cexpr(contract_simplify_compile(expr))) print()