Skip to content

SMEFT 4-fermion corrections #306

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Draft
wants to merge 8 commits into
base: master
Choose a base branch
from
Draft
Changes from 6 commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
97 changes: 93 additions & 4 deletions src/yadism/coefficient_functions/coupling_constants.py
Original file line number Diff line number Diff line change
Expand Up @@ -46,6 +46,24 @@ def __init__(self, theory_config, obs_config):
# neutrinos
for pid in [12, 14, 16]:
self.weak_isospin_3[pid] = 1 / 2

# BSM couplings -------------------------------------------------------
# Temporary couplings to Olq3
self.BSM_coupling_vectorial = {21: 0}
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I wonder is "BSM" too generic? i.e. should we be more specific, like e.g. "4F"*? please tell me! because from my non-existent BSM experience I only know there are many BSM models 🙃

*Python variables can not start with a number, but we can solve that problem 🙃

Copy link
Collaborator Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Indeed you are right! I realise I have already started mixing up BSM and 4F labels, 4F is a more precise one and I'll try to stick with it

Copy link
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I'd even be more specific and use Olq3, as you can have multiple 4F operators right?

self.BSM_coupling_axial = {21: 0}
for q in range(1, 7):
self.BSM_coupling_vectorial[q] = -1 if q % 2 == 0 else 1 # u if stmt else d
self.BSM_coupling_axial[q] = -1 if q % 2 == 0 else 1 # u if stmt else d
# leptons: 11 = e-(!)
for pid in [11, 13, 15]:
self.BSM_coupling_vectorial[pid] = 1
self.BSM_coupling_axial[pid] = 1
# neutrinos
# CC not yet implemented
for pid in [12, 14, 16]:
self.BSM_coupling_vectorial[pid] = 0
self.BSM_coupling_axial[pid] = 0

self.log()

def log(self):
Expand Down Expand Up @@ -92,9 +110,16 @@ def leptonic_coupling(self, mode, quark_coupling_type):
# load Z coupling
projectile_v = 0.0
projectile_a = 0.0
if mode in ["phZ", "ZZ"]:
projectile_BSM_v = 0.0
projectile_BSM_a = 0.0
if mode in ["phZ", "ZZ", "Z4F"]:
projectile_v = self.vectorial_coupling(abs(projectile_pid))
projectile_a = self.weak_isospin_3[abs(projectile_pid)]
projectile_BSM_v = self.BSM_coupling_vectorial[abs(projectile_pid)]
Copy link
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Maybe for consistency you should also rename projectile_v to projectile_Z_v and same for _a.
Or maybe this has to be even generalized more if you have multiple operators...

Copy link
Collaborator Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Hi Giacomo,
This makes good sense, thanks! I have implemented it

projectile_BSM_a = self.BSM_coupling_axial[abs(projectile_pid)]
if mode in ["ph4F"]:
projectile_BSM_v = self.BSM_coupling_vectorial[abs(projectile_pid)]
projectile_BSM_a = self.BSM_coupling_axial[abs(projectile_pid)]
# switch mode
if mode == "phph":
if quark_coupling_type in ["VV", "AA"]:
Expand All @@ -121,6 +146,36 @@ def leptonic_coupling(self, mode, quark_coupling_type):
return 2.0 * projectile_v * projectile_a + pol * (
projectile_v**2 + projectile_a**2
)
elif mode == "ph4F":
if quark_coupling_type in ["VV", "AA"]:
return self.electric_charge[abs(projectile_pid)] * (
projectile_BSM_v + pol * projectile_BSM_a
)
else:
return self.electric_charge[abs(projectile_pid)] * (
projectile_BSM_a + pol * projectile_BSM_v
)
elif mode == "Z4F":
if quark_coupling_type in ["VV", "AA"]:
return (
projectile_v * projectile_BSM_v
+ projectile_a * projectile_BSM_a
+ pol
* (
projectile_v * projectile_BSM_a
+ projectile_BSM_v * projectile_a
)
)
else:
return (
projectile_v * projectile_BSM_a
+ projectile_BSM_v * projectile_a
+ pol
* (
projectile_v * projectile_BSM_v
+ projectile_a * projectile_BSM_a
)
)
raise ValueError(f"Unknown mode: {mode}")

def nc_partonic_coupling(self, pid):
Expand All @@ -130,7 +185,12 @@ def nc_partonic_coupling(self, pid):
"A": 0, # axial coupling of the photon to the quark is not there of course
}
qZ = {"V": self.vectorial_coupling(pid), "A": self.weak_isospin_3[pid]}
return qph, qZ

q4F = {
"V": self.BSM_coupling_vectorial[pid],
"A": self.BSM_coupling_axial[pid],
}
return qph, qZ, q4F

def partonic_coupling(self, mode, pid, quark_coupling_type, cc_mask=None):
"""Compute the coupling of the boson to the parton.
Expand All @@ -157,7 +217,7 @@ def partonic_coupling(self, mode, pid, quark_coupling_type, cc_mask=None):
if mode == "WW":
return np.sum(self.theory_config["CKM"].masked(cc_mask)(pid))
# load couplings
qph, qZ = self.nc_partonic_coupling(pid)
qph, qZ, q4F = self.nc_partonic_coupling(pid)

first = quark_coupling_type[0]
second = quark_coupling_type[1]
Expand All @@ -168,6 +228,10 @@ def partonic_coupling(self, mode, pid, quark_coupling_type, cc_mask=None):
return qph[first] * qZ[second]
if mode == "ZZ":
return qZ[first] * qZ[second]
if mode == "ph4F":
return qph[first] * q4F[second]
if mode == "Z4F":
return qZ[first] * q4F[second]
raise ValueError(f"Unknown mode: {mode}")

def partonic_coupling_fl11(self, mode, pid, nf, quark_coupling_type):
Expand Down Expand Up @@ -247,10 +311,22 @@ def propagator_factor(self, mode, Q2):
* (1.0 - self.theory_config["sin2theta_weak"])
)
eta_phZ /= 1 - self.obs_config["propagatorCorrection"]

# Need to specify Wilson coefficient?
C_4F = 1 / 4
# Modify with more precise value
alpha = 1 / 137
# Should get it from param card
BSM_scale = 1000
eta_ph4F = ((C_4F) / (4 * np.pi * alpha)) * (Q2 / BSM_scale**2)
if mode == "phZ":
return eta_phZ
if mode == "ZZ":
return eta_phZ**2
if mode == "ph4F":
return eta_ph4F
if mode == "Z4F":
return eta_phZ * eta_ph4F
if mode == "WW":
eta_W = (
(eta_phZ / 2)
Expand Down Expand Up @@ -327,7 +403,20 @@ def get_weight(self, pid, Q2, quark_coupling_type, cc_mask=None):
* self.propagator_factor("ZZ", Q2)
* self.partonic_coupling("ZZ", pid, quark_coupling_type)
)
return w_phph + w_phZ + w_ZZ
# photon-4F interference
w_ph4F = (
2
* self.leptonic_coupling("ph4F", quark_coupling_type)
* self.propagator_factor("ph4F", Q2)
* self.partonic_coupling("ph4F", pid, quark_coupling_type)
)
# Z-4F interference
w_Z4F = (
self.leptonic_coupling("Z4F", quark_coupling_type)
* self.propagator_factor("Z4F", Q2)
* self.partonic_coupling("Z4F", pid, quark_coupling_type)
)
return w_phph + w_phZ + w_ZZ + w_ph4F + w_Z4F
raise ValueError(f"Unknown process: {self.obs_config['process']}")

def get_fl11_weight(self, pid, Q2, nf, quark_coupling_type):
Expand Down
Loading