{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Arithmetic and algebra\n",
    "\n",
    "This SageMath notebook does arithmetic and algebra for the paper *A theory of the dark matter*.\n",
    "\n",
    "The section numbering follows the paper.  Equation numbers refer to equations in the paper."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Preamble"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Code for displaying equations."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "%display latex\n",
    "LE = lambda latex_string: LatexExpr(latex_string);"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The numerical values of fundamental constants and basic physical quantities are stored in the dictionary value[].\n",
    "\n",
    "Physical units such as 's' and 'GeV' are defined as algebraic variables.\n",
    "\n",
    "The dictionary formula[] will contain the formulas for derived physical quantities.  For example, formula[kappa] = 8\\*pi\\*G.\n",
    "\n",
    "The function valof(x) substitutes in a formula x to obtain a numerical value.\n",
    "\n",
    "The function print_values(x1,x2,...) prints formulas xn with their numerical values.\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [],
   "source": [
    "value = {}\n",
    "formula={}"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [],
   "source": [
    "def valof(x):\n",
    "    SymRing = FractionField(PolynomialRing(RR,'s,GeV,J,m,kg,K,C'))\n",
    "    NumRing = FractionField(PolynomialRing(RDF,'s,GeV,J,m,kg,K,C'))\n",
    "    xval = SR(x)\n",
    "    for a in reversed(formula):\n",
    "        xval = xval.subs({a:formula[a]})\n",
    "    for key,val in value.items():\n",
    "        xval = xval.subs({key:val})\n",
    "    xval = NumRing(xval)\n",
    "    return xval\n",
    "def print_values(*args):\n",
    "    print_list=[]\n",
    "    for arg in args:\n",
    "        if arg in formula:\n",
    "            print_list += [arg,'=',formula[arg],'=',valof(arg),LE(r'\\qquad')]\n",
    "        else:\n",
    "            print_list += [arg,'=',valof(arg),LE(r'\\qquad')]\n",
    "    print_list.pop()\n",
    "    pretty_print(*print_list)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 1.3 Physical parameters"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Units and fundamental constants as variables"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "declare units as variables"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [],
   "source": [
    "s = var('s', domain='positive',latex_name='\\\\mathrm{s}'); assume(s,'real');\n",
    "GeV = var('GeV', domain='positive',latex_name='\\\\mathrm{GeV}'); assume(GeV,'real');\n",
    "J = var('J', domain='positive',latex_name='\\\\mathrm{J}'); assume(J,'real');\n",
    "m = var('m', domain='positive',latex_name='\\\\mathrm{m}'); assume(m,'real');\n",
    "kg = var('kg', domain='positive',latex_name='\\\\mathrm{kg}\\\\,'); assume(kg,'real');\n",
    "K = var('K', domain='positive',latex_name='\\\\mathrm{K}'); assume(K,'real');\n",
    "C = var('C', domain='positive',latex_name='\\\\mathrm{C}'); assume(C,'real');\n",
    "km = var('km', domain='positive',latex_name='\\\\mathrm{km}\\\\,'); assume(km,'real');\n",
    "Mpc = var('Mpc', domain='positive',latex_name='\\\\mathrm{Mpc}\\\\,'); assume(Mpc,'real');\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "declare fundamental constants and physical parameters as variables"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [],
   "source": [
    "# fundamental constants\n",
    "c = var('c',latex_name=r\"c\")\n",
    "c_mps = var('c_mps',latex_name=r\"c\")\n",
    "e_charge = var('e_charge',latex_name=r\"e\")\n",
    "hbar = var('hbar',latex_name=r\"\\hbar\")\n",
    "kB = var('kB',latex_name=r\"k_{B}\")\n",
    "G = var('G',latex_name=r\"G\")\n",
    "kappa = var('kappa',latex_name=r\"\\kappa\")\n",
    "#\n",
    "# numbered constants c_n for temporary use\n",
    "# c_ = var('c_',n=20,latex_name=r\"c\")\n",
    "#\n",
    "# Standard Model\n",
    "GFermi = var('GFermi',latex_name=r\"G_{\\mathrm{Fermi}}\")\n",
    "m_W = var('m_W',latex_name=r\"m_{W}\")\n",
    "m_Higgs = var('m_Higgs',latex_name=r\"m_{\\mathrm{Higgs}}\")\n",
    "hbar_v = var('hbar_v',latex_name=r\"\\hbar v\")\n",
    "v = var('v',latex_name=r\"v\")\n",
    "g = var('g',latex_name=r\"g\")\n",
    "lambdaH = var('lambdaH',latex_name=r\"\\lambda\")\n",
    "#\n",
    "# cosmological parameters\n",
    "H0 = var('H0',latex_name=r\"H_{0}\")\n",
    "Omega_curvature = var('Omega_curvature',latex_name=r\"\\Omega_{\\mathrm{curvature}}\")\n",
    "Omega_Lambda = var('Omega_Lambda',latex_name=r\"\\Omega_{\\Lambda}\")\n",
    "#\n",
    "# time and energy scales\n",
    "t_grav = var('t_grav',latex_name=r\"t_{\\mathrm{grav}}\")\n",
    "t_Higgs = var('t_Higgs',latex_name=r\"t_{\\mathrm{Higgs}}\")\n",
    "t_Hubble = var('t_Hubble',latex_name=r\"t_{\\mathrm{Hubble}}\")\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Values of the fundamental constants from NIST 2018"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [],
   "source": [
    "value[G]        = 6.67430e-11*m^3*kg^(-1)* s^(-2)\n",
    "value[e_charge] = 1.602176634e-19 * C\n",
    "value[c]=c_mps  = 2.99792458e8 * m * s^(-1)\n",
    "value[hbar]     = 1.054571817e-34*J*s\n",
    "value[kB]       = 1.380649e-23*J*K^(-1)\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<html>\\[\\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}\\]</html>"
      ],
      "text/latex": [
       "$$\\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}$$"
      ],
      "text/plain": [
       "G '=' (6.6743e-11*m^3)/(s^2*kg) \\qquad e_charge '=' 1.602176634e-19*C \\qquad c '=' 299792458.0*m/s"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "\n",
      "\n"
     ]
    },
    {
     "data": {
      "text/html": [
       "<html>\\[\\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}\\]</html>"
      ],
      "text/latex": [
       "$$\\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}$$"
      ],
      "text/plain": [
       "hbar '=' 1.054571817e-34*s*J \\qquad kB '=' (1.380649e-23*J)/K"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "print_values(G,e_charge,c)\n",
    "print('\\n')\n",
    "print_values(hbar,kB)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The constant $\\kappa = 8 \\pi G$."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<html>\\[\\newcommand{\\Bold}[1]{\\mathbf{#1}}{\\kappa} \\verb|=| 8 \\, \\pi {G} \\verb|=| \\frac{1.6774345478283484 \\times 10^{-09} m^{3}}{s^{2} \\mathit{kg}}\\]</html>"
      ],
      "text/latex": [
       "$$\\newcommand{\\Bold}[1]{\\mathbf{#1}}{\\kappa} \\verb|=| 8 \\, \\pi {G} \\verb|=| \\frac{1.6774345478283484 \\times 10^{-09} m^{3}}{s^{2} \\mathit{kg}}$$"
      ],
      "text/plain": [
       "kappa '=' 8*pi*G '=' (1.6774345478283484e-09*m^3)/(s^2*kg)"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "formula[kappa]  = 8*pi*G\n",
    "print_values(kappa)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### convert to c=1 units with unit of energy = GeV"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<html>\\[\\newcommand{\\Bold}[1]{\\mathbf{#1}}{c} \\verb|=| 1.0 \\qquad {\\mathrm{J}} \\verb|=| 6241509074.460764 \\mathit{GeV}\\]</html>"
      ],
      "text/latex": [
       "$$\\newcommand{\\Bold}[1]{\\mathbf{#1}}{c} \\verb|=| 1.0 \\qquad {\\mathrm{J}} \\verb|=| 6241509074.460764 \\mathit{GeV}$$"
      ],
      "text/plain": [
       "c '=' 1.0 \\qquad J '=' 6241509074.460764*GeV"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "\n",
      "\n"
     ]
    },
    {
     "data": {
      "text/html": [
       "<html>\\[\\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}\\]</html>"
      ],
      "text/latex": [
       "$$\\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}$$"
      ],
      "text/plain": [
       "m '=' 3.3356409519815204e-09*s \\qquad kg '=' 5.6095886038044526e+26*GeV"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "\n",
      "\n"
     ]
    },
    {
     "data": {
      "text/html": [
       "<html>\\[\\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}\\]</html>"
      ],
      "text/latex": [
       "$$\\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}$$"
      ],
      "text/plain": [
       "hbar '=' 6.582119565476076e-25*s*GeV \\qquad kB '=' (8.61733326214518e-14*GeV)/K"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "\n",
      "\n"
     ]
    },
    {
     "data": {
      "text/html": [
       "<html>\\[\\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}}\\]</html>"
      ],
      "text/latex": [
       "$$\\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}}$$"
      ],
      "text/plain": [
       "G '=' (4.415832614329417e-63*s)/GeV \\qquad kappa '=' 8*pi*G '=' (1.1098197840527606e-61*s)/GeV"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "conv={}\n",
    "conv[m]  = c^(-1)*m\n",
    "conv[kg] = J*s^2*m^(-2)\n",
    "conv[J]  = e_charge^(-1) * C * 1e-9 * GeV\n",
    "for a,b in conv.items():\n",
    "    value[a]= a.subs(conv).subs(value)\n",
    "#\n",
    "print_values(c,J)\n",
    "print('\\n')\n",
    "print_values(m,kg)\n",
    "print('\\n')\n",
    "print_values(hbar,kB)\n",
    "print('\\n')\n",
    "print_values(G,kappa)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Standard Model coupling constants from PDG (2020, 2021)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "measured quantities"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<html>\\[\\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}\\]</html>"
      ],
      "text/latex": [
       "$$\\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}$$"
      ],
      "text/plain": [
       "GFermi '=' (1.1663787e-05)/GeV^2 \\qquad m_W '=' 80.379*GeV \\qquad m_Higgs '=' 125.1*GeV"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "value[GFermi]  = 1.1663787e-5*GeV^(-2);\n",
    "value[m_W]     = 80.379*GeV;\n",
    "value[m_Higgs] = 125.10*GeV;\n",
    "#\n",
    "print_values(GFermi,m_W,m_Higgs)\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "derived quantities"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<html>\\[\\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}\\]</html>"
      ],
      "text/latex": [
       "$$\\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}$$"
      ],
      "text/plain": [
       "hbar_v '=' 1/2*2^(3/4)/sqrt(GFermi) '=' 246.2196507941374*GeV \\qquad v '=' 1/2*2^(3/4)/(sqrt(GFermi)*hbar) '=' (3.740735007087775e+26)/s"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "\n",
      "\n"
     ]
    },
    {
     "data": {
      "text/html": [
       "<html>\\[\\newcommand{\\Bold}[1]{\\mathbf{#1}}{g} \\verb|=| \\frac{2 \\, {m_{W}}}{{\\hbar v}} \\verb|=| 0.6529048330687819 \\qquad {g}^{2} \\verb|=| 0.4262847210445737\\]</html>"
      ],
      "text/latex": [
       "$$\\newcommand{\\Bold}[1]{\\mathbf{#1}}{g} \\verb|=| \\frac{2 \\, {m_{W}}}{{\\hbar v}} \\verb|=| 0.6529048330687819 \\qquad {g}^{2} \\verb|=| 0.4262847210445737$$"
      ],
      "text/plain": [
       "g '=' 2*m_W/hbar_v '=' 0.6529048330687819 \\qquad g^2 '=' 0.4262847210445737"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "\n",
      "\n"
     ]
    },
    {
     "data": {
      "text/html": [
       "<html>\\[\\newcommand{\\Bold}[1]{\\mathbf{#1}}{\\lambda} \\verb|=| \\frac{{m_{\\mathrm{Higgs}}}}{{\\hbar v}} \\verb|=| 0.5080829235055462 \\qquad {\\lambda}^{2} \\verb|=| 0.25814825715794265\\]</html>"
      ],
      "text/latex": [
       "$$\\newcommand{\\Bold}[1]{\\mathbf{#1}}{\\lambda} \\verb|=| \\frac{{m_{\\mathrm{Higgs}}}}{{\\hbar v}} \\verb|=| 0.5080829235055462 \\qquad {\\lambda}^{2} \\verb|=| 0.25814825715794265$$"
      ],
      "text/plain": [
       "lambdaH '=' m_Higgs/hbar_v '=' 0.5080829235055462 \\qquad lambdaH^2 '=' 0.25814825715794265"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "#\n",
    "formula[hbar_v] = 2^(-1/4)*GFermi^(-1/2)\n",
    "formula[v] = formula[hbar_v]/hbar\n",
    "formula[g] = 2*m_W/hbar_v;\n",
    "formula[lambdaH] = m_Higgs/hbar_v;\n",
    "#\n",
    "print_values(hbar_v,v)\n",
    "print('\\n')\n",
    "print_values(g,g^2)\n",
    "print('\\n')\n",
    "print_values(lambdaH,lambdaH^2)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Cosmological parameters"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "units of distance\n",
    "\n",
    "The value of 1 parsec in meters is taken from <I>IAU 2015 Resolution B2</I>, note 4 \n",
    "which references the definition as exactly 64000/$\\pi$ au and the definition of au from <I>IAU 2012 Resolution B2</I>.  "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<html>\\[\\newcommand{\\Bold}[1]{\\mathbf{#1}}{\\mathrm{km}\\,} \\verb|=| 1000.0 m \\qquad {\\mathrm{Mpc}\\,} \\verb|=| 3.0856775810000003 \\times 10^{22} m\\]</html>"
      ],
      "text/latex": [
       "$$\\newcommand{\\Bold}[1]{\\mathbf{#1}}{\\mathrm{km}\\,} \\verb|=| 1000.0 m \\qquad {\\mathrm{Mpc}\\,} \\verb|=| 3.0856775810000003 \\times 10^{22} m$$"
      ],
      "text/plain": [
       "km '=' 1000.0*m \\qquad Mpc '=' 3.0856775810000003e+22*m"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "value[km]       = 1e3 *m\n",
    "value[Mpc]      = 1e6 * 3.085677581 * 1e16 * m\n",
    "print_values(km,Mpc)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "\n",
    "From <I>Particle Data Group 2020 Particle Physics Booklet</I>"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<html>\\[\\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\\]</html>"
      ],
      "text/latex": [
       "$$\\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$$"
      ],
      "text/plain": [
       "H0 '=' 67.4000000000000*km/(Mpc*s) '=' (2.1842852414333302e-18)/s \\qquad Omega_Lambda '=' 0.685 \\qquad Omega_curvature '=' 0.001"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "formula[H0] = 6.74e1 *km * s^(-1) * Mpc^(-1)\n",
    "value[Omega_curvature] = 0.001\n",
    "value[Omega_Lambda] = 0.685\n",
    "#\n",
    "print_values(H0,Omega_Lambda,Omega_curvature)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### time and energy scales"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<html>\\[\\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}\\]</html>"
      ],
      "text/latex": [
       "$$\\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}$$"
      ],
      "text/plain": [
       "t_grav '=' sqrt(hbar*kappa) '=' 2.7027701557413476e-43*s \\qquad hbar/t_grav '=' 2.435323459338203e+18*GeV"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "\n",
      "\n"
     ]
    },
    {
     "data": {
      "text/html": [
       "<html>\\[\\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}\\]</html>"
      ],
      "text/latex": [
       "$$\\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}$$"
      ],
      "text/plain": [
       "t_Higgs '=' hbar/m_Higgs '=' 5.2614864632102925e-27*s \\qquad m_Higgs '=' 125.1*GeV"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "\n",
      "\n"
     ]
    },
    {
     "data": {
      "text/html": [
       "<html>\\[\\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}\\]</html>"
      ],
      "text/latex": [
       "$$\\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}$$"
      ],
      "text/plain": [
       "t_Hubble '=' 1/H0 '=' 4.5781566483679526e+17*s \\qquad hbar/t_Hubble '=' 1.4377226624218956e-42*GeV"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "formula[t_grav] = sqrt(kappa*hbar)\n",
    "formula[t_Higgs] = hbar/m_Higgs\n",
    "formula[t_Hubble] = 1/H0\n",
    "#\n",
    "print_values(t_grav,hbar/t_grav)\n",
    "print('\\n')\n",
    "print_values(t_Higgs,m_Higgs)\n",
    "print('\\n')\n",
    "print_values(t_Hubble,hbar/t_Hubble)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 2.2 Initial CGF energy"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "elliptic parameter $k^{2}_{\\mathrm{EW}}=\\frac12$ and elliptic integral of first kind $K_{\\mathrm{EW}}=K(k_{\\mathrm{EW}})=K'(k_{\\mathrm{EW}})$"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<html>\\[\\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\\]</html>"
      ],
      "text/latex": [
       "$$\\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$$"
      ],
      "text/plain": [
       "ksq_EW '=' 0.5 \\qquad K_EW '=' 1.8540746773013719 \\qquad Kp_EW '=' 1.8540746773013719"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "# elliptic modulus and integrals\n",
    "ksq_EW = var('ksq_EW',latex_name=r\"k^2_{\\mathrm{EW}}\")\n",
    "k_EW = var('k_EW',latex_name=r\"k_{\\mathrm{EW}}\")\n",
    "K_EW = var('K_EW',latex_name=r\"K_{\\mathrm{EW}}\")\n",
    "Kp_EW = var('Kp_EW',latex_name=r\"K'_{\\mathrm{EW}}\")\n",
    "#\n",
    "value[ksq_EW]= 1/2\n",
    "value[k_EW]= sqrt(1/2)\n",
    "value[K_EW]= elliptic_kc(0.5)\n",
    "value[Kp_EW]= elliptic_kc(0.5)\n",
    "print_values(ksq_EW,K_EW,Kp_EW)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 2.4 Start of the electroweak transition at $a=a_{\\mathrm{EW}}$"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 16,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<html>\\[\\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\\]</html>"
      ],
      "text/latex": [
       "$$\\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$$"
      ],
      "text/plain": [
       "1/4*sqrt(6)*sqrt(pi)/K_EW '=' 0.5854143283037644 \\qquad a_EW '=' 1/4*sqrt(6)*sqrt(pi)*t_Higgs/K_EW '=' 3.0801495637396022e-27*s"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "a_EW = var('a_EW',latex_name=r\"a_{\\mathrm{EW}}\")\n",
    "#\n",
    "a_EW_coeff = (6*pi)^(1/2)/(4*K_EW)\n",
    "formula[a_EW] = a_EW_coeff * t_Higgs\n",
    "print_values(a_EW_coeff,a_EW)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 17,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<html>\\[\\newcommand{\\Bold}[1]{\\mathbf{#1}}{\\hat a_{\\mathrm{EW}}} \\verb|=| \\frac{{a_{\\mathrm{EW}}}}{{t_{\\mathrm{Higgs}}}} \\verb|=| 0.5854143283037644\\]</html>"
      ],
      "text/latex": [
       "$$\\newcommand{\\Bold}[1]{\\mathbf{#1}}{\\hat a_{\\mathrm{EW}}} \\verb|=| \\frac{{a_{\\mathrm{EW}}}}{{t_{\\mathrm{Higgs}}}} \\verb|=| 0.5854143283037644$$"
      ],
      "text/plain": [
       "ahat_EW '=' a_EW/t_Higgs '=' 0.5854143283037644"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "ahat_EW = var('ahat_EW',latex_name=r\"\\hat a_{\\mathrm{EW}}\")\n",
    "formula[ahat_EW] = a_EW/t_Higgs\n",
    "print_values(ahat_EW)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 2.6 Realizing the electroweak transition"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 18,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<html>\\[\\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}\\]</html>"
      ],
      "text/latex": [
       "$$\\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}$$"
      ],
      "text/plain": [
       "asqHsq_EW '=' 1/24*t_grav^2*(a_EW^2/(lambdaH^2*t_Higgs^2) + 3*t_Higgs^2/(a_EW^2*g^2))/t_Higgs^2 '=' 2.4037614729587787e-33"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "\n",
      "\n"
     ]
    },
    {
     "data": {
      "text/html": [
       "<html>\\[\\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}\\]</html>"
      ],
      "text/latex": [
       "$$\\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}$$"
      ],
      "text/plain": [
       "1/12*sqrt(3)*t_grav^2/(g*lambdaH*t_Higgs^2) '=' 1.1481436122987478e-33"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "asqHsq_EW = var('asqHsq_EW',\\\n",
    "    latex_name=r\"a_{\\mathrm{EW}}^2 H_{\\mathrm{EW}}^2+\\epsilon^2\")\n",
    "formula[asqHsq_EW ] = (t_grav^2/t_Higgs^2) * ( (1/(24*lambdaH^2))\\\n",
    "    * (a_EW^2/t_Higgs^2) +(1/(8*g^2))* (t_Higgs^2/a_EW^2) )\n",
    "print_values(asqHsq_EW)\n",
    "print('\\n')\n",
    "expr = (t_grav^2/t_Higgs^2)/(4*sqrt(3)*g*lambdaH)\n",
    "print_values(expr)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 4.5 Parametrize the time evolution by $k^2$"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 19,
   "metadata": {},
   "outputs": [],
   "source": [
    "PSR.<k> = PowerSeriesRing(SR)\n",
    "K_ps = (pi/2)*(1+k^2/4 + 9*k^4/64)+O(k^6)\n",
    "E_ps = (pi/2)*(1-k^2/4 - 3*k^4/64)+O(k^6)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 20,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<html>\\[\\newcommand{\\Bold}[1]{\\mathbf{#1}}\\alpha^3 =  \\frac{3 \\, \\pi}{2 \\, {K_{\\mathrm{EW}}}} k^{2} + O(k^{4})\\]</html>"
      ],
      "text/latex": [
       "$$\\newcommand{\\Bold}[1]{\\mathbf{#1}}\\alpha^3 =  \\frac{3 \\, \\pi}{2 \\, {K_{\\mathrm{EW}}}} k^{2} + O(k^{4})$$"
      ],
      "text/plain": [
       "\\alpha^3 =  3/2*pi/K_EW*k^2 + O(k^4)"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "alpha3 = (2*(1-k^2)*K_ps + 2*(2*k^2-1)*E_ps)/K_EW\n",
    "pretty_print(LE(r\"\\alpha^3 = \"),alpha3+ O(k^4))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 21,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<html>\\[\\newcommand{\\Bold}[1]{\\mathbf{#1}}\\alpha^2 \\langle b^2\\rangle =  \\frac{1}{2} k^{2} + O(k^{4})\\]</html>"
      ],
      "text/latex": [
       "$$\\newcommand{\\Bold}[1]{\\mathbf{#1}}\\alpha^2 \\langle b^2\\rangle =  \\frac{1}{2} k^{2} + O(k^{4})$$"
      ],
      "text/plain": [
       "\\alpha^2 \\langle b^2\\rangle =  1/2*k^2 + O(k^4)"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "alpha2_b2 = k^2 -1 +E_ps/K_ps; pretty_print(LE(r\"\\alpha^2 \\langle b^2\\rangle = \"),alpha2_b2+O(k^4))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 22,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<html>\\[\\newcommand{\\Bold}[1]{\\mathbf{#1}}\\alpha^2 \\mu^2 =  1 + O(k^{2})\\]</html>"
      ],
      "text/latex": [
       "$$\\newcommand{\\Bold}[1]{\\mathbf{#1}}\\alpha^2 \\mu^2 =  1 + O(k^{2})$$"
      ],
      "text/plain": [
       "\\alpha^2 \\mu^2 =  1 + O(k^2)"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "alpha2_mu2 = 1-2*k^2; pretty_print(LE(r\"\\alpha^2 \\mu^2 = \"),alpha2_mu2+O(k^2))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 23,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<html>\\[\\newcommand{\\Bold}[1]{\\mathbf{#1}}\\alpha^4 E_{CGF} =  \\frac{1}{2} k^{2} + O(k^{4})\\]</html>"
      ],
      "text/latex": [
       "$$\\newcommand{\\Bold}[1]{\\mathbf{#1}}\\alpha^4 E_{CGF} =  \\frac{1}{2} k^{2} + O(k^{4})$$"
      ],
      "text/plain": [
       "\\alpha^4 E_{CGF} =  1/2*k^2 + O(k^4)"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "alpha4_E_CGF = k^2*(1-k^2)/2; pretty_print(LE(r\"\\alpha^4 E_{CGF} = \"),alpha4_E_CGF+O(k^4))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 24,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<html>\\[\\newcommand{\\Bold}[1]{\\mathbf{#1}}\\alpha^2 \\hat a^2=  \\frac{4 \\, {\\lambda}^{2}}{{g}^{2}} + O(k^{2})\\]</html>"
      ],
      "text/latex": [
       "$$\\newcommand{\\Bold}[1]{\\mathbf{#1}}\\alpha^2 \\hat a^2=  \\frac{4 \\, {\\lambda}^{2}}{{g}^{2}} + O(k^{2})$$"
      ],
      "text/plain": [
       "\\alpha^2 \\hat a^2=  4*lambdaH^2/g^2 + O(k^2)"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "alpha2_ahat2 = (3/2)* alpha2_b2 + (4*lambdaH^2/g^2)* alpha2_mu2\n",
    "pretty_print(LE(r\"\\alpha^2 \\hat a^2= \"),alpha2_ahat2+O(k^2))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 25,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<html>\\[\\newcommand{\\Bold}[1]{\\mathbf{#1}}\\hat \\rho_{CGF} =  \\frac{3 \\, {g}^{2}}{32 \\, {\\lambda}^{4}} k^{2} + O(k^{4})\\]</html>"
      ],
      "text/latex": [
       "$$\\newcommand{\\Bold}[1]{\\mathbf{#1}}\\hat \\rho_{CGF} =  \\frac{3 \\, {g}^{2}}{32 \\, {\\lambda}^{4}} k^{2} + O(k^{4})$$"
      ],
      "text/plain": [
       "\\hat \\rho_{CGF} =  3/32*g^2/lambdaH^4*k^2 + O(k^4)"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "rhohat_CGF = (3*alpha4_E_CGF/g^2 +9*alpha2_b2^2/(32*lambdaH^2))/alpha2_ahat2^2\n",
    "pretty_print(LE(r\"\\hat \\rho_{CGF} = \"),rhohat_CGF +O(k^4))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 26,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<html>\\[\\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})\\]</html>"
      ],
      "text/latex": [
       "$$\\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})$$"
      ],
      "text/plain": [
       "\\hat p_{CGF} =  (9/2048*g^4*(8/g^2 - 1/lambdaH^2)/lambdaH^4)*k^4 + O(k^6)"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "phat_CGF = ((alpha4_E_CGF-alpha2_mu2*alpha2_b2)/g^2 -9*alpha2_b2^2/(32*lambdaH^2))/alpha2_ahat2^2\n",
    "pretty_print(LE(r\"\\hat p_{CGF} = \"),phat_CGF +O(k^6))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 5.1 $\\Omega_{\\mathrm{CGF}} + \\Omega_{\\Lambda}=1$"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 27,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<html>\\[\\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}\\]</html>"
      ],
      "text/latex": [
       "$$\\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}$$"
      ],
      "text/plain": [
       "k0sq = 0.315 32*lambdaH^4*t_Higgs^4/(g^2*t_Hubble^2*t_grav^2) '=' 7.887392672844198e-56"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "k0sq = var('k0sq',latex_name=r\"k_{0}^{2}\")\n",
    "expr = t_Higgs^4/(t_grav^2*t_Hubble^2)*(32*lambdaH^4/g^2)\n",
    "formula[k0sq] = 0.315*expr\n",
    "pretty_print(k0sq,LE('= 0.315'),expr,'=',valof(k0sq))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 28,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<html>\\[\\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\\]</html>"
      ],
      "text/latex": [
       "$$\\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$$"
      ],
      "text/plain": [
       "a0 = 0.315^{-\\frac13} 1/6*6^(2/3)/(pi*g*lambdaH*t_Higgs/(K_EW*t_Hubble^2*t_grav^2))^(1/3) '=' 1.3991814678093661e-08*s"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "a0 = var('a0',latex_name=r\"a_{0}\")\n",
    "expr= ((6*pi*lambdaH*g/K_EW)*t_Higgs/(t_grav^2*t_Hubble^2))^(-1/3)\n",
    "formula[a0] = 0.315^(-1/3)*expr\n",
    "pretty_print(a0,LE(r\"= 0.315^{-\\frac13}\"),expr,'=',valof(a0))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 29,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<html>\\[\\newcommand{\\Bold}[1]{\\mathbf{#1}}\\frac{{a_{0}}}{{t_{\\mathrm{Higgs}}}} \\verb|=| 2.6592893046343716 \\times 10^{18}\\]</html>"
      ],
      "text/latex": [
       "$$\\newcommand{\\Bold}[1]{\\mathbf{#1}}\\frac{{a_{0}}}{{t_{\\mathrm{Higgs}}}} \\verb|=| 2.6592893046343716 \\times 10^{18}$$"
      ],
      "text/plain": [
       "a0/t_Higgs '=' 2.6592893046343716e+18"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "print_values(a0/t_Higgs)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 30,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<html>\\[\\newcommand{\\Bold}[1]{\\mathbf{#1}}\\frac{{a_{0}}}{{a_{\\mathrm{EW}}}} \\verb|=| 4.542576387461597 \\times 10^{18}\\]</html>"
      ],
      "text/latex": [
       "$$\\newcommand{\\Bold}[1]{\\mathbf{#1}}\\frac{{a_{0}}}{{a_{\\mathrm{EW}}}} \\verb|=| 4.542576387461597 \\times 10^{18}$$"
      ],
      "text/plain": [
       "a0/a_EW '=' 4.542576387461597e+18"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "print_values(a0/a_EW)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 31,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<html>\\[\\newcommand{\\Bold}[1]{\\mathbf{#1}}\\frac{{t_{\\mathrm{Hubble}}}^{2}}{{a_{0}}^{2}} \\verb|=| 1.0706147161725437 \\times 10^{51}\\]</html>"
      ],
      "text/latex": [
       "$$\\newcommand{\\Bold}[1]{\\mathbf{#1}}\\frac{{t_{\\mathrm{Hubble}}}^{2}}{{a_{0}}^{2}} \\verb|=| 1.0706147161725437 \\times 10^{51}$$"
      ],
      "text/plain": [
       "t_Hubble^2/a0^2 '=' 1.0706147161725437e+51"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "print_values(t_Hubble^2/a0^2)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 5.2 $w_{\\mathrm{CGF}} =0$"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 32,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<html>\\[\\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})\\]</html>"
      ],
      "text/latex": [
       "$$\\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})$$"
      ],
      "text/plain": [
       "w_{CGF} =  (3/64*g^2*(8/g^2 - 1/lambdaH^2))*k^2 + O(k^4)"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "w_CGF = phat_CGF/rhohat_CGF\n",
    "pretty_print(LE(r\"w_{CGF} = \"),w_CGF +O(k^4))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 33,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<html>\\[\\newcommand{\\Bold}[1]{\\mathbf{#1}}\\frac{{K_{\\mathrm{EW}}}}{2 \\, \\pi {g} {\\lambda}} \\verb|=| 0.8895346543935904\\]</html>"
      ],
      "text/latex": [
       "$$\\newcommand{\\Bold}[1]{\\mathbf{#1}}\\frac{{K_{\\mathrm{EW}}}}{2 \\, \\pi {g} {\\lambda}} \\verb|=| 0.8895346543935904$$"
      ],
      "text/plain": [
       "1/2*K_EW/(pi*g*lambdaH) '=' 0.8895346543935904"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "print_values(K_EW/(2*pi*g*lambdaH))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 5.3 The CGF in the present"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 34,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<html>\\[\\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}\\]</html>"
      ],
      "text/latex": [
       "$$\\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}$$"
      ],
      "text/plain": [
       "omega_W '=' 1/2*g/(lambdaH*t_Higgs) '=' (1.221171982678596e+26)/s \\qquad m_W/hbar '=' (1.221171982678596e+26)/s"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "omega_W = var('omega_W',latex_name=r\"\\omega_W\")\n",
    "formula[omega_W] = g/(2*lambdaH*t_Higgs)\n",
    "print_values(omega_W,m_W/hbar)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 35,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<html>\\[\\newcommand{\\Bold}[1]{\\mathbf{#1}}\\frac{2 \\, \\pi}{{\\omega_W}} \\verb|=| 5.145209189452289 \\times 10^{-26} s\\]</html>"
      ],
      "text/latex": [
       "$$\\newcommand{\\Bold}[1]{\\mathbf{#1}}\\frac{2 \\, \\pi}{{\\omega_W}} \\verb|=| 5.145209189452289 \\times 10^{-26} s$$"
      ],
      "text/plain": [
       "2*pi/omega_W '=' 5.145209189452289e-26*s"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "print_values(2*pi/omega_W)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 36,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<html>\\[\\newcommand{\\Bold}[1]{\\mathbf{#1}}{k_{0}} \\verb|=| \\sqrt{{k_{0}^{2}}} \\verb|=| 2.8084502261646374 \\times 10^{-28}\\]</html>"
      ],
      "text/latex": [
       "$$\\newcommand{\\Bold}[1]{\\mathbf{#1}}{k_{0}} \\verb|=| \\sqrt{{k_{0}^{2}}} \\verb|=| 2.8084502261646374 \\times 10^{-28}$$"
      ],
      "text/plain": [
       "k0 '=' sqrt(k0sq) '=' 2.8084502261646374e-28"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "\n",
      "\n"
     ]
    },
    {
     "data": {
      "text/html": [
       "<html>\\[\\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}\\]</html>"
      ],
      "text/latex": [
       "$$\\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}$$"
      ],
      "text/plain": [
       "k0*omega_W '=' 0.03429600730939622/s \\qquad k0*m_W '=' 2.257404207288874e-26*GeV"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "k0 = var('k0',latex_name=r\"k_{0}\")\n",
    "formula[k0]= sqrt(k0sq)\n",
    "print_values(k0)\n",
    "print('\\n')\n",
    "print_values(k0*omega_W,k0*m_W)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 5.4 CGF equation of state"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 37,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<html>\\[\\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}}\\]</html>"
      ],
      "text/latex": [
       "$$\\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}}$$"
      ],
      "text/plain": [
       "rho_b ' = ' hbar/t_Higgs^4 ' = ' (5.682491786274889e+28)*kg/m^3"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "rho_b = var('rho_b ',latex_name=r\"\\rho_{b}\")\n",
    "formula[rho_b] = hbar/t_Higgs^4\n",
    "pretty_print(rho_b,\" = \",formula[rho_b] ,\" = \", kg*valof(rho_b/kg)/c_mps^3)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 38,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<html>\\[\\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\\]</html>"
      ],
      "text/latex": [
       "$$\\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$$"
      ],
      "text/plain": [
       "rhohat_EW '=' 8/3*K_EW^4/(pi^2*g^2) + 1/8/lambdaH^2 '=' 7.97415392287409"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "rhohat_EW = var('rhohat_EW ',latex_name=r\"\\hat\\rho_{EW}\")\n",
    "formula[rhohat_EW ] = 8 *K_EW^4/(3*pi^2*g^2) + 1/(8*lambdaH^2)\n",
    "print_values(rhohat_EW)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 39,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<html>\\[\\newcommand{\\Bold}[1]{\\mathbf{#1}}c_{b} \\verb|=| \\frac{1}{2 \\, {\\lambda}^{2} {\\hat\\rho_{EW}}} \\verb|=| 0.24289366755876735\\]</html>"
      ],
      "text/latex": [
       "$$\\newcommand{\\Bold}[1]{\\mathbf{#1}}c_{b} \\verb|=| \\frac{1}{2 \\, {\\lambda}^{2} {\\hat\\rho_{EW}}} \\verb|=| 0.24289366755876735$$"
      ],
      "text/plain": [
       "c_b '=' 1/2/(lambdaH^2*rhohat_EW) '=' 0.24289366755876735"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "c_b = var('c_b')\n",
    "formula[c_b] = 1/(2*lambdaH^2*rhohat_EW)\n",
    "print_values(c_b)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 40,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<html>\\[\\newcommand{\\Bold}[1]{\\mathbf{#1}}c_{a} \\verb|=| -\\frac{{\\left({g}^{2} - 8 \\, {\\lambda}^{2}\\right)} {\\lambda}^{2}}{{g}^{2}} \\verb|=| 0.9924810876684255\\]</html>"
      ],
      "text/latex": [
       "$$\\newcommand{\\Bold}[1]{\\mathbf{#1}}c_{a} \\verb|=| -\\frac{{\\left({g}^{2} - 8 \\, {\\lambda}^{2}\\right)} {\\lambda}^{2}}{{g}^{2}} \\verb|=| 0.9924810876684255$$"
      ],
      "text/plain": [
       "c_a '=' -(g^2 - 8*lambdaH^2)*lambdaH^2/g^2 '=' 0.9924810876684255"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "0.9924810876684255\n"
     ]
    }
   ],
   "source": [
    "c_a = var('c_a')\n",
    "formula[c_a] = 2*lambdaH^2*(8*lambdaH^2-g^2)/(2*g^2)\n",
    "c_a_val = valof(c_a)\n",
    "print_values(c_a)\n",
    "print(c_a_val)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 5.5 Adiabatic condition for  $a\\ge a_{\\mathrm{EW}}$"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 41,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<html>\\[\\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}\\]</html>"
      ],
      "text/latex": [
       "$$\\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}$$"
      ],
      "text/plain": [
       "y '=' 4*K_EW*ahat_EW*sqrt(Omega_Lambda*t_Higgs^2/t_Hubble^2 + 1/3*rhohat_EW*t_grav^2/t_Higgs^2) '=' 3.6360755535373564e-16"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "y=var('y',latex_name=r\"4K_{\\mathrm{EW}}(a^2H^2+\\epsilon^2)^{1/2}_{\\mathrm{EW}}\")\n",
    "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))\n",
    "print_values(y)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 42,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<html>\\[\\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}\\]</html>"
      ],
      "text/latex": [
       "$$\\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}$$"
      ],
      "text/plain": [
       "1/2*pi^2*t_grav^2/(lambdaH^2*t_Higgs^2) '=' 5.044311146293211e-32 \\qquad 16*pi^2*Omega_Lambda*lambdaH^2*t_Higgs^2/(g^2*t_Hubble^2) '=' 8.651976825635332e-87"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "print_values((pi^2/(2*lambdaH^2))*(t_grav^2/t_Higgs^2),\\\n",
    "    (16*pi^2*lambdaH^2/g^2)*Omega_Lambda*(t_Higgs^2/t_Hubble^2))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 43,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<html>\\[\\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}\\]</html>"
      ],
      "text/latex": [
       "$$\\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}$$"
      ],
      "text/plain": [
       "4*H0*K0*a0*alpha0 '=' 1.1238604496607828e-43 \\qquad 4*pi*lambdaH*t_Higgs/(g*t_Hubble) '=' 1.1238604496607786e-43"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "K0 = var('K0',latex_name=r\"K_{0}\")\n",
    "formula[K0] = pi/2\n",
    "alpha0 = var('alpha0',latex_name=r\"\\alpha_{0}\")\n",
    "formula[alpha0] = (3*pi*k0sq/(2*K_EW))^(1/3)\n",
    "print_values(4*K0*alpha0*a0*H0,(2*pi*2*lambdaH/g)*(t_Higgs/t_Hubble))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 6.4 Temperature after $a_{\\mathrm{EW}}$"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 44,
   "metadata": {},
   "outputs": [],
   "source": [
    "binary_precision=668\n",
    "Reals = RealField(binary_precision)\n",
    "RealNumber = Reals\n",
    "myR = Reals"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 45,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<html>\\[\\newcommand{\\Bold}[1]{\\mathbf{#1}}{K'_{0}} \\verb|=| 64.82604415468899\\]</html>"
      ],
      "text/latex": [
       "$$\\newcommand{\\Bold}[1]{\\mathbf{#1}}{K'_{0}} \\verb|=| 64.82604415468899$$"
      ],
      "text/plain": [
       "Kp0 '=' 64.82604415468899"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "Kp0 = var('Kp0',latex_name=r\"K'_{0}\")\n",
    "value[Kp0]= elliptic_kc(1-myR(valof(k0sq)))\n",
    "print_values(Kp0)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 46,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<html>\\[\\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\\]</html>"
      ],
      "text/latex": [
       "$$\\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$$"
      ],
      "text/plain": [
       "T_CGF0 '=' 1/4*hbar/(Kp0*a0*alpha0*kB) '=' 3597163664137.4536*K"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "data": {
      "text/html": [
       "<html>\\[\\newcommand{\\Bold}[1]{\\mathbf{#1}}\\begin{array}{l}\n",
       "\\phantom{\\verb!x!}\\\\\n",
       "\\phantom{\\verb!x!}\\\\\n",
       "\\phantom{\\verb!x!}\n",
       "\\end{array} {T_{\\mathrm{CGF},0}} \\phantom{\\verb!x!}\\verb|=| \\verb|3.597164e+12| \\phantom{\\verb!x!}\\verb|K|\\]</html>"
      ],
      "text/latex": [
       "$$\\newcommand{\\Bold}[1]{\\mathbf{#1}}\\begin{array}{l}\n",
       "\\phantom{\\verb!x!}\\\\\n",
       "\\phantom{\\verb!x!}\\\\\n",
       "\\phantom{\\verb!x!}\n",
       "\\end{array} {T_{\\mathrm{CGF},0}} \\phantom{\\verb!x!}\\verb|=| \\verb|3.597164e+12| \\phantom{\\verb!x!}\\verb|K|$$"
      ],
      "text/plain": [
       "'\\n\\n' T_CGF0 ' = ' '3.597164e+12' ' K'"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "data": {
      "text/html": [
       "<html>\\[\\newcommand{\\Bold}[1]{\\mathbf{#1}}{T_{\\mathrm{CGF},0}} {k_{B}} \\verb|=| 0.3099795809235172 \\mathit{GeV}\\]</html>"
      ],
      "text/latex": [
       "$$\\newcommand{\\Bold}[1]{\\mathbf{#1}}{T_{\\mathrm{CGF},0}} {k_{B}} \\verb|=| 0.3099795809235172 \\mathit{GeV}$$"
      ],
      "text/plain": [
       "T_CGF0*kB '=' 0.3099795809235172*GeV"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "T_CGF0 = var('T_CGF0',latex_name=r\"T_{\\mathrm{CGF},0}\")\n",
    "formula[T_CGF0]=hbar/(kB*4*Kp0*alpha0*a0)\n",
    "value\n",
    "print_values(T_CGF0)\n",
    "pretty_print(\"\\n\\n\",T_CGF0,\" = \",\"%e\"%(valof(T_CGF0/K)),\" K\")\n",
    "print_values(T_CGF0*kB)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 47,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<html>\\[\\newcommand{\\Bold}[1]{\\mathbf{#1}}\\frac{{g} {\\hbar}}{8 \\, {k_{B}} {\\lambda} {t_{\\mathrm{Higgs}}}} \\verb|=| 233189890523018.44 K\\]</html>"
      ],
      "text/latex": [
       "$$\\newcommand{\\Bold}[1]{\\mathbf{#1}}\\frac{{g} {\\hbar}}{8 \\, {k_{B}} {\\lambda} {t_{\\mathrm{Higgs}}}} \\verb|=| 233189890523018.44 K$$"
      ],
      "text/plain": [
       "1/8*g*hbar/(kB*lambdaH*t_Higgs) '=' 233189890523018.44*K"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "print_values(hbar*g/(kB*8*lambdaH*t_Higgs))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 48,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<html>\\[\\newcommand{\\Bold}[1]{\\mathbf{#1}}\\frac{{g} {\\hbar}}{8 \\, {\\lambda} {t_{\\mathrm{Higgs}}}} \\verb|=| 20.09475 \\mathit{GeV}\\]</html>"
      ],
      "text/latex": [
       "$$\\newcommand{\\Bold}[1]{\\mathbf{#1}}\\frac{{g} {\\hbar}}{8 \\, {\\lambda} {t_{\\mathrm{Higgs}}}} \\verb|=| 20.09475 \\mathit{GeV}$$"
      ],
      "text/plain": [
       "1/8*g*hbar/(lambdaH*t_Higgs) '=' 20.09475*GeV"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "print_values(hbar*g/(8*lambdaH*t_Higgs))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 49,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<html>\\[\\newcommand{\\Bold}[1]{\\mathbf{#1}}\\frac{{g} {\\hbar}}{8 \\, {a_{0}} {k_{B}} {\\lambda}} \\verb|=| 8.768880095769793 \\times 10^{-05} K\\]</html>"
      ],
      "text/latex": [
       "$$\\newcommand{\\Bold}[1]{\\mathbf{#1}}\\frac{{g} {\\hbar}}{8 \\, {a_{0}} {k_{B}} {\\lambda}} \\verb|=| 8.768880095769793 \\times 10^{-05} K$$"
      ],
      "text/plain": [
       "1/8*g*hbar/(a0*kB*lambdaH) '=' 8.768880095769793e-05*K"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "print_values(hbar*g/(kB*8*lambdaH*a0))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 6.7 Semiclassical approximation"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 50,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<html>\\[\\newcommand{\\Bold}[1]{\\mathbf{#1}}\\frac{2 \\, \\pi {K_{\\mathrm{EW}}}}{{\\epsilon}^{3} {g}^{2}} \\verb|=| 2.7327966956656648 \\times 10^{82}\\]</html>"
      ],
      "text/latex": [
       "$$\\newcommand{\\Bold}[1]{\\mathbf{#1}}\\frac{2 \\, \\pi {K_{\\mathrm{EW}}}}{{\\epsilon}^{3} {g}^{2}} \\verb|=| 2.7327966956656648 \\times 10^{82}$$"
      ],
      "text/plain": [
       "2*pi*K_EW/(epsilon^3*g^2) '=' 2.7327966956656648e+82"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "epsilon = var('epsilon',latex_name=r\"\\epsilon\")\n",
    "value[epsilon] =1e-27\n",
    "print_values(2*pi^2*2*K_EW/(2*pi*epsilon^3*g^2))"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "SageMath 9.4",
   "language": "sage",
   "name": "sagemath"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.9.5"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}