Init¶

In [1]:

from pymule import *
from pymule.maths import Li2

Populating the interactive namespace from numpy and matplotlib

Let us normalise all results to the LO decay rate

In [2]:

r = Mel/Mmu
gamma0 = ((Mmu**5)/(192*pi**3))*(1 - 8*r**2 - 24*r**4*log(r) + 8*r**6 - r**8)

In [3]:

setup(folder='nlo-test/out.tar.bz2')
lo = scaleset(mergefks(sigma('m2eNN0')), 1/gamma0)
nlo = scaleset(mergefks(sigma('m2eNNF'), sigma('m2eNNR')), alpha/gamma0)

Equations from [2112.12395]¶

In [4]:

def tomalakLONue(ENue):
    """
    Compute the differentiate decay rate with respect to the anti-neutrino
    electron energy at LO (formula (37) of Tomalak's paper)
    qsqNue is the square of the transfered momentum from the muon to the
    electronic anti-neutrino
    ENue = electronic anti-neutrino energy
    """
    qsqNue = Mmu**2 - 2*Mmu*ENue
    return (1/(2*np.pi**3))*(((qsqNue-Mel**2)**2)/qsqNue)*ENue**2

In [5]:

def tomalakLONumu(ENumu):
    """
    Compute the differentiate decay rate with respect to the muonic neutrino
    energy at LO (formula (38) of Tomalak's paper)
    qsqNumu is the square of the transfered momentum from the muon to the
    muonic neutrino
    ENue = muonic neutrino energy
    """
    qsqNumu = Mmu**2 - 2*Mmu*ENumu
    a = (1/(np.pi**3))*(((qsqNumu-Mel**2)/qsqNumu)**2)*ENumu**2
    b = (qsqNumu - Mel**2)/12 + (qsqNumu + 2*Mel**2)*Mmu*(Mmu - ENumu)/(6*qsqNumu)
    return a*b

In [6]:

def tomalakNLONue(ENue):
    """
    Compute the differentiate decay rate with respect to the anti-neutrino
    electron energy at NLO (formula (55) of Tomalak's paper) for massless
    electron
    qsqNue is the square of the transfered momentum from the muon to the
    electronic anti-neutrino
    ENue = electronic anti-neutrino energy
    """
    qsqNue = Mmu**2 - 2*Mmu*ENue
    a = Li2(2*ENue/Mmu) + 0.5*np.log(1 - 2*ENue/Mmu)**2 + (np.pi**2)/3 - 19/24 + (5/24)*(Mmu/ENue)
    b = (2/3 + (1/3)*(Mmu/ENue) + (5/48)*((Mmu/ENue)**2))*np.log(1 - 2*ENue/Mmu)
    return (-alpha/np.pi)*(a + b)*tomalakLONue(ENue)

def tomalakNLONumu(ENumu):
    """
    Compute the differentiate decay rate with respect to the muonic neutrino
    energy at LO (formula (38) of Tomalak's paper) for massless electrons
    qsqNumu is the square of the transfered momentum from the muon to the
    muonic neutrino
    ENue = muonic neutrino energy
    """
    qsqNumu = Mmu**2 - 2*Mmu*ENumu
    a = Li2(2*ENumu/Mmu) + 0.5*np.log(1 - 2*ENumu/Mmu)**2 + (np.pi**2)/3
    b = ( (43/6)*(ENumu/Mmu) - 51/8 + (41/24)*(Mmu/ENumu) )/( 3 - 4*ENumu/Mmu)
    c = - ( (8/3)*(ENumu/Mmu) - 7/2 + (3/2)*(Mmu/ENumu) - (41/48)*((Mmu/ENumu)**2) )* np.log(1 - 2*ENumu/Mmu)/(3 - 4*ENumu/Mmu)
    return (-alpha/np.pi)*(a + b + c)*tomalakLONumu(ENumu)

In [7]:

wNu = (Mmu**2 - Mel**2)/(2*Mmu)

