{
"cells": [
{
"cell_type": "markdown",
"metadata": {},
"source": [
"\n",
"\n",
"\n",
"# [Polytropic TOV](https://en.wikipedia.org/wiki/Tolman%E2%80%93Oppenheimer%E2%80%93Volkoff_equation) Initial Data\n",
"\n",
"## Authors: Phil Chang, Zach Etienne, & Leo Werneck\n",
"### Formatting improvements courtesy Brandon Clark\n",
"\n",
"## This module sets up initial data for a [TOV](https://en.wikipedia.org/wiki/Tolman%E2%80%93Oppenheimer%E2%80%93Volkoff_equation) star in *spherical, isotropic coordinates*\n",
"\n",
"**Notebook Status:** Validated \n",
"\n",
"**Validation Notes:** This module has been validated to exhibit convergence to zero of the Hamiltonian constraint violation at the expected order to the exact solution (see [start-to-finish TOV module](Tutorial-Start_to_Finish-BSSNCurvilinear-Setting_up_TOV_initial_data.ipynb) for full test). Note that convergence at the surface of the star is lower order due to the sharp drop to zero in $T^{\\mu\\nu}$.\n",
"\n",
"### NRPy+ Source Code for this module: [TOV/TOV_Solver.py](../edit/TOV/TOV_Solver.py)\n",
"\n",
"[comment]: <> (Introduction: TODO)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"\n",
"\n",
"# Table of Contents\n",
"$$\\label{toc}$$\n",
"\n",
"This notebook is organized as follows:\n",
"\n",
"1. [Step 1](#initializenrpy): Initialize core Python/NRPy+ modules\n",
"1. [Step 2](#tov): The TOV Equations\n",
"1. [Step 3](#code_validation): Code Validation against `TOV.TOV_Solver` NRPy+ module\n",
"1. [Step 4](#latex_pdf_output): Output this notebook to $\\LaTeX$-formatted PDF file"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"\n",
"\n",
"# Step 1: Initialize core Python/NRPy+ modules \\[Back to [top](#toc)\\]\n",
"$$\\label{initializenrpy}$$"
]
},