Arithmetic and algebra

This SageMath notebook does arithmetic and algebra for the paper A theory of the dark matter.

The section numbering follows the paper. Equation numbers refer to equations in the paper.

Preamble

Code for displaying equations.

%display latex
LE = lambda latex_string: LatexExpr(latex_string);

The numerical values of fundamental constants and basic physical quantities are stored in the dictionary value[].

Physical units such as 's' and 'GeV' are defined as algebraic variables.

The dictionary formula[] will contain the formulas for derived physical quantities. For example, formula[kappa] = 8*pi*G.

The function valof(x) substitutes in a formula x to obtain a numerical value.

The function print_values(x1,x2,...) prints formulas xn with their numerical values.

value = {}
formula={}

def valof(x):
    SymRing = FractionField(PolynomialRing(RR,'s,GeV,J,m,kg,K,C'))
    NumRing = FractionField(PolynomialRing(RDF,'s,GeV,J,m,kg,K,C'))
    xval = SR(x)
    for a in reversed(formula):
        xval = xval.subs({a:formula[a]})
    for key,val in value.items():
        xval = xval.subs({key:val})
    xval = NumRing(xval)
    return xval
def print_values(*args):
    print_list=[]
    for arg in args:
        if arg in formula:
            print_list += [arg,'=',formula[arg],'=',valof(arg),LE(r'\qquad')]
        else:
            print_list += [arg,'=',valof(arg),LE(r'\qquad')]
    print_list.pop()
    pretty_print(*print_list)

1.3 Physical parameters

Units and fundamental constants as variables

declare units as variables

s = var('s', domain='positive',latex_name='\\mathrm{s}'); assume(s,'real');
GeV = var('GeV', domain='positive',latex_name='\\mathrm{GeV}'); assume(GeV,'real');
J = var('J', domain='positive',latex_name='\\mathrm{J}'); assume(J,'real');
m = var('m', domain='positive',latex_name='\\mathrm{m}'); assume(m,'real');
kg = var('kg', domain='positive',latex_name='\\mathrm{kg}\\,'); assume(kg,'real');
K = var('K', domain='positive',latex_name='\\mathrm{K}'); assume(K,'real');
C = var('C', domain='positive',latex_name='\\mathrm{C}'); assume(C,'real');
km = var('km', domain='positive',latex_name='\\mathrm{km}\\,'); assume(km,'real');
Mpc = var('Mpc', domain='positive',latex_name='\\mathrm{Mpc}\\,'); assume(Mpc,'real');

declare fundamental constants and physical parameters as variables

# fundamental constants
c = var('c',latex_name=r"c")
c_mps = var('c_mps',latex_name=r"c")
e_charge = var('e_charge',latex_name=r"e")
hbar = var('hbar',latex_name=r"\hbar")
kB = var('kB',latex_name=r"k_{B}")
G = var('G',latex_name=r"G")
kappa = var('kappa',latex_name=r"\kappa")
#
# numbered constants c_n for temporary use
# c_ = var('c_',n=20,latex_name=r"c")
#
# Standard Model
GFermi = var('GFermi',latex_name=r"G_{\mathrm{Fermi}}")
m_W = var('m_W',latex_name=r"m_{W}")
m_Higgs = var('m_Higgs',latex_name=r"m_{\mathrm{Higgs}}")
hbar_v = var('hbar_v',latex_name=r"\hbar v")
v = var('v',latex_name=r"v")
g = var('g',latex_name=r"g")
lambdaH = var('lambdaH',latex_name=r"\lambda")
#
# cosmological parameters
H0 = var('H0',latex_name=r"H_{0}")
Omega_curvature = var('Omega_curvature',latex_name=r"\Omega_{\mathrm{curvature}}")
Omega_Lambda = var('Omega_Lambda',latex_name=r"\Omega_{\Lambda}")
#
# time and energy scales
t_grav = var('t_grav',latex_name=r"t_{\mathrm{grav}}")
t_Higgs = var('t_Higgs',latex_name=r"t_{\mathrm{Higgs}}")
t_Hubble = var('t_Hubble',latex_name=r"t_{\mathrm{Hubble}}")