Compare with [2112.12395]¶

Compare $E(\nu_e)$¶

In [8]:

fig, (ax1,ax2) = kplot(
    {
        'lo': lo['EnuE'],
        'nlo': nlo['EnuE']
    },
    r'$E(\nu_e)$',
    labelsigma=r'$\D\mathcal{B}/\D E(\nu_e)$',
    show=[0],
    legend={}
)

Enu = np.linspace(0, wNu,1000)
ax1.plot(Enu, tomalakLONue(Enu)/gamma0, 'C2', linestyle='dashed')
ax2.plot(Enu, tomalakNLONue(Enu)/tomalakLONue(Enu), 'C0', linestyle='dashed')
ax2.set_ylim(-0.02,0.001)
ax2.legend(
    [
        matplotlib.lines.Line2D([0],[0],color='C0', linestyle=ls)
        for ls in ['solid', 'dashed']
    ],
    [r'\textsc{McMule} NLO', 'Tomalak 21']
)
watermark(fig)
fig.savefig("enue.pdf")

<ipython-input-6-f7a29defc315>:11: RuntimeWarning: divide by zero encountered in true_divide
  a = Li2(2*ENue/Mmu) + 0.5*np.log(1 - 2*ENue/Mmu)**2 + (np.pi**2)/3 - 19/24 + (5/24)*(Mmu/ENue)
<ipython-input-6-f7a29defc315>:12: RuntimeWarning: divide by zero encountered in true_divide
  b = (2/3 + (1/3)*(Mmu/ENue) + (5/48)*((Mmu/ENue)**2))*np.log(1 - 2*ENue/Mmu)
<ipython-input-6-f7a29defc315>:12: RuntimeWarning: invalid value encountered in multiply
  b = (2/3 + (1/3)*(Mmu/ENue) + (5/48)*((Mmu/ENue)**2))*np.log(1 - 2*ENue/Mmu)

Compare $E(\nu_\mu)$¶

In [9]:

fig, (ax1,ax2) = kplot(
    {
        'lo': lo['EnuMu'],
        'nlo': nlo['EnuMu']
    },
    r'$E(\nu_\mu)$',
    labelsigma=r'$\D\mathcal{B}/\D E(\nu_\mu)$',
    show=[0],
    legend={}
)

Enu = np.linspace(0, wNu,1000)
ax1.plot(Enu, tomalakLONumu(Enu)/gamma0, 'C2', linestyle='dashed')
ax2.plot(Enu, tomalakNLONumu(Enu)/tomalakLONumu(Enu), 'C0', linestyle='dashed')
ax2.set_ylim(-0.02,0.001)
ax2.legend(
    [
        matplotlib.lines.Line2D([0],[0],color='C0', linestyle=ls)
        for ls in ['solid', 'dashed']
    ],
    [r'\textsc{McMule} NLO', 'Tomalak 21']
)
watermark(fig)
fig.savefig("enumu.pdf")

<ipython-input-6-f7a29defc315>:25: RuntimeWarning: divide by zero encountered in true_divide
  b = ( (43/6)*(ENumu/Mmu) - 51/8 + (41/24)*(Mmu/ENumu) )/( 3 - 4*ENumu/Mmu)
<ipython-input-6-f7a29defc315>:26: RuntimeWarning: divide by zero encountered in true_divide
  c = - ( (8/3)*(ENumu/Mmu) - 7/2 + (3/2)*(Mmu/ENumu) - (41/48)*((Mmu/ENumu)**2) )* np.log(1 - 2*ENumu/Mmu)/(3 - 4*ENumu/Mmu)
<ipython-input-6-f7a29defc315>:26: RuntimeWarning: invalid value encountered in subtract
  c = - ( (8/3)*(ENumu/Mmu) - 7/2 + (3/2)*(Mmu/ENumu) - (41/48)*((Mmu/ENumu)**2) )* np.log(1 - 2*ENumu/Mmu)/(3 - 4*ENumu/Mmu)