
{
"cell_type": "code",
"execution_count": 1,
"metadata": {
"execution": {
"iopub.execute_input": "2021-03-07T17:06:26.375685Z",
"iopub.status.busy": "2021-03-07T17:06:26.374015Z",
"iopub.status.idle": "2021-03-07T17:06:26.609669Z",
"shell.execute_reply": "2021-03-07T17:06:26.610190Z"
}
},
"outputs": [],
"source": [
"# Step 1: Import needed Python/NRPy+ modules\n",
"import numpy as np # NumPy: A numerical methods module for Python\n",
"import scipy.integrate as si # SciPy: Python module for mathematics, science, and engineering applications\n",
"import math, sys # Standard Python modules for math; multiplatform OS-level functions\n",
"import TOV.Polytropic_EOSs as ppeos # NRPy+: Piecewise polytrope equation of state support"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"\n",
"\n",
"# Step 2: The TOV equations \\[Back to [top](#toc)\\]\n",
"$$\\label{tov}$$\n",
"\n",
"The [TOV line element](https://en.wikipedia.org/wiki/Tolman%E2%80%93Oppenheimer%E2%80%93Volkoff_equation) in terms of the *Schwarzschild coordinate* $r$ is written (in the $-+++$ form):\n",
"$$\n",
"ds^2 = - c^2 e^\\nu dt^2 + \\left(1 - \\frac{2Gm}{rc^2}\\right)^{-1} dr^2 + r^2 d\\Omega^2,\n",
"$$\n",
"where $m(r)$ is the mass-energy enclosed at a given $r$, and is equal to the total star's mass outside the stellar radius $r=R$.\n",
"\n",
"In terms of the *isotropic coordinate* $\\bar{r}$ with $G=c=1$ (i.e., the coordinate system and units we'd prefer to use), the ($-+++$ form) line element is written:\n",
"$$\n",
"ds^2 = - e^{\\nu} dt^2 + e^{4\\phi} \\left(d\\bar{r}^2 + \\bar{r}^2 d\\Omega^2\\right),\n",
"$$\n",
"where $\\phi$ here is the *conformal factor*.\n",
"\n",
"Setting components of the above line element equal to one another, we get (in $G=c=1$ units):\n",
"\n",
"\\begin{align}\n",
"r^2 &= e^{4\\phi} \\bar{r}^2 \\implies e^{4\\phi} = \\frac{r^2}{\\bar{r}^2} \\\\\n",
"\\left(1 - \\frac{2m}{r}\\right)^{-1} dr^2 &= e^{4\\phi} d\\bar{r}^2 \\\\\n",
"\\implies \\frac{d\\bar{r}(r)}{dr} &= \\left(1 - \\frac{2m}{r} \\right)^{-1/2} \\frac{\\bar{r}(r)}{r}.\n",
"\\end{align}\n",
"\n",
"The TOV equations provide radial ODEs for the pressure and $\\nu$ (from [the Wikipedia article on the TOV solution](https://en.wikipedia.org/wiki/Tolman%E2%80%93Oppenheimer%E2%80%93Volkoff_equation)):\n",
"\n",
"\\begin{align}\n",
"\\frac{dP}{dr} &= - \\frac{1}{r} \\left( \\frac{\\rho + P}{2} \\right) \\left(\\frac{2 m}{r} + 8 \\pi r^2 P\\right) \\left(1 - \\frac{2 m}{r}\\right)^{-1} \\\\\n",
"\\frac{d \\nu}{d r} &= \\frac{1}{r}\\left(1 - \\frac{2 m}{r}\\right)^{-1} \\left(\\frac{2 m}{r} + 8 \\pi r^2 P\\right) \\\\\n",
"\\end{align}\n",
"\n",
"Assuming a polytropic equation of state, which relates the pressure $P$ to the baryonic rest-mass density $\\rho_B$,\n",
"\n",
"$$\n",
"P(\\rho_B) = K \\rho_B^\\Gamma,\n",
"$$\n",
"the specific internal energy will be given by\n",
"$$\n",
"\\epsilon = \\frac{P}{\\rho_B (\\Gamma - 1)},\n",
"$$\n",
"\n",
"so the total mass-energy density $\\rho$ is given by\n",
"$$\n",
"\\rho = \\rho_B (1 + \\epsilon).\n",
"$$\n",
"\n",
"Given this, the mass-energy $m(r)$ density is the solution to the ODE:\n",
"$$\n",
"\\frac{dm(r)}{dr} = 4\\pi r^2 \\rho(r)\n",
"$$\n",
"\n",
"Thus the full set of ODEs that need to be solved is given by\n",
"\n",
"$$\n",
"\\boxed{\n",
"\\begin{matrix}\n",
"\\frac{dP}{dr} &=& - \\frac{1}{r} \\left( \\frac{\\rho + P}{2} \\right) \\left(\\frac{2 m}{r} + 8 \\pi r^2 P\\right) \\left(1 - \\frac{2 m}{r}\\right)^{-1} \\\\\n",
"\\frac{d \\nu}{d r} &=& \\frac{1}{r}\\left(1 - \\frac{2 m}{r}\\right)^{-1} \\left(\\frac{2 m}{r} + 8 \\pi r^2 P\\right) \\\\\n",
"\\frac{m(r)}{dr} &=& 4\\pi r^2 \\rho(r) \\\\\n",
"\\frac{d\\bar{r}(r)}{dr} &=& \\left(1 - \\frac{2m}{r} \\right)^{-1/2} \\frac{\\bar{r}(r)}{r}\n",
"\\end{matrix}\n",
"}\\ .\n",
"$$\n",
"\n",
"The following code solves these equations, and was largely written by Phil Chang."
]
},