Values of the fundamental constants from NIST 2018

value[G]        = 6.67430e-11*m^3*kg^(-1)* s^(-2)
value[e_charge] = 1.602176634e-19 * C
value[c]=c_mps  = 2.99792458e8 * m * s^(-1)
value[hbar]     = 1.054571817e-34*J*s
value[kB]       = 1.380649e-23*J*K^(-1)

print_values(G,e_charge,c)
print('\n')
print_values(hbar,kB)

\[\newcommand{\Bold}[1]{\mathbf{#1}}{G} \verb|=| \frac{6.6743 \times 10^{-11} m^{3}}{s^{2} \mathit{kg}} \qquad {e} \verb|=| 1.602176634 \times 10^{-19} C \qquad {c} \verb|=| \frac{299792458.0 m}{s}\]

\[\newcommand{\Bold}[1]{\mathbf{#1}}{\hbar} \verb|=| 1.054571817 \times 10^{-34} s J \qquad {k_{B}} \verb|=| \frac{1.380649 \times 10^{-23} J}{K}\]

The constant $\kappa = 8 \pi G$.

formula[kappa]  = 8*pi*G
print_values(kappa)

\[\newcommand{\Bold}[1]{\mathbf{#1}}{\kappa} \verb|=| 8 \, \pi {G} \verb|=| \frac{1.6774345478283484 \times 10^{-09} m^{3}}{s^{2} \mathit{kg}}\]

convert to c=1 units with unit of energy = GeV

conv={}
conv[m]  = c^(-1)*m
conv[kg] = J*s^2*m^(-2)
conv[J]  = e_charge^(-1) * C * 1e-9 * GeV
for a,b in conv.items():
    value[a]= a.subs(conv).subs(value)
#
print_values(c,J)
print('\n')
print_values(m,kg)
print('\n')
print_values(hbar,kB)
print('\n')
print_values(G,kappa)

\[\newcommand{\Bold}[1]{\mathbf{#1}}{c} \verb|=| 1.0 \qquad {\mathrm{J}} \verb|=| 6241509074.460764 \mathit{GeV}\]

\[\newcommand{\Bold}[1]{\mathbf{#1}}{\mathrm{m}} \verb|=| 3.3356409519815204 \times 10^{-09} s \qquad {\mathrm{kg}\,} \verb|=| 5.6095886038044526 \times 10^{26} \mathit{GeV}\]

\[\newcommand{\Bold}[1]{\mathbf{#1}}{\hbar} \verb|=| 6.582119565476076 \times 10^{-25} s \mathit{GeV} \qquad {k_{B}} \verb|=| \frac{8.61733326214518 \times 10^{-14} \mathit{GeV}}{K}\]

\[\newcommand{\Bold}[1]{\mathbf{#1}}{G} \verb|=| \frac{4.415832614329417 \times 10^{-63} s}{\mathit{GeV}} \qquad {\kappa} \verb|=| 8 \, \pi {G} \verb|=| \frac{1.1098197840527606 \times 10^{-61} s}{\mathit{GeV}}\]

Standard Model coupling constants from PDG (2020, 2021)

measured quantities

value[GFermi]  = 1.1663787e-5*GeV^(-2);
value[m_W]     = 80.379*GeV;
value[m_Higgs] = 125.10*GeV;
#
print_values(GFermi,m_W,m_Higgs)

\[\newcommand{\Bold}[1]{\mathbf{#1}}{G_{\mathrm{Fermi}}} \verb|=| \frac{1.1663787 \times 10^{-05}}{\mathit{GeV}^{2}} \qquad {m_{W}} \verb|=| 80.379 \mathit{GeV} \qquad {m_{\mathrm{Higgs}}} \verb|=| 125.1 \mathit{GeV}\]

derived quantities

#
formula[hbar_v] = 2^(-1/4)*GFermi^(-1/2)
formula[v] = formula[hbar_v]/hbar
formula[g] = 2*m_W/hbar_v;
formula[lambdaH] = m_Higgs/hbar_v;
#
print_values(hbar_v,v)
print('\n')
print_values(g,g^2)
print('\n')
print_values(lambdaH,lambdaH^2)

\[\newcommand{\Bold}[1]{\mathbf{#1}}{\hbar v} \verb|=| \frac{2^{\frac{3}{4}}}{2 \, \sqrt{{G_{\mathrm{Fermi}}}}} \verb|=| 246.2196507941374 \mathit{GeV} \qquad {v} \verb|=| \frac{2^{\frac{3}{4}}}{2 \, \sqrt{{G_{\mathrm{Fermi}}}} {\hbar}} \verb|=| \frac{3.740735007087775 \times 10^{26}}{s}\]

\[\newcommand{\Bold}[1]{\mathbf{#1}}{g} \verb|=| \frac{2 \, {m_{W}}}{{\hbar v}} \verb|=| 0.6529048330687819 \qquad {g}^{2} \verb|=| 0.4262847210445737\]

\[\newcommand{\Bold}[1]{\mathbf{#1}}{\lambda} \verb|=| \frac{{m_{\mathrm{Higgs}}}}{{\hbar v}} \verb|=| 0.5080829235055462 \qquad {\lambda}^{2} \verb|=| 0.25814825715794265\]

Cosmological parameters

units of distance

The value of 1 parsec in meters is taken from IAU 2015 Resolution B2, note 4 which references the definition as exactly 64000/$\pi$ au and the definition of au from IAU 2012 Resolution B2.

value[km]       = 1e3 *m
value[Mpc]      = 1e6 * 3.085677581 * 1e16 * m
print_values(km,Mpc)

\[\newcommand{\Bold}[1]{\mathbf{#1}}{\mathrm{km}\,} \verb|=| 1000.0 m \qquad {\mathrm{Mpc}\,} \verb|=| 3.0856775810000003 \times 10^{22} m\]

From Particle Data Group 2020 Particle Physics Booklet

formula[H0] = 6.74e1 *km * s^(-1) * Mpc^(-1)
value[Omega_curvature] = 0.001
value[Omega_Lambda] = 0.685
#
print_values(H0,Omega_Lambda,Omega_curvature)

\[\newcommand{\Bold}[1]{\mathbf{#1}}{H_{0}} \verb|=| \frac{67.4000000000000 \, {\mathrm{km}\,}}{{\mathrm{Mpc}\,} {\mathrm{s}}} \verb|=| \frac{2.1842852414333302 \times 10^{-18}}{s} \qquad {\Omega_{\Lambda}} \verb|=| 0.685 \qquad {\Omega_{\mathrm{curvature}}} \verb|=| 0.001\]

time and energy scales

formula[t_grav] = sqrt(kappa*hbar)
formula[t_Higgs] = hbar/m_Higgs
formula[t_Hubble] = 1/H0
#
print_values(t_grav,hbar/t_grav)
print('\n')
print_values(t_Higgs,m_Higgs)
print('\n')
print_values(t_Hubble,hbar/t_Hubble)

\[\newcommand{\Bold}[1]{\mathbf{#1}}{t_{\mathrm{grav}}} \verb|=| \sqrt{{\hbar} {\kappa}} \verb|=| 2.7027701557413476 \times 10^{-43} s \qquad \frac{{\hbar}}{{t_{\mathrm{grav}}}} \verb|=| 2.435323459338203 \times 10^{18} \mathit{GeV}\]

\[\newcommand{\Bold}[1]{\mathbf{#1}}{t_{\mathrm{Higgs}}} \verb|=| \frac{{\hbar}}{{m_{\mathrm{Higgs}}}} \verb|=| 5.2614864632102925 \times 10^{-27} s \qquad {m_{\mathrm{Higgs}}} \verb|=| 125.1 \mathit{GeV}\]

\[\newcommand{\Bold}[1]{\mathbf{#1}}{t_{\mathrm{Hubble}}} \verb|=| \frac{1}{{H_{0}}} \verb|=| 4.5781566483679526 \times 10^{17} s \qquad \frac{{\hbar}}{{t_{\mathrm{Hubble}}}} \verb|=| 1.4377226624218956 \times 10^{-42} \mathit{GeV}\]

2.2 Initial CGF energy

elliptic parameter $k^{2}_{\mathrm{EW}}=\frac12$ and elliptic integral of first kind $K_{\mathrm{EW}}=K(k_{\mathrm{EW}})=K'(k_{\mathrm{EW}})$

# elliptic modulus and integrals
ksq_EW = var('ksq_EW',latex_name=r"k^2_{\mathrm{EW}}")
k_EW = var('k_EW',latex_name=r"k_{\mathrm{EW}}")
K_EW = var('K_EW',latex_name=r"K_{\mathrm{EW}}")
Kp_EW = var('Kp_EW',latex_name=r"K'_{\mathrm{EW}}")
#
value[ksq_EW]= 1/2
value[k_EW]= sqrt(1/2)
value[K_EW]= elliptic_kc(0.5)
value[Kp_EW]= elliptic_kc(0.5)
print_values(ksq_EW,K_EW,Kp_EW)

\[\newcommand{\Bold}[1]{\mathbf{#1}}{k^2_{\mathrm{EW}}} \verb|=| 0.5 \qquad {K_{\mathrm{EW}}} \verb|=| 1.8540746773013719 \qquad {K'_{\mathrm{EW}}} \verb|=| 1.8540746773013719\]

2.4 Start of the electroweak transition at $a=a_{\mathrm{EW}}$

a_EW = var('a_EW',latex_name=r"a_{\mathrm{EW}}")
#
a_EW_coeff = (6*pi)^(1/2)/(4*K_EW)
formula[a_EW] = a_EW_coeff * t_Higgs
print_values(a_EW_coeff,a_EW)

\[\newcommand{\Bold}[1]{\mathbf{#1}}\frac{\sqrt{6} \sqrt{\pi}}{4 \, {K_{\mathrm{EW}}}} \verb|=| 0.5854143283037644 \qquad {a_{\mathrm{EW}}} \verb|=| \frac{\sqrt{6} \sqrt{\pi} {t_{\mathrm{Higgs}}}}{4 \, {K_{\mathrm{EW}}}} \verb|=| 3.0801495637396022 \times 10^{-27} s\]

ahat_EW = var('ahat_EW',latex_name=r"\hat a_{\mathrm{EW}}")
formula[ahat_EW] = a_EW/t_Higgs
print_values(ahat_EW)

\[\newcommand{\Bold}[1]{\mathbf{#1}}{\hat a_{\mathrm{EW}}} \verb|=| \frac{{a_{\mathrm{EW}}}}{{t_{\mathrm{Higgs}}}} \verb|=| 0.5854143283037644\]

2.6 Realizing the electroweak transition

asqHsq_EW = var('asqHsq_EW',\
    latex_name=r"a_{\mathrm{EW}}^2 H_{\mathrm{EW}}^2+\epsilon^2")
formula[asqHsq_EW ] = (t_grav^2/t_Higgs^2) * ( (1/(24*lambdaH^2))\
    * (a_EW^2/t_Higgs^2) +(1/(8*g^2))* (t_Higgs^2/a_EW^2) )
print_values(asqHsq_EW)
print('\n')
expr = (t_grav^2/t_Higgs^2)/(4*sqrt(3)*g*lambdaH)
print_values(expr)

\[\newcommand{\Bold}[1]{\mathbf{#1}}{a_{\mathrm{EW}}^2 H_{\mathrm{EW}}^2+\epsilon^2} \verb|=| \frac{{t_{\mathrm{grav}}}^{2} {\left(\frac{{a_{\mathrm{EW}}}^{2}}{{\lambda}^{2} {t_{\mathrm{Higgs}}}^{2}} + \frac{3 \, {t_{\mathrm{Higgs}}}^{2}}{{a_{\mathrm{EW}}}^{2} {g}^{2}}\right)}}{24 \, {t_{\mathrm{Higgs}}}^{2}} \verb|=| 2.4037614729587787 \times 10^{-33}\]

\[\newcommand{\Bold}[1]{\mathbf{#1}}\frac{\sqrt{3} {t_{\mathrm{grav}}}^{2}}{12 \, {g} {\lambda} {t_{\mathrm{Higgs}}}^{2}} \verb|=| 1.1481436122987478 \times 10^{-33}\]

4.5 Parametrize the time evolution by $k^2$

PSR.<k> = PowerSeriesRing(SR)
K_ps = (pi/2)*(1+k^2/4 + 9*k^4/64)+O(k^6)
E_ps = (pi/2)*(1-k^2/4 - 3*k^4/64)+O(k^6)

alpha3 = (2*(1-k^2)*K_ps + 2*(2*k^2-1)*E_ps)/K_EW
pretty_print(LE(r"\alpha^3 = "),alpha3+ O(k^4))

\[\newcommand{\Bold}[1]{\mathbf{#1}}\alpha^3 = \frac{3 \, \pi}{2 \, {K_{\mathrm{EW}}}} k^{2} + O(k^{4})\]

alpha2_b2 = k^2 -1 +E_ps/K_ps; pretty_print(LE(r"\alpha^2 \langle b^2\rangle = "),alpha2_b2+O(k^4))

\[\newcommand{\Bold}[1]{\mathbf{#1}}\alpha^2 \langle b^2\rangle = \frac{1}{2} k^{2} + O(k^{4})\]

alpha2_mu2 = 1-2*k^2; pretty_print(LE(r"\alpha^2 \mu^2 = "),alpha2_mu2+O(k^2))

\[\newcommand{\Bold}[1]{\mathbf{#1}}\alpha^2 \mu^2 = 1 + O(k^{2})\]

alpha4_E_CGF = k^2*(1-k^2)/2; pretty_print(LE(r"\alpha^4 E_{CGF} = "),alpha4_E_CGF+O(k^4))

\[\newcommand{\Bold}[1]{\mathbf{#1}}\alpha^4 E_{CGF} = \frac{1}{2} k^{2} + O(k^{4})\]

alpha2_ahat2 = (3/2)* alpha2_b2 + (4*lambdaH^2/g^2)* alpha2_mu2
pretty_print(LE(r"\alpha^2 \hat a^2= "),alpha2_ahat2+O(k^2))

\[\newcommand{\Bold}[1]{\mathbf{#1}}\alpha^2 \hat a^2= \frac{4 \, {\lambda}^{2}}{{g}^{2}} + O(k^{2})\]

rhohat_CGF = (3*alpha4_E_CGF/g^2 +9*alpha2_b2^2/(32*lambdaH^2))/alpha2_ahat2^2
pretty_print(LE(r"\hat \rho_{CGF} = "),rhohat_CGF +O(k^4))

\[\newcommand{\Bold}[1]{\mathbf{#1}}\hat \rho_{CGF} = \frac{3 \, {g}^{2}}{32 \, {\lambda}^{4}} k^{2} + O(k^{4})\]

phat_CGF = ((alpha4_E_CGF-alpha2_mu2*alpha2_b2)/g^2 -9*alpha2_b2^2/(32*lambdaH^2))/alpha2_ahat2^2
pretty_print(LE(r"\hat p_{CGF} = "),phat_CGF +O(k^6))

\[\newcommand{\Bold}[1]{\mathbf{#1}}\hat p_{CGF} = \left(\frac{9 \, {g}^{4} {\left(\frac{8}{{g}^{2}} - \frac{1}{{\lambda}^{2}}\right)}}{2048 \, {\lambda}^{4}}\right) k^{4} + O(k^{6})\]

5.1 $\Omega_{\mathrm{CGF}} + \Omega_{\Lambda}=1$

k0sq = var('k0sq',latex_name=r"k_{0}^{2}")
expr = t_Higgs^4/(t_grav^2*t_Hubble^2)*(32*lambdaH^4/g^2)
formula[k0sq] = 0.315*expr
pretty_print(k0sq,LE('= 0.315'),expr,'=',valof(k0sq))

\[\newcommand{\Bold}[1]{\mathbf{#1}}{k_{0}^{2}} = 0.315 \frac{32 \, {\lambda}^{4} {t_{\mathrm{Higgs}}}^{4}}{{g}^{2} {t_{\mathrm{Hubble}}}^{2} {t_{\mathrm{grav}}}^{2}} \verb|=| 7.887392672844198 \times 10^{-56}\]

a0 = var('a0',latex_name=r"a_{0}")
expr= ((6*pi*lambdaH*g/K_EW)*t_Higgs/(t_grav^2*t_Hubble^2))^(-1/3)
formula[a0] = 0.315^(-1/3)*expr
pretty_print(a0,LE(r"= 0.315^{-\frac13}"),expr,'=',valof(a0))

\[\newcommand{\Bold}[1]{\mathbf{#1}}{a_{0}} = 0.315^{-\frac13} \frac{6^{\frac{2}{3}}}{6 \, \left(\frac{\pi {g} {\lambda} {t_{\mathrm{Higgs}}}}{{K_{\mathrm{EW}}} {t_{\mathrm{Hubble}}}^{2} {t_{\mathrm{grav}}}^{2}}\right)^{\frac{1}{3}}} \verb|=| 1.3991814678093661 \times 10^{-08} s\]

print_values(a0/t_Higgs)

\[\newcommand{\Bold}[1]{\mathbf{#1}}\frac{{a_{0}}}{{t_{\mathrm{Higgs}}}} \verb|=| 2.6592893046343716 \times 10^{18}\]

print_values(a0/a_EW)

\[\newcommand{\Bold}[1]{\mathbf{#1}}\frac{{a_{0}}}{{a_{\mathrm{EW}}}} \verb|=| 4.542576387461597 \times 10^{18}\]

print_values(t_Hubble^2/a0^2)

\[\newcommand{\Bold}[1]{\mathbf{#1}}\frac{{t_{\mathrm{Hubble}}}^{2}}{{a_{0}}^{2}} \verb|=| 1.0706147161725437 \times 10^{51}\]

5.2 $w_{\mathrm{CGF}} =0$

w_CGF = phat_CGF/rhohat_CGF
pretty_print(LE(r"w_{CGF} = "),w_CGF +O(k^4))

\[\newcommand{\Bold}[1]{\mathbf{#1}}w_{CGF} = \left(\frac{3}{64} \, {g}^{2} {\left(\frac{8}{{g}^{2}} - \frac{1}{{\lambda}^{2}}\right)}\right) k^{2} + O(k^{4})\]

print_values(K_EW/(2*pi*g*lambdaH))

\[\newcommand{\Bold}[1]{\mathbf{#1}}\frac{{K_{\mathrm{EW}}}}{2 \, \pi {g} {\lambda}} \verb|=| 0.8895346543935904\]

5.3 The CGF in the present

omega_W = var('omega_W',latex_name=r"\omega_W")
formula[omega_W] = g/(2*lambdaH*t_Higgs)
print_values(omega_W,m_W/hbar)

\[\newcommand{\Bold}[1]{\mathbf{#1}}{\omega_W} \verb|=| \frac{{g}}{2 \, {\lambda} {t_{\mathrm{Higgs}}}} \verb|=| \frac{1.221171982678596 \times 10^{26}}{s} \qquad \frac{{m_{W}}}{{\hbar}} \verb|=| \frac{1.221171982678596 \times 10^{26}}{s}\]

print_values(2*pi/omega_W)

\[\newcommand{\Bold}[1]{\mathbf{#1}}\frac{2 \, \pi}{{\omega_W}} \verb|=| 5.145209189452289 \times 10^{-26} s\]

k0 = var('k0',latex_name=r"k_{0}")
formula[k0]= sqrt(k0sq)
print_values(k0)
print('\n')
print_values(k0*omega_W,k0*m_W)

\[\newcommand{\Bold}[1]{\mathbf{#1}}{k_{0}} \verb|=| \sqrt{{k_{0}^{2}}} \verb|=| 2.8084502261646374 \times 10^{-28}\]

\[\newcommand{\Bold}[1]{\mathbf{#1}}{k_{0}} {\omega_W} \verb|=| \frac{0.03429600730939622}{s} \qquad {k_{0}} {m_{W}} \verb|=| 2.257404207288874 \times 10^{-26} \mathit{GeV}\]

5.4 CGF equation of state

rho_b = var('rho_b ',latex_name=r"\rho_{b}")
formula[rho_b] = hbar/t_Higgs^4
pretty_print(rho_b," = ",formula[rho_b] ," = ", kg*valof(rho_b/kg)/c_mps^3)

\[\newcommand{\Bold}[1]{\mathbf{#1}}{\rho_{b}} \phantom{\verb!x!}\verb|=| \frac{{\hbar}}{{t_{\mathrm{Higgs}}}^{4}} \phantom{\verb!x!}\verb|=| \frac{5.682491786274889 \times 10^{28} \, {\mathrm{kg}\,}}{{\mathrm{m}}^{3}}\]

rhohat_EW = var('rhohat_EW ',latex_name=r"\hat\rho_{EW}")
formula[rhohat_EW ] = 8 *K_EW^4/(3*pi^2*g^2) + 1/(8*lambdaH^2)
print_values(rhohat_EW)

\[\newcommand{\Bold}[1]{\mathbf{#1}}{\hat\rho_{EW}} \verb|=| \frac{8 \, {K_{\mathrm{EW}}}^{4}}{3 \, \pi^{2} {g}^{2}} + \frac{1}{8 \, {\lambda}^{2}} \verb|=| 7.97415392287409\]

c_b = var('c_b')
formula[c_b] = 1/(2*lambdaH^2*rhohat_EW)
print_values(c_b)

\[\newcommand{\Bold}[1]{\mathbf{#1}}c_{b} \verb|=| \frac{1}{2 \, {\lambda}^{2} {\hat\rho_{EW}}} \verb|=| 0.24289366755876735\]

c_a = var('c_a')
formula[c_a] = 2*lambdaH^2*(8*lambdaH^2-g^2)/(2*g^2)
c_a_val = valof(c_a)
print_values(c_a)
print(c_a_val)

\[\newcommand{\Bold}[1]{\mathbf{#1}}c_{a} \verb|=| -\frac{{\left({g}^{2} - 8 \, {\lambda}^{2}\right)} {\lambda}^{2}}{{g}^{2}} \verb|=| 0.9924810876684255\]

0.9924810876684255

5.5 Adiabatic condition for $a\ge a_{\mathrm{EW}}$

y=var('y',latex_name=r"4K_{\mathrm{EW}}(a^2H^2+\epsilon^2)^{1/2}_{\mathrm{EW}}")
formula[y] = 4*K_EW*ahat_EW*sqrt((t_grav^2/t_Higgs^2)*rhohat_EW/3+ Omega_Lambda*(t_Higgs^2/t_Hubble^2))
print_values(y)

\[\newcommand{\Bold}[1]{\mathbf{#1}}{4K_{\mathrm{EW}}(a^2H^2+\epsilon^2)^{1/2}_{\mathrm{EW}}} \verb|=| 4 \, {K_{\mathrm{EW}}} {\hat a_{\mathrm{EW}}} \sqrt{\frac{{\Omega_{\Lambda}} {t_{\mathrm{Higgs}}}^{2}}{{t_{\mathrm{Hubble}}}^{2}} + \frac{{\hat\rho_{EW}} {t_{\mathrm{grav}}}^{2}}{3 \, {t_{\mathrm{Higgs}}}^{2}}} \verb|=| 3.6360755535373564 \times 10^{-16}\]

print_values((pi^2/(2*lambdaH^2))*(t_grav^2/t_Higgs^2),\
    (16*pi^2*lambdaH^2/g^2)*Omega_Lambda*(t_Higgs^2/t_Hubble^2))

\[\newcommand{\Bold}[1]{\mathbf{#1}}\frac{\pi^{2} {t_{\mathrm{grav}}}^{2}}{2 \, {\lambda}^{2} {t_{\mathrm{Higgs}}}^{2}} \verb|=| 5.044311146293211 \times 10^{-32} \qquad \frac{16 \, \pi^{2} {\Omega_{\Lambda}} {\lambda}^{2} {t_{\mathrm{Higgs}}}^{2}}{{g}^{2} {t_{\mathrm{Hubble}}}^{2}} \verb|=| 8.651976825635332 \times 10^{-87}\]

K0 = var('K0',latex_name=r"K_{0}")
formula[K0] = pi/2
alpha0 = var('alpha0',latex_name=r"\alpha_{0}")
formula[alpha0] = (3*pi*k0sq/(2*K_EW))^(1/3)
print_values(4*K0*alpha0*a0*H0,(2*pi*2*lambdaH/g)*(t_Higgs/t_Hubble))

\[\newcommand{\Bold}[1]{\mathbf{#1}}4 \, {H_{0}} {K_{0}} {a_{0}} {\alpha_{0}} \verb|=| 1.1238604496607828 \times 10^{-43} \qquad \frac{4 \, \pi {\lambda} {t_{\mathrm{Higgs}}}}{{g} {t_{\mathrm{Hubble}}}} \verb|=| 1.1238604496607786 \times 10^{-43}\]

6.4 Temperature after $a_{\mathrm{EW}}$

binary_precision=668
Reals = RealField(binary_precision)
RealNumber = Reals
myR = Reals

Kp0 = var('Kp0',latex_name=r"K'_{0}")
value[Kp0]= elliptic_kc(1-myR(valof(k0sq)))
print_values(Kp0)

\[\newcommand{\Bold}[1]{\mathbf{#1}}{K'_{0}} \verb|=| 64.82604415468899\]

T_CGF0 = var('T_CGF0',latex_name=r"T_{\mathrm{CGF},0}")
formula[T_CGF0]=hbar/(kB*4*Kp0*alpha0*a0)
value
print_values(T_CGF0)
pretty_print("\n\n",T_CGF0," = ","%e"%(valof(T_CGF0/K))," K")
print_values(T_CGF0*kB)

\[\newcommand{\Bold}[1]{\mathbf{#1}}{T_{\mathrm{CGF},0}} \verb|=| \frac{{\hbar}}{4 \, {K'_{0}} {a_{0}} {\alpha_{0}} {k_{B}}} \verb|=| 3597163664137.4536 K\]

\[\newcommand{\Bold}[1]{\mathbf{#1}}\begin{array}{l} \phantom{\verb!x!}\\ \phantom{\verb!x!}\\ \phantom{\verb!x!} \end{array} {T_{\mathrm{CGF},0}} \phantom{\verb!x!}\verb|=| \verb|3.597164e+12| \phantom{\verb!x!}\verb|K|\]

\[\newcommand{\Bold}[1]{\mathbf{#1}}{T_{\mathrm{CGF},0}} {k_{B}} \verb|=| 0.3099795809235172 \mathit{GeV}\]

print_values(hbar*g/(kB*8*lambdaH*t_Higgs))

\[\newcommand{\Bold}[1]{\mathbf{#1}}\frac{{g} {\hbar}}{8 \, {k_{B}} {\lambda} {t_{\mathrm{Higgs}}}} \verb|=| 233189890523018.44 K\]

print_values(hbar*g/(8*lambdaH*t_Higgs))

\[\newcommand{\Bold}[1]{\mathbf{#1}}\frac{{g} {\hbar}}{8 \, {\lambda} {t_{\mathrm{Higgs}}}} \verb|=| 20.09475 \mathit{GeV}\]

print_values(hbar*g/(kB*8*lambdaH*a0))

\[\newcommand{\Bold}[1]{\mathbf{#1}}\frac{{g} {\hbar}}{8 \, {a_{0}} {k_{B}} {\lambda}} \verb|=| 8.768880095769793 \times 10^{-05} K\]

6.7 Semiclassical approximation

epsilon = var('epsilon',latex_name=r"\epsilon")
value[epsilon] =1e-27
print_values(2*pi^2*2*K_EW/(2*pi*epsilon^3*g^2))

\[\newcommand{\Bold}[1]{\mathbf{#1}}\frac{2 \, \pi {K_{\mathrm{EW}}}}{{\epsilon}^{3} {g}^{2}} \verb|=| 2.7327966956656648 \times 10^{82}\]