{ "cells": [ { "cell_type": "markdown", "metadata": {}, "source": [ "## Tolman-Oppenheimer-Volkoff equations" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "*This worksheet illustrates some features of [SageManifolds](http://sagemanifolds.obspm.fr/) (v0.8) on the derivation of the Tolman-Oppenheimer-Volkoff equations (spherically symmetric, stationary solution of general relativity).*" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "We will calculate the Einstein equations \n", "$$R_{\\mu\\nu} - \\frac{1}{2}Rg_{\\mu\\nu} = T_{\\mu\\nu}$$\n", "for a corresponding spherically symmetric, stationary metric $g$. In the above, $R_{\\mu\\nu}$ is the Ricci tensor, $R=R^\\mu_\\mu$ is the Ricci scalar, and $T_{\\mu\\nu}$ is the energy-momentum tensor (left side of Einstein's equations describe the geometry of spacetime, and the right side the matter in the spacetime). " ] }, { "cell_type": "code", "execution_count": 1, "metadata": { "collapsed": true }, "outputs": [], "source": [ "%display text latex\n", "set_nproc()\n", "omit_function_args(True)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "We first declare the spacetime `M` as a general 4-dimensional manifold," ] }, { "cell_type": "code", "execution_count": 2, "metadata": { "collapsed": true }, "outputs": [], "source": [ "M = Manifold(4, 'M')" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "with the standard spherical coordinates (`X` denotes the coordinate chart on `M`):" ] }, { "cell_type": "code", "execution_count": 3, "metadata": { "collapsed": false }, "outputs": [], "source": [ "X. = M.chart(r't r:(0,+oo) th:(0,pi):\\theta phi:(0,2*pi):\\phi')" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "In order to define a general spherically symmetric, stationary metric one needs a few auxiliary functions of the radial coordinate $r$ - metric functions $\\nu(r)$ and $\\lambda(r)$, matter pressure $p(r)$ and energy density $\\rho(r)$, as well as the mass $m(r)$ enclosed within the sphere of the radius $r$: " ] }, { "cell_type": "code", "execution_count": 4, "metadata": { "collapsed": true }, "outputs": [], "source": [ "# metric functions: \n", "nu = function(\"nu\", r, latex_name=r\"\\nu\")\n", "lam = function(\"lambda\", r, latex_name=r\"\\lambda\")\n", "# density and pressure: \n", "rho = function(\"rho\", r, latex_name=r\"\\rho\")\n", "p = function(\"P\", r)\n", "# mass enclosed in sphere of radius r: \n", "m = function(\"m\", r)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "In general, such metric reads as follows: " ] }, { "cell_type": "code", "execution_count": 5, "metadata": { "collapsed": false }, "outputs": [ { "data": { "text/html": [ " $\\newcommand{\\Bold}[1]{\\mathbf{#1}}g = -e^{\\left(2 \\, \\nu\\right)} \\mathrm{d} t\\otimes \\mathrm{d} t + e^{\\left(2 \\, \\lambda\\right)} \\mathrm{d} r\\otimes \\mathrm{d} r + r^{2} \\mathrm{d} {\\theta}\\otimes \\mathrm{d} {\\theta} + r^{2} \\sin\\left({\\theta}\\right)^{2} \\mathrm{d} {\\phi}\\otimes \\mathrm{d} {\\phi}$ " ], "text/plain": [ "g = -e^(2*nu) dt*dt + e^(2*lambda) dr*dr + r^2 dth*dth + r^2*sin(th)^2 dphi*dphi" ] }, "execution_count": 5, "metadata": {}, "output_type": "execute_result" } ], "source": [ "g = M.lorentz_metric('g')\n", "g[0,0] = -exp(2*nu)\n", "g[1,1] = exp(2*lam)\n", "g[2,2], g[3,3] = r^2, r^2*sin(th)^2\n", "g.display()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "which works assuming that the physical constants $G=c=1$. Let's introduce $G$ and $c$ as variables to obtain the dimensional version of the equations:" ] }, { "cell_type": "code", "execution_count": 6, "metadata": { "collapsed": true }, "outputs": [], "source": [ "var('G c pi'); assume(G>0); assume(c>0)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "From the Newtonian weak field limit considerations (Newtonian force far from the source) one may simplify the above expression and replace $\\lambda(r)$ with $\\frac{1-2Gm}{rc^2}$, as well as explicitly put $c^2$ in front of $g_{tt}$. Then" ] }, { "cell_type": "code", "execution_count": 7, "metadata": { "collapsed": false }, "outputs": [ { "data": { "text/html": [ " $\\newcommand{\\Bold}[1]{\\mathbf{#1}}g = -c^{2} e^{\\left(2 \\, \\nu\\right)} \\mathrm{d} t\\otimes \\mathrm{d} t + \\left( \\frac{c^{2} r}{c^{2} r - 2 \\, G m} \\right) \\mathrm{d} r\\otimes \\mathrm{d} r + r^{2} \\mathrm{d} {\\theta}\\otimes \\mathrm{d} {\\theta} + r^{2} \\sin\\left({\\theta}\\right)^{2} \\mathrm{d} {\\phi}\\otimes \\mathrm{d} {\\phi}$ " ], "text/plain": [ "g = -c^2*e^(2*nu) dt*dt + c^2*r/(c^2*r - 2*G*m) dr*dr + r^2 dth*dth + r^2*sin(th)^2 dphi*dphi" ] }, "execution_count": 7, "metadata": {}, "output_type": "execute_result" } ], "source": [ "g[0,0] = -c^2*exp(2*nu)\n", "g[1,1] = 1/(1-2*G*m/(r*c^2))\n", "g.display()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "The Ricci tensor is a result of a method `ricci()` acting on the metric `g`:" ] }, { "cell_type": "code", "execution_count": 8, "metadata": { "collapsed": false }, "outputs": [ { "data": { "text/html": [ " $\\newcommand{\\Bold}[1]{\\mathbf{#1}}\\mathrm{Ric}\\left(g\\right) = \\left( \\frac{c^{2} r^{2} e^{\\left(2 \\, \\nu\\right)} \\left(\\frac{\\partial\\,\\nu}{\\partial r}\\right)^{2} + c^{2} r^{2} e^{\\left(2 \\, \\nu\\right)} \\frac{\\partial^2\\,\\nu}{\\partial r^2} + 2 \\, c^{2} r e^{\\left(2 \\, \\nu\\right)} \\frac{\\partial\\,\\nu}{\\partial r} - {\\left(2 \\, r e^{\\left(2 \\, \\nu\\right)} m \\left(\\frac{\\partial\\,\\nu}{\\partial r}\\right)^{2} + 2 \\, r e^{\\left(2 \\, \\nu\\right)} m \\frac{\\partial^2\\,\\nu}{\\partial r^2} + {\\left(r e^{\\left(2 \\, \\nu\\right)} \\frac{\\partial\\,m}{\\partial r} + 3 \\, e^{\\left(2 \\, \\nu\\right)} m\\right)} \\frac{\\partial\\,\\nu}{\\partial r}\\right)} G}{r^{2}} \\right) \\mathrm{d} t\\otimes \\mathrm{d} t + \\left( -\\frac{c^{2} r^{3} \\left(\\frac{\\partial\\,\\nu}{\\partial r}\\right)^{2} + c^{2} r^{3} \\frac{\\partial^2\\,\\nu}{\\partial r^2} - {\\left(2 \\, r^{2} m \\left(\\frac{\\partial\\,\\nu}{\\partial r}\\right)^{2} + 2 \\, r^{2} m \\frac{\\partial^2\\,\\nu}{\\partial r^2} + 2 \\, r \\frac{\\partial\\,m}{\\partial r} + {\\left(r^{2} \\frac{\\partial\\,m}{\\partial r} - r m\\right)} \\frac{\\partial\\,\\nu}{\\partial r} - 2 \\, m\\right)} G}{c^{2} r^{3} - 2 \\, G r^{2} m} \\right) \\mathrm{d} r\\otimes \\mathrm{d} r + \\left( -\\frac{c^{2} r^{2} \\frac{\\partial\\,\\nu}{\\partial r} - {\\left(2 \\, r m \\frac{\\partial\\,\\nu}{\\partial r} + r \\frac{\\partial\\,m}{\\partial r} + m\\right)} G}{c^{2} r} \\right) \\mathrm{d} {\\theta}\\otimes \\mathrm{d} {\\theta} + \\left( -\\frac{c^{2} r^{2} \\sin\\left({\\theta}\\right)^{2} \\frac{\\partial\\,\\nu}{\\partial r} - {\\left(2 \\, r m \\frac{\\partial\\,\\nu}{\\partial r} + r \\frac{\\partial\\,m}{\\partial r} + m\\right)} G \\sin\\left({\\theta}\\right)^{2}}{c^{2} r} \\right) \\mathrm{d} {\\phi}\\otimes \\mathrm{d} {\\phi}$ " ], "text/plain": [ "Ric(g) = (c^2*r^2*e^(2*nu)*(d(nu)/dr)^2 + c^2*r^2*e^(2*nu)*d^2(nu)/dr^2 + 2*c^2*r*e^(2*nu)*d(nu)/dr - (2*r*e^(2*nu)*m*(d(nu)/dr)^2 + 2*r*e^(2*nu)*m*d^2(nu)/dr^2 + (r*e^(2*nu)*d(m)/dr + 3*e^(2*nu)*m)*d(nu)/dr)*G)/r^2 dt*dt - (c^2*r^3*(d(nu)/dr)^2 + c^2*r^3*d^2(nu)/dr^2 - (2*r^2*m*(d(nu)/dr)^2 + 2*r^2*m*d^2(nu)/dr^2 + 2*r*d(m)/dr + (r^2*d(m)/dr - r*m)*d(nu)/dr - 2*m)*G)/(c^2*r^3 - 2*G*r^2*m) dr*dr - (c^2*r^2*d(nu)/dr - (2*r*m*d(nu)/dr + r*d(m)/dr + m)*G)/(c^2*r) dth*dth - (c^2*r^2*sin(th)^2*d(nu)/dr - (2*r*m*d(nu)/dr + r*d(m)/dr + m)*G*sin(th)^2)/(c^2*r) dphi*dphi" ] }, "execution_count": 8, "metadata": {}, "output_type": "execute_result" } ], "source": [ "Ricci = g.ricci(); Ricci.display()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "For example, the $R_{tt}$ component is " ] }, { "cell_type": "code", "execution_count": 9, "metadata": { "collapsed": false }, "outputs": [ { "data": { "text/html": [ " $\\newcommand{\\Bold}[1]{\\mathbf{#1}}\\frac{c^{2} r^{2} e^{\\left(2 \\, \\nu\\right)} \\left(\\frac{\\partial\\,\\nu}{\\partial r}\\right)^{2} + c^{2} r^{2} e^{\\left(2 \\, \\nu\\right)} \\frac{\\partial^2\\,\\nu}{\\partial r^2} + 2 \\, c^{2} r e^{\\left(2 \\, \\nu\\right)} \\frac{\\partial\\,\\nu}{\\partial r} - {\\left(2 \\, r e^{\\left(2 \\, \\nu\\right)} m \\left(\\frac{\\partial\\,\\nu}{\\partial r}\\right)^{2} + 2 \\, r e^{\\left(2 \\, \\nu\\right)} m \\frac{\\partial^2\\,\\nu}{\\partial r^2} + {\\left(r e^{\\left(2 \\, \\nu\\right)} \\frac{\\partial\\,m}{\\partial r} + 3 \\, e^{\\left(2 \\, \\nu\\right)} m\\right)} \\frac{\\partial\\,\\nu}{\\partial r}\\right)} G}{r^{2}}$ " ], "text/plain": [ "(c^2*r^2*e^(2*nu)*(d(nu)/dr)^2 + c^2*r^2*e^(2*nu)*d^2(nu)/dr^2 + 2*c^2*r*e^(2*nu)*d(nu)/dr - (2*r*e^(2*nu)*m*(d(nu)/dr)^2 + 2*r*e^(2*nu)*m*d^2(nu)/dr^2 + (r*e^(2*nu)*d(m)/dr + 3*e^(2*nu)*m)*d(nu)/dr)*G)/r^2" ] }, "execution_count": 9, "metadata": {}, "output_type": "execute_result" } ], "source": [ "Ricci[0,0]" ] }, { "cell_type": "code", "execution_count": 10, "metadata": { "collapsed": false }, "outputs": [ { "data": { "text/html": [ " $\\newcommand{\\Bold}[1]{\\mathbf{#1}}{\\left(c^{2} \\left(\\frac{\\partial\\,\\nu}{\\partial r}\\right)^{2} - \\frac{2 \\, G m \\left(\\frac{\\partial\\,\\nu}{\\partial r}\\right)^{2}}{r} + c^{2} \\frac{\\partial^2\\,\\nu}{\\partial r^2} + \\frac{2 \\, c^{2} \\frac{\\partial\\,\\nu}{\\partial r}}{r} - \\frac{G \\frac{\\partial\\,m}{\\partial r} \\frac{\\partial\\,\\nu}{\\partial r}}{r} - \\frac{2 \\, G m \\frac{\\partial^2\\,\\nu}{\\partial r^2}}{r} - \\frac{3 \\, G m \\frac{\\partial\\,\\nu}{\\partial r}}{r^{2}}\\right)} e^{\\left(2 \\, \\nu\\right)}$ " ], "text/plain": [ "(c^2*(d(nu)/dr)^2 - 2*G*m*(d(nu)/dr)^2/r + c^2*d^2(nu)/dr^2 + 2*c^2*d(nu)/dr/r - G*d(m)/dr*d(nu)/dr/r - 2*G*m*d^2(nu)/dr^2/r - 3*G*m*d(nu)/dr/r^2)*e^(2*nu)" ] }, "execution_count": 10, "metadata": {}, "output_type": "execute_result" } ], "source": [ "Ricci[0,0].expand().collect(nu.diff(r)).collect(nu.diff(r,r)).collect(c^2*exp(2*nu)).collect_common_factors()" ] }, { "cell_type": "code", "execution_count": 11, "metadata": { "collapsed": false }, "outputs": [ { "data": { "text/html": [ " $\\newcommand{\\Bold}[1]{\\mathbf{#1}}-\\frac{c^{2} r^{3} \\left(\\frac{\\partial\\,\\nu}{\\partial r}\\right)^{2} - 2 \\, G r^{2} m \\left(\\frac{\\partial\\,\\nu}{\\partial r}\\right)^{2} + c^{2} r^{3} \\frac{\\partial^2\\,\\nu}{\\partial r^2} - G r^{2} \\frac{\\partial\\,m}{\\partial r} \\frac{\\partial\\,\\nu}{\\partial r} - 2 \\, G r^{2} m \\frac{\\partial^2\\,\\nu}{\\partial r^2} + G r m \\frac{\\partial\\,\\nu}{\\partial r} - 2 \\, G r \\frac{\\partial\\,m}{\\partial r} + 2 \\, G m}{{\\left(c^{2} r - 2 \\, G m\\right)} r^{2}}$ " ], "text/plain": [ "-(c^2*r^3*(d(nu)/dr)^2 - 2*G*r^2*m*(d(nu)/dr)^2 + c^2*r^3*d^2(nu)/dr^2 - G*r^2*d(m)/dr*d(nu)/dr - 2*G*r^2*m*d^2(nu)/dr^2 + G*r*m*d(nu)/dr - 2*G*r*d(m)/dr + 2*G*m)/((c^2*r - 2*G*m)*r^2)" ] }, "execution_count": 11, "metadata": {}, "output_type": "execute_result" } ], "source": [ "Ricci[1,1].expand().collect_common_factors()" ] }, { "cell_type": "code", "execution_count": 12, "metadata": { "collapsed": false }, "outputs": [ { "data": { "text/html": [ " $\\newcommand{\\Bold}[1]{\\mathbf{#1}}-r \\frac{\\partial\\,\\nu}{\\partial r} + \\frac{2 \\, G m \\frac{\\partial\\,\\nu}{\\partial r}}{c^{2}} + \\frac{G \\frac{\\partial\\,m}{\\partial r}}{c^{2}} + \\frac{G m}{c^{2} r}$ " ], "text/plain": [ "-r*d(nu)/dr + 2*G*m*d(nu)/dr/c^2 + G*d(m)/dr/c^2 + G*m/(c^2*r)" ] }, "execution_count": 12, "metadata": {}, "output_type": "execute_result" } ], "source": [ "Ricci[2,2].expand().collect(nu.diff(r)).collect(nu.diff(r,r)).collect(c^2*exp(2*nu))" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Ricci scalar is obtained by the `ricci_scalar()` method acting on `g`:" ] }, { "cell_type": "code", "execution_count": 13, "metadata": { "collapsed": false }, "outputs": [ { "data": { "text/html": [ " $\\newcommand{\\Bold}[1]{\\mathbf{#1}}-\\frac{2 \\, {\\left(c^{2} r^{2} \\left(\\frac{\\partial\\,\\nu}{\\partial r}\\right)^{2} + c^{2} r^{2} \\frac{\\partial^2\\,\\nu}{\\partial r^2} + 2 \\, c^{2} r \\frac{\\partial\\,\\nu}{\\partial r} - {\\left(2 \\, r m \\left(\\frac{\\partial\\,\\nu}{\\partial r}\\right)^{2} + 2 \\, r m \\frac{\\partial^2\\,\\nu}{\\partial r^2} + {\\left(r \\frac{\\partial\\,m}{\\partial r} + 3 \\, m\\right)} \\frac{\\partial\\,\\nu}{\\partial r} + 2 \\, \\frac{\\partial\\,m}{\\partial r}\\right)} G\\right)}}{c^{2} r^{2}}$ " ], "text/plain": [ "-2*(c^2*r^2*(d(nu)/dr)^2 + c^2*r^2*d^2(nu)/dr^2 + 2*c^2*r*d(nu)/dr - (2*r*m*(d(nu)/dr)^2 + 2*r*m*d^2(nu)/dr^2 + (r*d(m)/dr + 3*m)*d(nu)/dr + 2*d(m)/dr)*G)/(c^2*r^2)" ] }, "execution_count": 13, "metadata": {}, "output_type": "execute_result" } ], "source": [ "Ric_scalar = g.ricci_scalar()\n", "(Ric_scalar.function_chart(X)).collect(nu.diff(r)).collect(nu.diff(r,r)).collect(c^2*exp(2*nu))" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "It is the trace of the Ricci tensor, $R = R_\\mu^\\mu$: " ] }, { "cell_type": "code", "execution_count": 14, "metadata": { "collapsed": false }, "outputs": [ { "data": { "text/html": [ " $\\newcommand{\\Bold}[1]{\\mathbf{#1}}\\mathrm{True}$ " ], "text/plain": [ "True" ] }, "execution_count": 14, "metadata": {}, "output_type": "execute_result" } ], "source": [ "Ric_scalar == Ricci.up(g, 1).trace(0, 1)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Left side of the Einstein equations is" ] }, { "cell_type": "code", "execution_count": 15, "metadata": { "collapsed": false }, "outputs": [ { "data": { "text/html": [ " $\\newcommand{\\Bold}[1]{\\mathbf{#1}}\\mathrm{Ric}\\left(g\\right)-\\mbox{unnamed metric} = \\frac{2 \\, G e^{\\left(2 \\, \\nu\\right)} \\frac{\\partial\\,m}{\\partial r}}{r^{2}} \\mathrm{d} t\\otimes \\mathrm{d} t + \\left( \\frac{2 \\, {\\left(c^{2} r^{2} \\frac{\\partial\\,\\nu}{\\partial r} - {\\left(2 \\, r m \\frac{\\partial\\,\\nu}{\\partial r} + m\\right)} G\\right)}}{c^{2} r^{3} - 2 \\, G r^{2} m} \\right) \\mathrm{d} r\\otimes \\mathrm{d} r + \\left( \\frac{c^{2} r^{3} \\left(\\frac{\\partial\\,\\nu}{\\partial r}\\right)^{2} + c^{2} r^{3} \\frac{\\partial^2\\,\\nu}{\\partial r^2} + c^{2} r^{2} \\frac{\\partial\\,\\nu}{\\partial r} - {\\left(2 \\, r^{2} m \\left(\\frac{\\partial\\,\\nu}{\\partial r}\\right)^{2} + 2 \\, r^{2} m \\frac{\\partial^2\\,\\nu}{\\partial r^2} + r \\frac{\\partial\\,m}{\\partial r} + {\\left(r^{2} \\frac{\\partial\\,m}{\\partial r} + r m\\right)} \\frac{\\partial\\,\\nu}{\\partial r} - m\\right)} G}{c^{2} r} \\right) \\mathrm{d} {\\theta}\\otimes \\mathrm{d} {\\theta} + \\left( -\\frac{{\\left(2 \\, r^{2} m \\left(\\frac{\\partial\\,\\nu}{\\partial r}\\right)^{2} + 2 \\, r^{2} m \\frac{\\partial^2\\,\\nu}{\\partial r^2} + r \\frac{\\partial\\,m}{\\partial r} + {\\left(r^{2} \\frac{\\partial\\,m}{\\partial r} + r m\\right)} \\frac{\\partial\\,\\nu}{\\partial r} - m\\right)} G \\sin\\left({\\theta}\\right)^{2} - {\\left(c^{2} r^{3} \\left(\\frac{\\partial\\,\\nu}{\\partial r}\\right)^{2} + c^{2} r^{3} \\frac{\\partial^2\\,\\nu}{\\partial r^2} + c^{2} r^{2} \\frac{\\partial\\,\\nu}{\\partial r}\\right)} \\sin\\left({\\theta}\\right)^{2}}{c^{2} r} \\right) \\mathrm{d} {\\phi}\\otimes \\mathrm{d} {\\phi}$ " ], "text/plain": [ "Ric(g)-unnamed metric = 2*G*e^(2*nu)*d(m)/dr/r^2 dt*dt + 2*(c^2*r^2*d(nu)/dr - (2*r*m*d(nu)/dr + m)*G)/(c^2*r^3 - 2*G*r^2*m) dr*dr + (c^2*r^3*(d(nu)/dr)^2 + c^2*r^3*d^2(nu)/dr^2 + c^2*r^2*d(nu)/dr - (2*r^2*m*(d(nu)/dr)^2 + 2*r^2*m*d^2(nu)/dr^2 + r*d(m)/dr + (r^2*d(m)/dr + r*m)*d(nu)/dr - m)*G)/(c^2*r) dth*dth - ((2*r^2*m*(d(nu)/dr)^2 + 2*r^2*m*d^2(nu)/dr^2 + r*d(m)/dr + (r^2*d(m)/dr + r*m)*d(nu)/dr - m)*G*sin(th)^2 - (c^2*r^3*(d(nu)/dr)^2 + c^2*r^3*d^2(nu)/dr^2 + c^2*r^2*d(nu)/dr)*sin(th)^2)/(c^2*r) dphi*dphi" ] }, "execution_count": 15, "metadata": {}, "output_type": "execute_result" } ], "source": [ "E = Ricci - (Ric_scalar*g)/2; E.display()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Now for the energy-momentum tensor, $T_{\\mu\\nu}$: " ] }, { "cell_type": "code", "execution_count": 16, "metadata": { "collapsed": false }, "outputs": [ { "data": { "text/html": [ " $\\newcommand{\\Bold}[1]{\\mathbf{#1}}u = e^{\\left(-\\nu\\right)} \\frac{\\partial}{\\partial t }$ " ], "text/plain": [ "u = e^(-nu) d/dt" ] }, "execution_count": 16, "metadata": {}, "output_type": "execute_result" } ], "source": [ "u = M.vector_field('u')\n", "u[0] = exp(-nu)\n", "u.display()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "We can check if it is indeed the timelike 4-vector by checking $u_\\mu u^\\mu = -c^2$ by contracting it with the metric `g` using a method `contract()`: " ] }, { "cell_type": "code", "execution_count": 17, "metadata": { "collapsed": false }, "outputs": [ { "data": { "text/html": [ " $\\newcommand{\\Bold}[1]{\\mathbf{#1}}\\mathrm{True}$ " ], "text/plain": [ "True" ] }, "execution_count": 17, "metadata": {}, "output_type": "execute_result" } ], "source": [ "umuumu = g.contract(0,u,0).contract(0,u,0).function_chart(X)\n", "umuumu == -c^2" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "The product $u_\\mu u^\\mu$ can be also obtained in much a simpler way, by just invoking" ] }, { "cell_type": "code", "execution_count": 18, "metadata": { "collapsed": false }, "outputs": [ { "data": { "text/html": [ " $\\newcommand{\\Bold}[1]{\\mathbf{#1}}\\mathrm{True}$ " ], "text/plain": [ "True" ] }, "execution_count": 18, "metadata": {}, "output_type": "execute_result" } ], "source": [ "umuumu = g(u,u)\n", "umuumu == -c^2" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Let's now addopt $T_{\\mu\\nu}$ in perfect fluid form: " ] }, { "cell_type": "code", "execution_count": 19, "metadata": { "collapsed": false }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "field of symmetric bilinear forms 'T' on the 4-dimensional manifold 'M'\n" ] }, { "data": { "text/html": [ " $\\newcommand{\\Bold}[1]{\\mathbf{#1}}T = c^{4} e^{\\left(2 \\, \\nu\\right)} \\rho \\mathrm{d} t\\otimes \\mathrm{d} t + \\left( \\frac{c^{2} r P}{c^{2} r - 2 \\, G m} \\right) \\mathrm{d} r\\otimes \\mathrm{d} r + r^{2} P \\mathrm{d} {\\theta}\\otimes \\mathrm{d} {\\theta} + r^{2} P \\sin\\left({\\theta}\\right)^{2} \\mathrm{d} {\\phi}\\otimes \\mathrm{d} {\\phi}$ " ], "text/plain": [ "T = c^4*e^(2*nu)*rho dt*dt + c^2*r*P/(c^2*r - 2*G*m) dr*dr + r^2*P dth*dth + r^2*P*sin(th)^2 dphi*dphi" ] }, "execution_count": 19, "metadata": {}, "output_type": "execute_result" } ], "source": [ "u_form = u.down(g)\n", "T = (rho + p/c^2)*(u_form*u_form) + p*g\n", "T.set_name('T')\n", "print T\n", "T.display()" ] }, { "cell_type": "code", "execution_count": 20, "metadata": { "collapsed": false }, "outputs": [ { "data": { "text/html": [ " $\\newcommand{\\Bold}[1]{\\mathbf{#1}}\\begin{array}{llcl} & M & \\longrightarrow & \\mathbb{R} \\\\ & \\left(t, r, {\\theta}, {\\phi}\\right) & \\longmapsto & -c^{2} \\rho + 3 \\, P \\end{array}$ " ], "text/plain": [ "M --> R\n", "(t, r, th, phi) |--> -c^2*rho + 3*P" ] }, "execution_count": 20, "metadata": {}, "output_type": "execute_result" } ], "source": [ "Ttrace = (T.up(g, 0)).trace(0, 1)\n", "Ttrace.display()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Three components of the Einstein equations are as follows - the $E_{tt}$ one: " ] }, { "cell_type": "code", "execution_count": 21, "metadata": { "collapsed": false }, "outputs": [], "source": [ "E0=(E[0,0] - (8*pi*G/c^4)*T[0,0]).expr() == 0" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "A small reorganization of the first equation, using the function solve() to solve for $dm/dr$: " ] }, { "cell_type": "code", "execution_count": 22, "metadata": { "collapsed": false }, "outputs": [], "source": [ "E0 = solve((E0*(-r^2/(G*exp(2*nu))/2)).expand().simplify(), m.diff(r))[0]" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Using SageManifolds ExpressionNice to display the derivatives in textbook form: " ] }, { "cell_type": "code", "execution_count": 23, "metadata": { "collapsed": false }, "outputs": [], "source": [ "from sage.geometry.manifolds.utilities import ExpressionNice" ] }, { "cell_type": "code", "execution_count": 24, "metadata": { "collapsed": false }, "outputs": [ { "data": { "text/html": [ " $\\newcommand{\\Bold}[1]{\\mathbf{#1}}\\frac{\\partial\\,m}{\\partial r} = 4 \\, \\pi r^{2} \\rho$ " ], "text/plain": [ "d(m)/dr == 4*pi*r^2*rho" ] }, "execution_count": 24, "metadata": {}, "output_type": "execute_result" } ], "source": [ "ExpressionNice(E0)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Radial component of Einstein's equations, $E_{rr}$: " ] }, { "cell_type": "code", "execution_count": 25, "metadata": { "collapsed": false }, "outputs": [], "source": [ "E1 = (E[1,1] - (8*pi*G/c^4)*T[1,1]).expr() == 0" ] }, { "cell_type": "code", "execution_count": 26, "metadata": { "collapsed": false }, "outputs": [ { "data": { "text/html": [ " $\\newcommand{\\Bold}[1]{\\mathbf{#1}}\\frac{\\partial\\,\\nu}{\\partial r} = \\frac{{\\left(4 \\, \\pi r^{3} P + c^{2} m\\right)} G}{c^{4} r^{2} - 2 \\, G c^{2} r m}$ " ], "text/plain": [ "d(nu)/dr == (4*pi*r^3*P + c^2*m)*G/(c^4*r^2 - 2*G*c^2*r*m)" ] }, "execution_count": 26, "metadata": {}, "output_type": "execute_result" } ], "source": [ "E1 = solve((E1*(c^4*r^3 - 2*G*c^2*r^2*m)/2).expand().simplify_full(), nu.diff(r))[0]\n", "ExpressionNice(E1)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "For the third equation we use radial part of the energy-momentum conservation equation $\\nabla_\\mu T^{\\mu\\nu}$. First, to define the energy-momentum tensor $T^{\\mu\\nu}$ itself:" ] }, { "cell_type": "code", "execution_count": 27, "metadata": { "collapsed": false }, "outputs": [ { "data": { "text/html": [ " $\\newcommand{\\Bold}[1]{\\mathbf{#1}}\\left(\\begin{array}{rrrr}\n", "e^{\\left(-2 \\, \\nu\\left(r\\right)\\right)} \\rho\\left(r\\right) & 0 & 0 & 0 \\\\\n", "0 & \\frac{c^{2} r P\\left(r\\right) - 2 \\, G P\\left(r\\right) m\\left(r\\right)}{c^{2} r} & 0 & 0 \\\\\n", "0 & 0 & \\frac{P\\left(r\\right)}{r^{2}} & 0 \\\\\n", "0 & 0 & 0 & \\frac{P\\left(r\\right)}{r^{2} \\sin\\left({\\theta}\\right)^{2}}\n", "\\end{array}\\right)$ " ], "text/plain": [ "[ e^(-2*nu(r))*rho(r) 0 0 0]\n", "[ 0 (c^2*r*P(r) - 2*G*P(r)*m(r))/(c^2*r) 0 0]\n", "[ 0 0 P(r)/r^2 0]\n", "[ 0 0 0 P(r)/(r^2*sin(th)^2)]" ] }, "execution_count": 27, "metadata": {}, "output_type": "execute_result" } ], "source": [ "Tup = T.up(g,0).up(g,1)\n", "Tup[:]" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Connection ${\\tt nab}$ for the covariant derivative, and the printout of the non-vanishing Christoffel symbols: " ] }, { "cell_type": "code", "execution_count": 28, "metadata": { "collapsed": false }, "outputs": [ { "data": { "text/html": [ " $\\newcommand{\\Bold}[1]{\\mathbf{#1}}\\begin{array}{lcl} \\Gamma_{ \\phantom{\\, t } \\, t \\, r }^{ \\, t \\phantom{\\, t } \\phantom{\\, r } } & = & \\frac{\\partial\\,\\nu}{\\partial r} \\\\ \\Gamma_{ \\phantom{\\, t } \\, r \\, t }^{ \\, t \\phantom{\\, r } \\phantom{\\, t } } & = & \\frac{\\partial\\,\\nu}{\\partial r} \\\\ \\Gamma_{ \\phantom{\\, r } \\, t \\, t }^{ \\, r \\phantom{\\, t } \\phantom{\\, t } } & = & \\frac{c^{2} r e^{\\left(2 \\, \\nu\\right)} \\frac{\\partial\\,\\nu}{\\partial r} - 2 \\, G e^{\\left(2 \\, \\nu\\right)} m \\frac{\\partial\\,\\nu}{\\partial r}}{r} \\\\ \\Gamma_{ \\phantom{\\, r } \\, r \\, r }^{ \\, r \\phantom{\\, r } \\phantom{\\, r } } & = & \\frac{{\\left(r \\frac{\\partial\\,m}{\\partial r} - m\\right)} G}{c^{2} r^{2} - 2 \\, G r m} \\\\ \\Gamma_{ \\phantom{\\, r } \\, {\\theta} \\, {\\theta} }^{ \\, r \\phantom{\\, {\\theta} } \\phantom{\\, {\\theta} } } & = & -\\frac{c^{2} r - 2 \\, G m}{c^{2}} \\\\ \\Gamma_{ \\phantom{\\, r } \\, {\\phi} \\, {\\phi} }^{ \\, r \\phantom{\\, {\\phi} } \\phantom{\\, {\\phi} } } & = & -\\frac{c^{2} r \\sin\\left({\\theta}\\right)^{2} - 2 \\, G m \\sin\\left({\\theta}\\right)^{2}}{c^{2}} \\\\ \\Gamma_{ \\phantom{\\, {\\theta} } \\, r \\, {\\theta} }^{ \\, {\\theta} \\phantom{\\, r } \\phantom{\\, {\\theta} } } & = & \\frac{1}{r} \\\\ \\Gamma_{ \\phantom{\\, {\\theta} } \\, {\\theta} \\, r }^{ \\, {\\theta} \\phantom{\\, {\\theta} } \\phantom{\\, r } } & = & \\frac{1}{r} \\\\ \\Gamma_{ \\phantom{\\, {\\theta} } \\, {\\phi} \\, {\\phi} }^{ \\, {\\theta} \\phantom{\\, {\\phi} } \\phantom{\\, {\\phi} } } & = & -\\cos\\left({\\theta}\\right) \\sin\\left({\\theta}\\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)} \\\\ \\Gamma_{ \\phantom{\\, {\\phi} } \\, {\\phi} \\, r }^{ \\, {\\phi} \\phantom{\\, {\\phi} } \\phantom{\\, r } } & = & \\frac{1}{r} \\\\ \\Gamma_{ \\phantom{\\, {\\phi} } \\, {\\phi} \\, {\\theta} }^{ \\, {\\phi} \\phantom{\\, {\\phi} } \\phantom{\\, {\\theta} } } & = & \\frac{\\cos\\left({\\theta}\\right)}{\\sin\\left({\\theta}\\right)} \\end{array}$ " ], "text/plain": [ "Gam^t_t,r = d(nu)/dr \n", "Gam^t_r,t = d(nu)/dr \n", "Gam^r_t,t = (c^2*r*e^(2*nu)*d(nu)/dr - 2*G*e^(2*nu)*m*d(nu)/dr)/r \n", "Gam^r_r,r = (r*d(m)/dr - m)*G/(c^2*r^2 - 2*G*r*m) \n", "Gam^r_th,th = -(c^2*r - 2*G*m)/c^2 \n", "Gam^r_phi,phi = -(c^2*r*sin(th)^2 - 2*G*m*sin(th)^2)/c^2 \n", "Gam^th_r,th = 1/r \n", "Gam^th_th,r = 1/r \n", "Gam^th_phi,phi = -cos(th)*sin(th) \n", "Gam^phi_r,phi = 1/r \n", "Gam^phi_th,phi = cos(th)/sin(th) \n", "Gam^phi_phi,r = 1/r \n", "Gam^phi_phi,th = cos(th)/sin(th) " ] }, "execution_count": 28, "metadata": {}, "output_type": "execute_result" } ], "source": [ "nab = g.connection()\n", "nab.display()" ] }, { "cell_type": "code", "execution_count": 29, "metadata": { "collapsed": false }, "outputs": [], "source": [ "co = nab(Tup)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "The following calculates the radial component of $\\nabla_\\mu T^{\\mu\\nu}$: " ] }, { "cell_type": "code", "execution_count": 30, "metadata": { "collapsed": false }, "outputs": [ { "data": { "text/html": [ " $\\newcommand{\\Bold}[1]{\\mathbf{#1}}\\frac{c^{2} r \\frac{\\partial\\,P}{\\partial r} - 2 \\, {\\left(m \\frac{\\partial\\,P}{\\partial r} + {\\left(c^{2} m \\rho + P m\\right)} \\frac{\\partial\\,\\nu}{\\partial r}\\right)} G + {\\left(c^{4} r \\rho + c^{2} r P\\right)} \\frac{\\partial\\,\\nu}{\\partial r}}{c^{2} r}$ " ], "text/plain": [ "(c^2*r*d(P)/dr - 2*(m*d(P)/dr + (c^2*m*rho + P*m)*d(nu)/dr)*G + (c^4*r*rho + c^2*r*P)*d(nu)/dr)/(c^2*r)" ] }, "execution_count": 30, "metadata": {}, "output_type": "execute_result" } ], "source": [ "cosum = 0\n", "# radial component of the covariant derivative: \n", "for i in M.irange():\n", " cosum += co[i,1,i]\n", "cosum" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Solve for $dp/dr$: " ] }, { "cell_type": "code", "execution_count": 31, "metadata": { "collapsed": false }, "outputs": [ { "data": { "text/html": [ " $\\newcommand{\\Bold}[1]{\\mathbf{#1}}\\frac{\\partial\\,P}{\\partial r} = -{\\left(c^{2} \\rho + P\\right)} \\frac{\\partial\\,\\nu}{\\partial r}$ " ], "text/plain": [ "d(P)/dr == -(c^2*rho + P)*d(nu)/dr" ] }, "execution_count": 31, "metadata": {}, "output_type": "execute_result" } ], "source": [ "E2 = solve(cosum.expr(), p.diff(r))[0]\n", "ExpressionNice(E2)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Finally, the three TOV equations: " ] }, { "cell_type": "code", "execution_count": 32, "metadata": { "collapsed": false }, "outputs": [ { "data": { "text/html": [ " $\\newcommand{\\Bold}[1]{\\mathbf{#1}}\\left(\\frac{\\partial\\,m}{\\partial r} = 4 \\, \\pi r^{2} \\rho, \\frac{\\partial\\,\\nu}{\\partial r} = \\frac{{\\left(4 \\, \\pi r^{3} P + c^{2} m\\right)} G}{c^{4} r^{2} - 2 \\, G c^{2} r m}, \\frac{\\partial\\,P}{\\partial r} = -{\\left(c^{2} \\rho + P\\right)} \\frac{\\partial\\,\\nu}{\\partial r}\\right)$ " ], "text/plain": [ "(d(m)/dr == 4*pi*r^2*rho,\n", " d(nu)/dr == (4*pi*r^3*P + c^2*m)*G/(c^4*r^2 - 2*G*c^2*r*m),\n", " d(P)/dr == -(c^2*rho + P)*d(nu)/dr)" ] }, "execution_count": 32, "metadata": {}, "output_type": "execute_result" } ], "source": [ "ExpressionNice(E0), ExpressionNice(E1), ExpressionNice(E2)" ] } ], "metadata": { "kernelspec": { "display_name": "Sage 6.7", "language": "", "name": "sage_6_7" }, "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.8" } }, "nbformat": 4, "nbformat_minor": 0 }