
{
"cell_type": "code",
"execution_count": 2,
"metadata": {
"execution": {
"iopub.execute_input": "2021-03-07T17:06:26.774454Z",
"iopub.status.busy": "2021-03-07T17:06:26.639256Z",
"iopub.status.idle": "2021-03-07T17:06:26.791454Z",
"shell.execute_reply": "2021-03-07T17:06:26.791961Z"
}
},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"1256 1256 1256 1256 1256 1256\n",
"Just generated a TOV star with R_Schw = 9.566044579232513e-01 , M = 1.405030336771405e-01 , M/R_Schw = 1.468768334847266e-01 .\n"
]
}
],
"source": [
"# Step 2: The TOV equations\n",
"\n",
"## TOV SOLVER FOR SINGLE AND PIECEWISE POLYTROPES\n",
"## Authors: Phil Chang, Zachariah B. Etienne, Leo Werneck\n",
"\n",
"# Full documentation for this module may be found in the NRPy+ tutorial Jupyter notebook:\n",
"# Tutorial-Start_to_Finish-BSSNCurvilinear-Setting_up_TOV_initial_data.ipynb\n",
"\n",
"\n",
"# Inputs:\n",
"# * Output data file name\n",
"# * rho_baryon_central, the central density of the TOV star.\n",
"# * n, the polytropic equation of state index. n=1 models cold, degenerate neutron star matter.\n",
"# * K_Polytrope, the polytropic constant.\n",
"# * Verbose output toggle (default = True)\n",
"\n",
"# Output: An initial data file (default file name = \"outputTOVpolytrope.txt\") that well\n",
"# samples the (spherically symmetric) solution both inside and outside the star.\n",
"# It is up to the initial data module to perform the 1D interpolation to generate\n",
"# the solution at arbitrary radius. The file has the following columns:\n",
"# Column 1: Schwarzschild radius\n",
"# Column 2: rho(r), *total* mass-energy density (as opposed to baryonic rest-mass density)\n",
"# Column 3: P(r), Pressure\n",
"# Column 4: m(r), mass enclosed\n",
"# Column 5: e^{nu(r)}, g_{tt}(r)\n",
"# Column 6: e^{4 phi(r)}, conformal factor g_{rr}(r)\n",
"# Column 7: rbar(r), Isotropic radius\n",
"\n",
"# rbar refers to the isotropic radius, and\n",
"# R_Schw refers to the Schwarzschild radius\n",
"\n",
"def TOV_Solver(eos,\n",
" outfile = \"outputTOVpolytrope.txt\",\n",
" rho_baryon_central = 0.129285,\n",
" verbose = True,\n",
" return_M_and_RSchw = False,\n",
" accuracy = \"medium\",\n",
" integrator_type = \"default\",\n",
" no_output_File = False,\n",
" pressure_renormalization=1): # reset the pressure to stellar oscillations studies\n",
"\n",
" def TOV_rhs(r_Schw, y) :\n",
" # In \\tilde units\n",
" #\n",
" P = y[0]\n",
" m = y[1]\n",
" # nu = y[2] # nu is not needed as input into TOV_rhs\n",
" rbar = y[3]\n",
"\n",
" # Compute rho_b and eps_cold, to be used below\n",
" # to compute rho_(total)\n",
" rho_baryon, eps_cold = ppeos.Polytrope_EOS__compute_rhob_and_eps_cold_from_P_cold(eos,P)\n",
"\n",
"# with open(\"rhob_P_cold_and_eps_cold.dat\",\"a+\") as file:\n",
"# file.write(str(r_Schw).format(\"%.15e\")+\" \"+str(rho_baryon).format(\"%.15e\")+\" \"+str(P).format(\"%.15e\")+\" \"+str(eps_cold).format(\"%.15e\")+\"\\n\")\n",
"\n",
" # Compute rho, the *total* mass-energy density:\n",
" # .------------------------------.\n",
" # | rho = (1 + eps)*rho_(baryon) |\n",
" # .------------------------------.\n",
" # with eps = eps_cold, for the initial data.\n",
" rho = (1.0 + eps_cold)*rho_baryon\n",
"\n",
"# m = 4*math.pi/3. * rho*r_Schw**3\n",
" if( r_Schw < 1e-4 or m <= 0.):\n",
" # From https://github.com/natj/tov/blob/master/tov.py#L33:\n",
" # dPdr = -cgs.G*(eden + P/cgs.c**2)*(m + 4.0*pi*r**3*P/cgs.c**2)\n",
" # dPdr = dPdr/(r*(r - 2.0*cgs.G*m/cgs.c**2))\n",
" dPdrSchw = -(rho + P)*(4.*math.pi/3.*r_Schw*rho + 4.*math.pi*r_Schw*P)/(1.-8.*math.pi*rho*r_Schw*r_Schw)\n",
" drbardrSchw = 1./(1. - 8.*math.pi*rho*r_Schw*r_Schw)**0.5\n",
" else:\n",
" dPdrSchw = -(rho + P)*(m + 4.*math.pi*r_Schw**3*P)/(r_Schw*r_Schw*(1.-2.*m/r_Schw))\n",
" drbardrSchw = 1./(1. - 2.*m/r_Schw)**0.5*rbar/r_Schw\n",
"\n",
" dmdrSchw = 4.*math.pi*r_Schw*r_Schw*rho\n",
" dnudrSchw = -2./(P + rho)*dPdrSchw\n",
" return [dPdrSchw, dmdrSchw, dnudrSchw, drbardrSchw]\n",
"\n",
" def integrateStar( eos, P, dumpData = False ):\n",
"\n",
" if accuracy == \"medium\":\n",
" min_step_size = 1e-5\n",
" max_step_size = 1e-2\n",
" integrator = 'dop853'\n",
" elif accuracy == \"low\":\n",
" min_step_size = 1e-3\n",
" max_step_size = 1e-1\n",
" integrator = 'dopri5'\n",
" elif accuracy == \"verylow\":\n",
" min_step_size = 1e-1\n",
" max_step_size = 5e-1\n",
" integrator = 'dopri5'\n",
" elif accuracy == \"high\":\n",
" min_step_size = 1e-5\n",
" max_step_size = 1e-5\n",
" integrator = 'dop853'\n",
" elif accuracy == \"veryhigh\":\n",
" min_step_size = 1e-7\n",
" max_step_size = 1e-6\n",
" integrator = 'dop853'\n",
" else:\n",
" print(\"Unknown accuracy option: \"+str(accuracy))\n",
"\n",
" if integrator_type == \"default\":\n",
" pass\n",
" else:\n",
" integrator = integrator_type\n",
"\n",
" integrator = si.ode(TOV_rhs).set_integrator(integrator)#,rtol=1e-4,atol=1e-4)\n",
"# integrator = si.ode(TOV_rhs).set_integrator('dopri5',rtol=1e-4)\n",
" y0 = [P, 0., 0., 0.]\n",
" r_Schw = 0.\n",
" integrator.set_initial_value(y0,r_Schw)\n",
" dr_Schw = min_step_size\n",
" P = y0[0]\n",
"\n",
" PArr = []\n",
" r_SchwArr = []\n",
" mArr = []\n",
" nuArr = []\n",
" rbarArr = []\n",
"\n",
" while integrator.successful() and P > 1e-19*y0[0] :\n",
" P, m, nu, rbar = integrator.integrate(r_Schw + dr_Schw)\n",
" # Update the value of r_Schw to the latest integrated value\n",
" r_Schw += dr_Schw\n",
"\n",
" dPdrSchw, dmdrSchw, dnudrSchw, drbardrSchw = TOV_rhs( r_Schw, [P,m,nu,rbar])\n",
" dr_Schw = 0.1*min(abs(P/dPdrSchw), abs(m/dmdrSchw))\n",
" dr_Schw = min(dr_Schw, max_step_size)\n",
" PArr.append(P)\n",
" r_SchwArr.append(r_Schw)\n",
" mArr.append(m)\n",
" nuArr.append(nu)\n",
" rbarArr.append(rbar)\n",
"\n",
" M = mArr[-1]\n",
" R_Schw = r_SchwArr[-1]\n",
"\n",
" if no_output_File == True:\n",
" return R_Schw, M\n",
"\n",
" # Apply integration constant to ensure rbar is continuous across TOV surface\n",
" for ii in range(len(rbarArr)):\n",
" rbarArr[ii] *= 0.5*(np.sqrt(R_Schw*(R_Schw - 2.0*M)) + R_Schw - M) / rbarArr[-1]\n",
"\n",
" nuArr_np = np.array(nuArr)\n",
" # Rescale solution to nu so that it satisfies BC: exp(nu(R))=exp(nutilde-nu(r=R)) * (1 - 2m(R)/R)\n",
" # Thus, nu(R) = (nutilde - nu(r=R)) + log(1 - 2*m(R)/R)\n",
" nuArr_np = nuArr_np - nuArr_np[-1] + math.log(1.-2.*mArr[-1]/r_SchwArr[-1])\n",
"\n",
" r_SchwArrExtend_np = 10.**(np.arange(0.01,5.0,0.01))*r_SchwArr[-1]\n",
"\n",
" r_SchwArr.extend(r_SchwArrExtend_np)\n",
" mArr.extend(r_SchwArrExtend_np*0. + M)\n",
" PArr.extend(r_SchwArrExtend_np*0.)\n",
" exp2phiArr_np = np.append( np.exp(nuArr_np), 1. - 2.*M/r_SchwArrExtend_np)\n",
" nuArr.extend(np.log(1. - 2.*M/r_SchwArrExtend_np))\n",
" rbarArr.extend( 0.5*(np.sqrt(r_SchwArrExtend_np**2 - 2.*M*r_SchwArrExtend_np) + r_SchwArrExtend_np - M) )\n",
"\n",
" #phiArr_np = np.append( np.exp(nuArr_np), 1. - 2.*M/r_SchwArrExtend_np)\n",
"\n",
" # Appending to a Python array does what one would reasonably expect.\n",
" # Appending to a numpy array allocates space for a new array with size+1,\n",
" # then copies the data over... over and over... super inefficient.\n",
" r_SchwArr_np = np.array(r_SchwArr)\n",
" PArr_np = np.array(PArr)\n",
" rho_baryonArr_np = np.array(PArr) # This is just to initialize the array\n",
" for j in range(len(PArr_np)):\n",
" # Compute rho_b from P\n",
" rho_baryonArr_np[j] = ppeos.Polytrope_EOS__compute_rhob_from_P_cold(eos,PArr_np[j])\n",
"\n",
" mArr_np = np.array(mArr)\n",
" rbarArr_np = np.array(rbarArr)\n",
" confFactor_exp4phi_np = (r_SchwArr_np/rbarArr_np)**2\n",
"\n",
" # Compute the *total* mass-energy density (as opposed to the *baryonic* mass density)\n",
" rhoArr_np = []\n",
" for i in range(len(PArr)):\n",
" rho_baryon, eps_cold = ppeos.Polytrope_EOS__compute_rhob_and_eps_cold_from_P_cold(eos,PArr[i])\n",
" rho = (1.0 + eps_cold ) * rho_baryon\n",
" rhoArr_np.append(rho)\n",
"\n",
" if verbose:\n",
" print(len(r_SchwArr_np),len(rhoArr_np),len(rho_baryonArr_np),len(PArr_np),len(mArr_np),len(exp2phiArr_np))\n",
"\n",
" PArr_np *= pressure_renormalization # set for pressure renormalization studies\n",
"\n",
" # Special thanks to Leonardo Werneck for pointing out this issue with zip()\n",
" if sys.version_info[0] < 3:\n",
" np.savetxt(outfile, zip(r_SchwArr_np,rhoArr_np,rho_baryonArr_np,PArr_np,mArr_np,exp2phiArr_np,confFactor_exp4phi_np,rbarArr_np),\n",
" fmt=\"%.15e\")\n",
" else:\n",
" np.savetxt(outfile, list(zip(r_SchwArr_np,rhoArr_np,rho_baryonArr_np,PArr_np,mArr_np,exp2phiArr_np,confFactor_exp4phi_np,rbarArr_np)),\n",
" fmt=\"%.15e\")\n",
"\n",
" return R_Schw, M\n",
"\n",
" # Set initial condition from rho_baryon_central\n",
" P_initial_condition = ppeos.Polytrope_EOS__compute_P_cold_from_rhob(eos, rho_baryon_central)\n",
"\n",
" # Integrate the initial condition\n",
" R_Schw_TOV, M_TOV = integrateStar(eos, P_initial_condition, True)\n",
" if verbose:\n",
" print(\"Just generated a TOV star with R_Schw = %.15e , M = %.15e , M/R_Schw = %.15e .\" %(R_Schw_TOV,M_TOV,(M_TOV / R_Schw_TOV)))\n",
"\n",
" if return_M_and_RSchw:\n",
" return M_TOV, R_Schw_TOV\n",
"\n",
"############################\n",
"# Single polytrope example #\n",
"############################\n",
"# Set neos = 1 (single polytrope)\n",
"neos = 1\n",
"\n",
"# Set rho_poly_tab (not needed for a single polytrope)\n",
"rho_poly_tab = []\n",
"\n",
"# Set Gamma_poly_tab\n",
"Gamma_poly_tab = [2.0]\n",
"\n",
"# Set K_poly_tab0\n",
"K_poly_tab0 = 1. # ZACH NOTES: CHANGED FROM 100.\n",
"\n",
"# Set the eos quantities\n",
"eos = ppeos.set_up_EOS_parameters__complete_set_of_input_variables(neos,rho_poly_tab,Gamma_poly_tab,K_poly_tab0)\n",
"\n",
"# Set initial condition (Pressure computed from central density)\n",
"rho_baryon_central = 0.129285\n",
"\n",
"M_TOV, R_Schw_TOV = TOV_Solver(eos,outfile=\"outputTOVpolytrope.txt\",rho_baryon_central=0.129285,verbose = True,\n",
" return_M_and_RSchw=True,accuracy=\"medium\",integrator_type=\"default\",\n",
" pressure_renormalization=1.0)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"\n",
"\n",
"# Step 3: Code Validation \\[Back to [top](#toc)\\]\n",
"$$\\label{code_validation}$$\n",
"\n",
"Here, as a code validation check, we verify agreement in the SymPy expressions for these TOV initial data between\n",
"\n",
"1. this tutorial and \n",
"2. the NRPy+ [TOV.TOV_Solver](../edit/TOV/TOV_Solver.py) module."
]
},