{ "cells": [ { "cell_type": "markdown", "metadata": {}, "source": [ "# Friedmann equations\n", "\n", "This notebook demonstrates a few capabilities of SageMath in computations regarding cosmological spacetimes with Friedmann-Lemaître-Robertson-Walker (FLRW) metrics.\n", "The corresponding tools have been developed within the [SageManifolds](https://sagemanifolds.obspm.fr) project (version 1.3, as included in SageMath 8.3).\n", "\n", "Click [here](https://raw.githubusercontent.com/sagemanifolds/SageManifolds/master/Worksheets/v1.3/SM_Friedmann_equations.ipynb) to download the notebook file (ipynb format). To run it, you must start SageMath within the Jupyter notebook, via the command `sage -n jupyter`" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "*NB:* a version of SageMath at least equal to 8.2 is required to run this notebook:" ] }, { "cell_type": "code", "execution_count": 1, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "'SageMath version 8.3, Release Date: 2018-08-03'" ] }, "execution_count": 1, "metadata": {}, "output_type": "execute_result" } ], "source": [ "version()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "First we set up the notebook to display mathematical objects using LaTeX formatting:" ] }, { "cell_type": "code", "execution_count": 2, "metadata": {}, "outputs": [], "source": [ "%display latex" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "We declare the spacetime M as a 4-dimensional Lorentzian manifold:" ] }, { "cell_type": "code", "execution_count": 3, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "4-dimensional Lorentzian manifold M\n" ] } ], "source": [ "M = Manifold(4, 'M', structure='Lorentzian')\n", "print(M)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "We introduce the standard FLRW coordinates, via the method `chart()`, the argument of which is a string expressing the coordinates names, their ranges (the default is $(-\\infty,+\\infty)$) and their LaTeX symbols:" ] }, { "cell_type": "code", "execution_count": 4, "metadata": {}, "outputs": [ { "data": { "text/html": [ " $\\newcommand{\\Bold}[1]{\\mathbf{#1}}\\left(M,(t, r, {\\theta}, {\\phi})\\right)$ " ], "text/plain": [ "Chart (M, (t, r, th, ph))" ] }, "execution_count": 4, "metadata": {}, "output_type": "execute_result" } ], "source": [ "fr. = M.chart(r't r:[0,+oo) th:[0,pi]:\\theta ph:[0,2*pi):\\phi')\n", "fr" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "

Assuming that the speed of light c=1, let us define a few variables: Newton's constant $G$, the cosmological constant $\\Lambda$, the spatial curvature constant $k$, the scale factor $a(t)$, the fluid proper density $\\rho(t)$ and the fluid pressure $p(t)$:

" ] }, { "cell_type": "code", "execution_count": 5, "metadata": {}, "outputs": [], "source": [ "var('G, Lambda, k', domain='real')\n", "a = M.scalar_field(function('a')(t), name='a')\n", "rho = M.scalar_field(function('rho')(t), name='rho')\n", "p = M.scalar_field(function('p')(t), name='p')" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "The FLRW metric is defined by its components in the manifold's default frame, i.e. the frame associated with the FLRW coordinates:" ] }, { "cell_type": "code", "execution_count": 6, "metadata": {}, "outputs": [ { "data": { "text/html": [ " $\\newcommand{\\Bold}[1]{\\mathbf{#1}}g = -\\mathrm{d} t\\otimes \\mathrm{d} t + \\left( -\\frac{a\\left(t\\right)^{2}}{k r^{2} - 1} \\right) \\mathrm{d} r\\otimes \\mathrm{d} r + r^{2} a\\left(t\\right)^{2} \\mathrm{d} {\\theta}\\otimes \\mathrm{d} {\\theta} + r^{2} a\\left(t\\right)^{2} \\sin\\left({\\theta}\\right)^{2} \\mathrm{d} {\\phi}\\otimes \\mathrm{d} {\\phi}$ " ], "text/plain": [ "g = -dt*dt - a(t)^2/(k*r^2 - 1) dr*dr + r^2*a(t)^2 dth*dth + r^2*a(t)^2*sin(th)^2 dph*dph" ] }, "execution_count": 6, "metadata": {}, "output_type": "execute_result" } ], "source": [ "g = M.metric()\n", "g[0,0] = -1\n", "g[1,1] = a*a/(1 - k*r^2)\n", "g[2,2] = a*a*r^2\n", "g[3,3] = a*a*(r*sin(th))^2\n", "g.display()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "

A matrix view of the metric components:

" ] }, { "cell_type": "code", "execution_count": 7, "metadata": {}, "outputs": [ { "data": { "text/html": [ " $\\newcommand{\\Bold}[1]{\\mathbf{#1}}\\left(\\begin{array}{rrrr}\n", "-1 & 0 & 0 & 0 \\\\\n", "0 & -\\frac{a\\left(t\\right)^{2}}{k r^{2} - 1} & 0 & 0 \\\\\n", "0 & 0 & r^{2} a\\left(t\\right)^{2} & 0 \\\\\n", "0 & 0 & 0 & r^{2} a\\left(t\\right)^{2} \\sin\\left({\\theta}\\right)^{2}\n", "\\end{array}\\right)$ " ], "text/plain": [ "[ -1 0 0 0]\n", "[ 0 -a(t)^2/(k*r^2 - 1) 0 0]\n", "[ 0 0 r^2*a(t)^2 0]\n", "[ 0 0 0 r^2*a(t)^2*sin(th)^2]" ] }, "execution_count": 7, "metadata": {}, "output_type": "execute_result" } ], "source": [ "g[:]" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "

The Levi-Civita connection associated with the metric is computed:

" ] }, { "cell_type": "code", "execution_count": 8, "metadata": {}, "outputs": [ { "data": { "text/html": [ " $\\newcommand{\\Bold}[1]{\\mathbf{#1}}\\begin{array}{lcl} \\Gamma_{ \\phantom{\\, t} \\, r \\, r }^{ \\, t \\phantom{\\, r} \\phantom{\\, r} } & = & -\\frac{a\\left(t\\right) \\frac{\\partial\\,a}{\\partial t}}{k r^{2} - 1} \\\\ \\Gamma_{ \\phantom{\\, t} \\, {\\theta} \\, {\\theta} }^{ \\, t \\phantom{\\, {\\theta}} \\phantom{\\, {\\theta}} } & = & r^{2} a\\left(t\\right) \\frac{\\partial\\,a}{\\partial t} \\\\ \\Gamma_{ \\phantom{\\, t} \\, {\\phi} \\, {\\phi} }^{ \\, t \\phantom{\\, {\\phi}} \\phantom{\\, {\\phi}} } & = & r^{2} a\\left(t\\right) \\sin\\left({\\theta}\\right)^{2} \\frac{\\partial\\,a}{\\partial t} \\\\ \\Gamma_{ \\phantom{\\, r} \\, t \\, r }^{ \\, r \\phantom{\\, t} \\phantom{\\, r} } & = & \\frac{\\frac{\\partial\\,a}{\\partial t}}{a\\left(t\\right)} \\\\ \\Gamma_{ \\phantom{\\, r} \\, r \\, r }^{ \\, r \\phantom{\\, r} \\phantom{\\, r} } & = & -\\frac{k r}{k r^{2} - 1} \\\\ \\Gamma_{ \\phantom{\\, r} \\, {\\theta} \\, {\\theta} }^{ \\, r \\phantom{\\, {\\theta}} \\phantom{\\, {\\theta}} } & = & k r^{3} - r \\\\ \\Gamma_{ \\phantom{\\, r} \\, {\\phi} \\, {\\phi} }^{ \\, r \\phantom{\\, {\\phi}} \\phantom{\\, {\\phi}} } & = & {\\left(k r^{3} - r\\right)} \\sin\\left({\\theta}\\right)^{2} \\\\ \\Gamma_{ \\phantom{\\, {\\theta}} \\, t \\, {\\theta} }^{ \\, {\\theta} \\phantom{\\, t} \\phantom{\\, {\\theta}} } & = & \\frac{\\frac{\\partial\\,a}{\\partial t}}{a\\left(t\\right)} \\\\ \\Gamma_{ \\phantom{\\, {\\theta}} \\, r \\, {\\theta} }^{ \\, {\\theta} \\phantom{\\, r} \\phantom{\\, {\\theta}} } & = & \\frac{1}{r} \\\\ \\Gamma_{ \\phantom{\\, {\\theta}} \\, {\\phi} \\, {\\phi} }^{ \\, {\\theta} \\phantom{\\, {\\phi}} \\phantom{\\, {\\phi}} } & = & -\\cos\\left({\\theta}\\right) \\sin\\left({\\theta}\\right) \\\\ \\Gamma_{ \\phantom{\\, {\\phi}} \\, t \\, {\\phi} }^{ \\, {\\phi} \\phantom{\\, t} \\phantom{\\, {\\phi}} } & = & \\frac{\\frac{\\partial\\,a}{\\partial t}}{a\\left(t\\right)} \\\\ \\Gamma_{ \\phantom{\\, {\\phi}} \\, r \\, {\\phi} }^{ \\, {\\phi} \\phantom{\\, r} \\phantom{\\, {\\phi}} } & = & \\frac{1}{r} \\\\ \\Gamma_{ \\phantom{\\, {\\phi}} \\, {\\theta} \\, {\\phi} }^{ \\, {\\phi} \\phantom{\\, {\\theta}} \\phantom{\\, {\\phi}} } & = & \\frac{\\cos\\left({\\theta}\\right)}{\\sin\\left({\\theta}\\right)} \\end{array}$ " ], "text/plain": [ "Gam^t_r,r = -a(t)*d(a)/dt/(k*r^2 - 1) \n", "Gam^t_th,th = r^2*a(t)*d(a)/dt \n", "Gam^t_ph,ph = r^2*a(t)*sin(th)^2*d(a)/dt \n", "Gam^r_t,r = d(a)/dt/a(t) \n", "Gam^r_r,r = -k*r/(k*r^2 - 1) \n", "Gam^r_th,th = k*r^3 - r \n", "Gam^r_ph,ph = (k*r^3 - r)*sin(th)^2 \n", "Gam^th_t,th = d(a)/dt/a(t) \n", "Gam^th_r,th = 1/r \n", "Gam^th_ph,ph = -cos(th)*sin(th) \n", "Gam^ph_t,ph = d(a)/dt/a(t) \n", "Gam^ph_r,ph = 1/r \n", "Gam^ph_th,ph = cos(th)/sin(th) " ] }, "execution_count": 8, "metadata": {}, "output_type": "execute_result" } ], "source": [ "nabla = g.connection()\n", "g.christoffel_symbols_display()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Ricci tensor:" ] }, { "cell_type": "code", "execution_count": 9, "metadata": {}, "outputs": [ { "data": { "text/html": [ " $\\newcommand{\\Bold}[1]{\\mathbf{#1}}\\mathrm{Ric}\\left(g\\right) = -\\frac{3 \\, \\frac{\\partial^2\\,a}{\\partial t ^ 2}}{a\\left(t\\right)} \\mathrm{d} t\\otimes \\mathrm{d} t + \\left( -\\frac{2 \\, \\left(\\frac{\\partial\\,a}{\\partial t}\\right)^{2} + a\\left(t\\right) \\frac{\\partial^2\\,a}{\\partial t ^ 2} + 2 \\, k}{k r^{2} - 1} \\right) \\mathrm{d} r\\otimes \\mathrm{d} r + \\left( 2 \\, r^{2} \\left(\\frac{\\partial\\,a}{\\partial t}\\right)^{2} + r^{2} a\\left(t\\right) \\frac{\\partial^2\\,a}{\\partial t ^ 2} + 2 \\, k r^{2} \\right) \\mathrm{d} {\\theta}\\otimes \\mathrm{d} {\\theta} + {\\left(2 \\, r^{2} \\left(\\frac{\\partial\\,a}{\\partial t}\\right)^{2} + r^{2} a\\left(t\\right) \\frac{\\partial^2\\,a}{\\partial t ^ 2} + 2 \\, k r^{2}\\right)} \\sin\\left({\\theta}\\right)^{2} \\mathrm{d} {\\phi}\\otimes \\mathrm{d} {\\phi}$ " ], "text/plain": [ "Ric(g) = -3*d^2(a)/dt^2/a(t) dt*dt - (2*(d(a)/dt)^2 + a(t)*d^2(a)/dt^2 + 2*k)/(k*r^2 - 1) dr*dr + (2*r^2*(d(a)/dt)^2 + r^2*a(t)*d^2(a)/dt^2 + 2*k*r^2) dth*dth + (2*r^2*(d(a)/dt)^2 + r^2*a(t)*d^2(a)/dt^2 + 2*k*r^2)*sin(th)^2 dph*dph" ] }, "execution_count": 9, "metadata": {}, "output_type": "execute_result" } ], "source": [ "Ricci = nabla.ricci()\n", "Ricci.display()" ] }, { "cell_type": "code", "execution_count": 10, "metadata": {}, "outputs": [ { "data": { "text/html": [ " $\\newcommand{\\Bold}[1]{\\mathbf{#1}}\\begin{array}{lcl} \\mathrm{Ric}\\left(g\\right)_{ \\, t \\, t }^{ \\phantom{\\, t}\\phantom{\\, t} } & = & -\\frac{3 \\, \\frac{\\partial^2\\,a}{\\partial t ^ 2}}{a\\left(t\\right)} \\\\ \\mathrm{Ric}\\left(g\\right)_{ \\, r \\, r }^{ \\phantom{\\, r}\\phantom{\\, r} } & = & -\\frac{2 \\, \\left(\\frac{\\partial\\,a}{\\partial t}\\right)^{2} + a\\left(t\\right) \\frac{\\partial^2\\,a}{\\partial t ^ 2} + 2 \\, k}{k r^{2} - 1} \\\\ \\mathrm{Ric}\\left(g\\right)_{ \\, {\\theta} \\, {\\theta} }^{ \\phantom{\\, {\\theta}}\\phantom{\\, {\\theta}} } & = & 2 \\, r^{2} \\left(\\frac{\\partial\\,a}{\\partial t}\\right)^{2} + r^{2} a\\left(t\\right) \\frac{\\partial^2\\,a}{\\partial t ^ 2} + 2 \\, k r^{2} \\\\ \\mathrm{Ric}\\left(g\\right)_{ \\, {\\phi} \\, {\\phi} }^{ \\phantom{\\, {\\phi}}\\phantom{\\, {\\phi}} } & = & {\\left(2 \\, r^{2} \\left(\\frac{\\partial\\,a}{\\partial t}\\right)^{2} + r^{2} a\\left(t\\right) \\frac{\\partial^2\\,a}{\\partial t ^ 2} + 2 \\, k r^{2}\\right)} \\sin\\left({\\theta}\\right)^{2} \\end{array}$ " ], "text/plain": [ "Ric(g)_t,t = -3*d^2(a)/dt^2/a(t) \n", "Ric(g)_r,r = -(2*(d(a)/dt)^2 + a(t)*d^2(a)/dt^2 + 2*k)/(k*r^2 - 1) \n", "Ric(g)_th,th = 2*r^2*(d(a)/dt)^2 + r^2*a(t)*d^2(a)/dt^2 + 2*k*r^2 \n", "Ric(g)_ph,ph = (2*r^2*(d(a)/dt)^2 + r^2*a(t)*d^2(a)/dt^2 + 2*k*r^2)*sin(th)^2 " ] }, "execution_count": 10, "metadata": {}, "output_type": "execute_result" } ], "source": [ "Ricci.display_comp()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "

Ricci scalar ($R^\\mu_{\\ \\, \\mu}$):

" ] }, { "cell_type": "code", "execution_count": 11, "metadata": {}, "outputs": [ { "data": { "text/html": [ " $\\newcommand{\\Bold}[1]{\\mathbf{#1}}\\begin{array}{llcl} \\mathrm{r}\\left(g\\right):& M & \\longrightarrow & \\mathbb{R} \\\\ & \\left(t, r, {\\theta}, {\\phi}\\right) & \\longmapsto & \\frac{6 \\, {\\left(\\left(\\frac{\\partial\\,a}{\\partial t}\\right)^{2} + a\\left(t\\right) \\frac{\\partial^2\\,a}{\\partial t ^ 2} + k\\right)}}{a\\left(t\\right)^{2}} \\end{array}$ " ], "text/plain": [ "r(g): M --> R\n", " (t, r, th, ph) |--> 6*((d(a)/dt)^2 + a(t)*d^2(a)/dt^2 + k)/a(t)^2" ] }, "execution_count": 11, "metadata": {}, "output_type": "execute_result" } ], "source": [ "Ricci_scalar = g.ricci_scalar()\n", "Ricci_scalar.display()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "

The fluid 4-velocity:

" ] }, { "cell_type": "code", "execution_count": 12, "metadata": {}, "outputs": [ { "data": { "text/html": [ " $\\newcommand{\\Bold}[1]{\\mathbf{#1}}u = \\frac{\\partial}{\\partial t }$ " ], "text/plain": [ "u = d/dt" ] }, "execution_count": 12, "metadata": {}, "output_type": "execute_result" } ], "source": [ "u = M.vector_field('u')\n", "u[0] = 1\n", "u.display()" ] }, { "cell_type": "code", "execution_count": 13, "metadata": {}, "outputs": [ { "data": { "text/html": [ " $\\newcommand{\\Bold}[1]{\\mathbf{#1}}-1$ " ], "text/plain": [ "-1" ] }, "execution_count": 13, "metadata": {}, "output_type": "execute_result" } ], "source": [ "g(u,u).expr()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "`u.dot(u)` is equivalent to `g(u,u)`:" ] }, { "cell_type": "code", "execution_count": 14, "metadata": {}, "outputs": [ { "data": { "text/html": [ " $\\newcommand{\\Bold}[1]{\\mathbf{#1}}-1$ " ], "text/plain": [ "-1" ] }, "execution_count": 14, "metadata": {}, "output_type": "execute_result" } ], "source": [ "u.dot(u).expr()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "

Perfect fluid energy-momentum tensor $T$:

" ] }, { "cell_type": "code", "execution_count": 15, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "Field of symmetric bilinear forms T on the 4-dimensional Lorentzian manifold M\n" ] }, { "data": { "text/html": [ " $\\newcommand{\\Bold}[1]{\\mathbf{#1}}T = \\rho\\left(t\\right) \\mathrm{d} t\\otimes \\mathrm{d} t + \\left( -\\frac{a\\left(t\\right)^{2} p\\left(t\\right)}{k r^{2} - 1} \\right) \\mathrm{d} r\\otimes \\mathrm{d} r + r^{2} a\\left(t\\right)^{2} p\\left(t\\right) \\mathrm{d} {\\theta}\\otimes \\mathrm{d} {\\theta} + r^{2} a\\left(t\\right)^{2} p\\left(t\\right) \\sin\\left({\\theta}\\right)^{2} \\mathrm{d} {\\phi}\\otimes \\mathrm{d} {\\phi}$ " ], "text/plain": [ "T = rho(t) dt*dt - a(t)^2*p(t)/(k*r^2 - 1) dr*dr + r^2*a(t)^2*p(t) dth*dth + r^2*a(t)^2*p(t)*sin(th)^2 dph*dph" ] }, "execution_count": 15, "metadata": {}, "output_type": "execute_result" } ], "source": [ "u_form = u.down(g) # the 1-form associated to u by metric duality\n", "T = (rho+p)*(u_form*u_form) + p*g\n", "T.set_name('T')\n", "print(T)\n", "T.display()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "

The trace of $T$ (we use index notation to denote the double contraction $g^{ab} T_{ab}$):

" ] }, { "cell_type": "code", "execution_count": 16, "metadata": {}, "outputs": [ { "data": { "text/html": [ " $\\newcommand{\\Bold}[1]{\\mathbf{#1}}\\begin{array}{llcl} & M & \\longrightarrow & \\mathbb{R} \\\\ & \\left(t, r, {\\theta}, {\\phi}\\right) & \\longmapsto & 3 \\, p\\left(t\\right) - \\rho\\left(t\\right) \\end{array}$ " ], "text/plain": [ "M --> R\n", "(t, r, th, ph) |--> 3*p(t) - rho(t)" ] }, "execution_count": 16, "metadata": {}, "output_type": "execute_result" } ], "source": [ "Ttrace = g.inverse()['^ab']*T['_ab']\n", "Ttrace.display()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "

Einstein equation: $R_{\\mu \\nu} - {1 \\over 2} R g_{\\mu \\nu} + \\Lambda g_{\\mu \\nu} = {8 \\pi G} T_{\\mu \\nu}$

" ] }, { "cell_type": "code", "execution_count": 17, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "First Friedmann equation:\n", "\n" ] }, { "data": { "text/html": [ " $\\newcommand{\\Bold}[1]{\\mathbf{#1}}-8 \\, \\pi G \\rho\\left(t\\right) - \\Lambda + \\frac{3 \\, \\frac{\\partial}{\\partial t}a\\left(t\\right)^{2}}{a\\left(t\\right)^{2}} + \\frac{3 \\, k}{a\\left(t\\right)^{2}} = 0$ " ], "text/plain": [ "-8*pi*G*rho(t) - Lambda + 3*diff(a(t), t)^2/a(t)^2 + 3*k/a(t)^2 == 0" ] }, "execution_count": 17, "metadata": {}, "output_type": "execute_result" } ], "source": [ "E1 = Ricci - Ricci_scalar/2*g + Lambda*g - (8*pi*G)*T\n", "print(\"First Friedmann equation:\\n\")\n", "E1[0,0].expr().expand() == 0" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "

Trace-reversed version of the Einstein equation: $R_{\\mu \\nu} - \\Lambda g_{\\mu \\nu} = {8 \\pi G} \\left(T_{\\mu \\nu} - {1 \\over 2}T\\,g_{\\mu \\nu}\\right)$

" ] }, { "cell_type": "code", "execution_count": 18, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "Second Friedmann equation:\n", "\n" ] }, { "data": { "text/html": [ " $\\newcommand{\\Bold}[1]{\\mathbf{#1}}-12 \\, \\pi G p\\left(t\\right) - 4 \\, \\pi G \\rho\\left(t\\right) + \\Lambda - \\frac{3 \\, \\frac{\\partial^{2}}{(\\partial t)^{2}}a\\left(t\\right)}{a\\left(t\\right)} = 0$ " ], "text/plain": [ "-12*pi*G*p(t) - 4*pi*G*rho(t) + Lambda - 3*diff(a(t), t, t)/a(t) == 0" ] }, "execution_count": 18, "metadata": {}, "output_type": "execute_result" } ], "source": [ "E2 = Ricci - Lambda*g - (8*pi*G)*(T - Ttrace/2*g)\n", "print(\"Second Friedmann equation:\\n\")\n", "E2[0,0].expr().expand() == 0" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [] } ], "metadata": { "kernelspec": { "display_name": "SageMath 8.3", "language": "", "name": "sagemath" }, "language": "python", "language_info": { "codemirror_mode": { "name": "ipython", "version": 2 }, "file_extension": ".py", "mimetype": "text/x-python", "name": "python", "nbconvert_exporter": "python", "pygments_lexer": "ipython2", "version": "2.7.15" } }, "nbformat": 4, "nbformat_minor": 1 }