
{
"cell_type": "code",
"execution_count": 3,
"metadata": {
"execution": {
"iopub.execute_input": "2021-03-07T17:06:26.937141Z",
"iopub.status.busy": "2021-03-07T17:06:26.797807Z",
"iopub.status.idle": "2021-03-07T17:06:26.954030Z",
"shell.execute_reply": "2021-03-07T17:06:26.954518Z"
}
},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"1256 1256 1256 1256 1256 1256\n",
"Just generated a TOV star with\n",
"* M = 1.405030336771405e-01 ,\n",
"* R_Schw = 9.566044579232513e-01 ,\n",
"* R_iso = 8.100085557410308e-01 ,\n",
"* M/R_Schw = 1.468768334847266e-01 \n",
"\n",
"TOV initial data test PASSED.\n"
]
}
],
"source": [
"# Step 3: Code Validation against TOV.TOV_Solver module\n",
"\n",
"import filecmp\n",
"\n",
"import TOV.TOV_Solver as TOV\n",
"TOV.TOV_Solver(eos,\n",
" outfile=\"outputTOVpolytrope-validation.txt\",\n",
" rho_baryon_central=0.129285,\n",
" verbose = True,\n",
" accuracy=\"medium\",\n",
" integrator_type=\"default\",\n",
" no_output_File = False)\n",
"\n",
"if filecmp.cmp('outputTOVpolytrope.txt',\n",
" 'outputTOVpolytrope-validation.txt') == False:\n",
" print(\"ERROR: TOV initial data test FAILED!\")\n",
" sys.exit(1)\n",
"else:\n",
" print(\"TOV initial data test PASSED.\")"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"\n",
"\n",
"# Step 4: Output this notebook to $\\LaTeX$-formatted PDF file \\[Back to [top](#toc)\\]\n",
"$$\\label{latex_pdf_output}$$\n",
"\n",
"The following code cell converts this Jupyter notebook into a proper, clickable $\\LaTeX$-formatted PDF file. After the cell is successfully run, the generated PDF may be found in the root NRPy+ tutorial directory, with filename\n",
"[Tutorial-ADM_Initial_Data-TOV](Tutorial-ADM_Initial_Data-TOV.pdf) (Note that clicking on this link may not work; you may need to open the PDF file through another means.)"
]
},
{
"cell_type": "code",
"execution_count": 4,
"metadata": {
"execution": {
"iopub.execute_input": "2021-03-07T17:06:26.960083Z",
"iopub.status.busy": "2021-03-07T17:06:26.958715Z",
"iopub.status.idle": "2021-03-07T17:06:30.249332Z",
"shell.execute_reply": "2021-03-07T17:06:30.250235Z"
}
},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"Created Tutorial-ADM_Initial_Data-TOV.tex, and compiled LaTeX file to PDF\n",
" file Tutorial-ADM_Initial_Data-TOV.pdf\n"
]
}
],
"source": [
"import cmdline_helper as cmd # NRPy+: Multi-platform Python command-line interface\n",
"cmd.output_Jupyter_notebook_to_LaTeXed_PDF(\"Tutorial-ADM_Initial_Data-TOV\")"
]
}
],
"metadata": {
"kernelspec": {
"display_name": "Python 3",
"language": "python",
"name": "python3"
},
"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.8.2"
}
},
"nbformat": 4,
"nbformat_minor": 2
}