{ "cells": [ { "cell_type": "markdown", "metadata": {}, "source": [ "\n", "\n", "\n", "# $H_{\\rm Orb, NS}$, up to and including third post-Newtonian order\n", "\n", "## Author: Zach Etienne\n", "\n", "## This notebook constructs the orbital, nonspinning Hamiltonian up to third post-Newtonian order, as summarized in [Buonanno, Chen, and Damour (2006)](https://arxiv.org/abs/gr-qc/0508067) (see references therein for sources).\n", "\n", "**Notebook Status:** Validated \n", "\n", "**Validation Notes:** All expressions in this notebook were transcribed twice by hand on separate occasions, and expressions were corrected as needed to ensure consistency with published PN expressions. In addition, this tutorial notebook has been confirmed to be self-consistent with its corresponding NRPy+ module, as documented [below](#code_validation). **Additional validation tests may have been performed, but are as yet, undocumented.**\n", "\n", "### This notebook exists as the following Python module:\n", "1. [PN_Hamiltonian_NS.py](../../edit/NRPyPN/PN_Hamiltonian_NS.py)\n", "\n", "### This notebook & corresponding Python module depend on the following NRPy+/NRPyPN Python modules:\n", "1. [indexedexp.py](../../edit/indexedexp.py): [**documentation+tutorial**](../Tutorial-Indexed_Expressions.ipynb)\n", "1. [NRPyPN_shortcuts.py](../../edit/NRPyPN/NRPyPN_shortcuts.py): [**documentation**](NRPyPN_shortcuts.ipynb)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "\n", "\n", "# Table of Contents\n", "$$\\label{toc}$$\n", "\n", "1. Part 1: [$H_{\\rm Orb,\\ NS,\\ Newt}+H_{\\rm Orb,\\ NS,\\ 1PN}+H_{\\rm Orb,\\ NS,\\ 2PN}$](#uptotwopn): Contributions up to and including second post-Newtonian order\n", "1. Part 2: [$H_{\\rm Orb,\\ NS,\\ 3PN}$](#threepn): Third post-Newtonian contribution\n", "1. Part 3: [Validation against second transcription and corresponding Python module](#code_validation)\n", "1. Part 4: [LaTeX PDF output](#latex_pdf_output): $\\LaTeX$ PDF Output" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "\n", "\n", "# Part 1: $H_{\\rm Orb,\\ NS,\\ Newt}+H_{\\rm Orb,\\ NS,\\ 1PN}+H_{\\rm Orb,\\ NS,\\ 2PN}$ \\[Back to [top](#toc)\\]\n", "$$\\label{uptotwopn}$$ \n", "\n", "As detailed in [Buonanno, Chen, and Damour (2006)](https://arxiv.org/abs/gr-qc/0508067) (henceforth BCD2006),\n", "this model assumes two point masses of mass $m_1$ and $m_2$ with corresponding momentum vectors $\\mathbf{P}_1$ and $\\mathbf{P}_2$, and displacement vectors $\\mathbf{X}_1$ and $\\mathbf{X}_2$ with respect to the center of mass.\n", "\n", "Following [BCD2006](https://arxiv.org/abs/gr-qc/0508067), we define the following quantities\n", "\n", "\\begin{align}\n", "\\mu &= m_1 m_2 / (m_1+m_2)\\\\\n", "\\eta &= m_1 m_2 / (m_1+m_2)^2\\\\\n", "\\mathbf{p} &= \\mathbf{P}_1/\\mu = -\\mathbf{P}_2/\\mu\\\\\n", "\\mathbf{q} &= (\\mathbf{X}_1-\\mathbf{X}_2)/M\\\\\n", "q &= |\\mathbf{q}|\\\\\n", "\\mathbf{n} &= \\frac{\\mathbf{q}}{q}\n", "\\end{align}\n", "\n", "Then the Hamiltonian up to and including second PN order is given by (to reduce the possibility of copying error, these equations are taken directly from Eqs 2.2-4 of the LaTeX source code of [BCD2006](https://arxiv.org/abs/gr-qc/0508067), and only mildly formatted to (1) improve presentation in Jupyter notebooks and (2) to ensure some degree of consistency in notation across different terms in other Hamiltonian notebooks):\n", "\n", "\\begin{align}\n", "H_{\\rm Newt}\\left({\\bf q},{\\bf p}\\right) &= \\mu \\left[\\frac{{\\bf p}^2}{2} -\n", "\\frac{1}{q}\\right]\\,, \\\\\n", "H_{\\rm 1PN}\\left({\\bf q},{\\bf p}\\right) &= \\mu\\left[\\frac{1}{8}(3\\eta-1)({\\bf p}^2)^2\n", "- \\frac{1}{2}\\left[(3+\\eta){\\bf p}^2+\\eta({\\bf n}\\cdot{\\bf p})^2\\right]\\frac{1}{q} + \\frac{1}{2q^2}\\right]\\,,\\\\\n", "H_{\\rm 2PN}\\left({\\bf q},{\\bf p}\\right)\n", "&= \\mu\\left[\\frac{1}{16}\\left(1-5\\eta+5\\eta^2\\right)({\\bf p}^2)^3\n", "+ \\frac{1}{8} \\left[\n", "\\left(5-20\\eta-3\\eta^2\\right)({\\bf p}^2)^2-2\\eta^2({\\bf n}\\cdot{\\bf p})^2{\\bf p}^2-3\\eta^2({\\bf n}\\cdot{\\bf p})^4 \\right]\\frac{1}{q}\\right.\n", "\\\\\n", "&\\quad\\quad\\quad \\left.+ \\frac{1}{2} \\left[(5+8\\eta){\\bf p}^2+3\\eta({\\bf n}\\cdot{\\bf p})^2\\right]\\frac{1}{q^2}\n", "- \\frac{1}{4}(1+3\\eta)\\frac{1}{q^3}\\right]\\,,\n", "\\end{align}" ] }, { "cell_type": "code", "execution_count": 1, "metadata": { "execution": { "iopub.execute_input": "2021-03-07T17:35:21.733138Z", "iopub.status.busy": "2021-03-07T17:35:21.732281Z", "iopub.status.idle": "2021-03-07T17:35:22.058067Z", "shell.execute_reply": "2021-03-07T17:35:22.057510Z" } }, "outputs": [], "source": [ "# Step 0: Add NRPy's directory to the path\n", "# https://stackoverflow.com/questions/16780014/import-file-from-parent-directory\n", "import sympy as sp # SymPy: The Python computer algebra package upon which NRPy+ depends\n", "import indexedexpNRPyPN as ixp # NRPy+: Symbolic indexed expression (e.g., tensors, vectors, etc.) support\n", "from NRPyPN_shortcuts import div,dot # NRPyPN: shortcuts for e.g., vector operations\n", "\n", "def f_H_Newt__H_NS_1PN__H_NS_2PN(m1,m2, PU, nU, q):\n", " mu = m1*m2 / (m1+m2)\n", " eta = m1*m2 / (m1+m2)**2\n", " pU = ixp.zerorank1()\n", " for i in range(3):\n", " pU[i] = PU[i]/mu\n", "\n", " global H_Newt, H_NS_1PN, H_NS_2PN\n", " H_Newt = mu*(+div(1,2)*dot(pU,pU) - 1/q)\n", "\n", " H_NS_1PN = mu*(+div(1,8)*(3*eta-1)*dot(pU,pU)**2\n", " -div(1,2)*((3+eta)*dot(pU,pU) + eta*dot(nU,pU)**2)/q\n", " +div(1,2)/q**2)\n", "\n", " H_NS_2PN = mu*(+div(1,16)*(1 - 5*eta + 5*eta**2)*dot(pU,pU)**3\n", " +div(1,8)*(+(5 - 20*eta - 3*eta**2)*dot(pU,pU)**2\n", " -2*eta**2*dot(nU,pU)**2*dot(pU,pU)\n", " -3*eta**2*dot(nU,pU)**4)/q\n", " +div(1,2)*((5+8*eta)*dot(pU,pU) + 3*eta*dot(nU,pU)**2)/q**2\n", " -div(1,4)*(1+3*eta)/q**3)" ] }, { "cell_type": "code", "execution_count": 2, "metadata": { "execution": { "iopub.execute_input": "2021-03-07T17:35:22.070031Z", "iopub.status.busy": "2021-03-07T17:35:22.068878Z", "iopub.status.idle": "2021-03-07T17:35:22.072056Z", "shell.execute_reply": "2021-03-07T17:35:22.071450Z" } }, "outputs": [], "source": [ "# Second version. This one was mostly a search+replace version\n", "# of the original TeX'ed up equations in the paper.\n", "# Used for validation purposes only.\n", "def f_H_Newt__H_NS_1PN__H_NS_2PNv2(m1,m2, PU, nU, q):\n", " mu = m1*m2/(m1+m2)\n", " eta = m1*m2/(m1+m2)**2\n", " pU = ixp.zerorank1()\n", " for i in range(3):\n", " pU[i] = PU[i]/mu\n", " p_dot_p = dot(pU,pU)\n", " n_dot_p = dot(nU,pU)\n", "\n", " # H_{\\rm Newt} = \\frac{p^i p^i}{2} - \\frac{1}{q}\n", " global H_Newtv2, H_NS_1PNv2, H_NS_2PNv2\n", " H_Newtv2 = mu*(div(1,2)*p_dot_p - 1/q)\n", " H_NS_1PNv2 = mu*(div(1,8)*(3*eta-1)*p_dot_p**2 - \\\n", " div(1,2)*((3+eta)*p_dot_p + eta*n_dot_p**2)/q + 1/(2*q**2))\n", " H_NS_2PNv2 = mu*(div(1,16)*(1 - 5*eta + 5*eta**2)*p_dot_p**3 +\n", " div(1,8)*((5 - 20*eta - 3*eta**2)*p_dot_p**2\n", " - 2*eta**2*n_dot_p**2*p_dot_p - 3*eta**2*n_dot_p**4)/q +\n", " div(1,2)*((5 + 8*eta)*p_dot_p + 3*eta*n_dot_p**2)/q**2 -\n", " div(1,4)*(1 + 3*eta)/q**3)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "\n", "\n", "# Part 2: $H_{\\rm Orb,\\ NS,\\ 3PN}$ Third post-Newtonian contribution \\[Back to [top](#toc)\\]\n", "$$\\label{threepn}$$ \n", "\n", "To reduce the possibility of copying errors, equations are taken directly from the LaTeX source code of Eqs 2.2-4 in [BCD2006](https://arxiv.org/abs/gr-qc/0508067), and only mildly formatted to (1) improve presentation in Jupyter notebooks and (2) to ensure some degree of consistency in notation across different terms in other Hamiltonian notebooks:\n", "\n", "\\begin{align}\n", "H_{\\rm 3PN}\\left({\\bf q},{\\bf p}\\right)\n", "&= \\mu\\left\\{\\frac{1}{128}\\left(-5+35\\eta-70\\eta^2+35\\eta^3\\right)({\\bf p}^2)^4\\right.\n", " \\\\\n", "&\\quad\\quad + \\frac{1}{16}\\left[\n", "\\left(-7+42\\eta-53\\eta^2-5\\eta^3\\right)({\\bf p}^2)^3\n", "+ (2-3\\eta)\\eta^2({\\bf n}\\cdot{\\bf p})^2({\\bf p}^2)^2\n", "+ 3(1-\\eta)\\eta^2({\\bf n}\\cdot{\\bf p})^4{\\bf p}^2 - 5\\eta^3({\\bf n}\\cdot{\\bf p})^6\n", "\\right]\\frac{1}{q}\n", " \\\\\n", "&\\quad\\quad +\\left[ \\frac{1}{16}\\left(-27+136\\eta+109\\eta^2\\right)({\\bf p}^2)^2\n", "+ \\frac{1}{16}(17+30\\eta)\\eta({\\bf n}\\cdot{\\bf p})^2{\\bf p}^2 + \\frac{1}{12}(5+43\\eta)\\eta({\\bf n}\\cdot{\\bf p})^4\n", "\\right]\\frac{1}{q^2} \\\\\n", "&\\quad\\quad +\\left\\{ \\left[ -\\frac{25}{8} + \\left(\\frac{1}{64}\\pi^2-\\frac{335}{48}\\right)\\eta \n", "- \\frac{23}{8}\\eta^2 \\right]{\\bf p}^2\n", "+ \\left(-\\frac{85}{16}-\\frac{3}{64}\\pi^2-\\frac{7}{4}\\eta\\right)\\eta({\\bf n}\\cdot{\\bf p})^2 \n", "\\right\\}\\frac{1}{q^3}\n", "\\\\\n", "&\\quad\\quad\\left. + \\left[ \\frac{1}{8} + \\left(\\frac{109}{12}-\\frac{21}{32}\\pi^2\\right)\\eta \n", "\\right]\\frac{1}{q^4}\\right\\}\\,,\n", "\\end{align}" ] }, { "cell_type": "code", "execution_count": 3, "metadata": { "execution": { "iopub.execute_input": "2021-03-07T17:35:22.086537Z", "iopub.status.busy": "2021-03-07T17:35:22.085568Z", "iopub.status.idle": "2021-03-07T17:35:22.088646Z", "shell.execute_reply": "2021-03-07T17:35:22.088002Z" } }, "outputs": [], "source": [ "def f_H_NS_3PN(m1,m2, PU, nU, q):\n", " mu = m1*m2 / (m1+m2)\n", " eta = m1*m2 / (m1+m2)**2\n", " pU = ixp.zerorank1()\n", " for i in range(3):\n", " pU[i] = PU[i]/mu\n", "\n", " global H_NS_3PN\n", " H_NS_3PN = mu*(+div(1,128)*(-5 + 35*eta - 70*eta**2 + 35*eta**3)*dot(pU,pU)**4\n", " +div(1, 16)*(+(-7 + 42*eta - 53*eta**2 - 5*eta**3)*dot(pU,pU)**3\n", " +(2-3*eta)*eta**2*dot(nU,pU)**2*dot(pU,pU)**2\n", " +3*(1-eta)*eta**2*dot(nU,pU)**4*dot(pU,pU) - 5*eta**3*dot(nU,pU)**6)/q\n", " +(+div(1,16)*(-27 + 136*eta + 109*eta**2)*dot(pU,pU)**2\n", " +div(1,16)*(+17 + 30*eta)*eta*dot(nU,pU)**2*dot(pU,pU)\n", " +div(1,12)*(+ 5 + 43*eta)*eta*dot(nU,pU)**4)/q**2\n", " +(+(-div(25, 8) + (div(1,64)*sp.pi**2 - div(335,48))*eta - div(23,8)*eta**2)*dot(pU,pU)\n", " +(-div(85,16) - div(3,64)*sp.pi**2 - div(7,4)*eta)*eta*dot(nU,pU)**2)/q**3\n", " +(+div(1,8)+(div(109,12) - div(21,32)*sp.pi**2)*eta)/q**4)" ] }, { "cell_type": "code", "execution_count": 4, "metadata": { "execution": { "iopub.execute_input": "2021-03-07T17:35:22.100331Z", "iopub.status.busy": "2021-03-07T17:35:22.099621Z", "iopub.status.idle": "2021-03-07T17:35:22.101944Z", "shell.execute_reply": "2021-03-07T17:35:22.102467Z" } }, "outputs": [], "source": [ "# Second version. This one was mostly a search+replace version\n", "# of the original TeX'ed up equations in the paper.\n", "# Used for validation purposes only.\n", "def f_H_NS_3PNv2(m1,m2, pU, nU, q):\n", " mu = m1*m2/(m1+m2)\n", " eta = m1*m2/(m1+m2)**2\n", " PU = ixp.zerorank1()\n", " for i in range(3):\n", " PU[i] = pU[i]/mu\n", " P_dot_P = dot(PU,PU)\n", " n_dot_P = dot(nU,PU)\n", "\n", " global H_NS_3PNv2\n", " # The following is simply by-hand search/replaced from the above LaTeX to minimize error\n", " H_NS_3PNv2 = \\\n", " mu*( div(1,128)*(-5+35*eta-70*eta**2+35*eta**3)*P_dot_P**4 +\n", " div(1,16)* ( (-7+42*eta-53*eta**2-5*eta**3)*P_dot_P**3\n", " +(2-3*eta)*eta**2*n_dot_P**2*P_dot_P**2 +\n", " +3*(1-eta)*eta**2*n_dot_P**4*P_dot_P - 5*eta**3*n_dot_P**6 )/(q) +\n", " ( div(1,16)*(-27+136*eta+109*eta**2)*P_dot_P**2\n", " + div(1,16)*(17+30*eta)*eta*n_dot_P**2*P_dot_P + div(1,12)*(5+43*eta)*eta*n_dot_P**4)/(q**2) +\n", " ( ( -div(25,8) + (div(1,64)*sp.pi**2-div(335,48))*eta\n", " - div(23,8)*eta**2 )*P_dot_P\n", " + (-div(85,16)-div(3,64)*sp.pi**2-div(7,4)*eta)*eta*n_dot_P**2)/(q**3) +\n", " ( div(1,8) + (div(109,12)-div(21,32)*sp.pi**2)*eta)/(q**4) )" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "\n", "\n", "# Part 3: Validation against second transcription and corresponding Python module \\[Back to [top](#toc)\\]\n", "$$\\label{code_validation}$$ \n", "\n", "As a code validation check, we verify agreement between \n", "* the SymPy expressions transcribed from the cited published work on two separate occasions, and\n", "* the SymPy expressions generated in this notebook, and the corresponding Python module." ] }, { "cell_type": "code", "execution_count": 5, "metadata": { "execution": { "iopub.execute_input": "2021-03-07T17:35:22.112417Z", "iopub.status.busy": "2021-03-07T17:35:22.111695Z", "iopub.status.idle": "2021-03-07T17:35:22.287409Z", "shell.execute_reply": "2021-03-07T17:35:22.287855Z" } }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "ALL TESTS PASS\n" ] } ], "source": [ "from NRPyPN_shortcuts import m1,m2,pU,nU,q # NRPyPN: import needed input variables.\n", "\n", "f_H_Newt__H_NS_1PN__H_NS_2PN(m1,m2, pU, nU, q)\n", "f_H_NS_3PN(m1,m2, pU, nU, q)\n", "\n", "def error(varname):\n", " print(\"ERROR: When comparing Python module & notebook, \"+varname+\" was found not to match.\")\n", " sys.exit(1)\n", "\n", "# Validation against second transcription of the expressions.\n", "f_H_Newt__H_NS_1PN__H_NS_2PNv2(m1,m2, pU, nU, q)\n", "f_H_NS_3PNv2(m1,m2, pU, nU, q)\n", "if sp.simplify(H_Newt - H_Newtv2) != 0: error(\"H_Newtv2\")\n", "if sp.simplify(H_NS_1PN - H_NS_1PNv2) != 0: error(\"H_NS_1PNv2\")\n", "if sp.simplify(H_NS_2PN - H_NS_2PNv2) != 0: error(\"H_NS_2PNv2\")\n", "if sp.simplify(H_NS_3PN - H_NS_3PNv2) != 0: error(\"H_NS_3PNv2\")\n", "\n", "# Validation against corresponding Python module\n", "import PN_Hamiltonian_NS as HNS\n", "HNS.f_H_Newt__H_NS_1PN__H_NS_2PN(m1,m2, pU, nU, q)\n", "HNS.f_H_NS_3PN(m1,m2, pU, nU, q)\n", "if sp.simplify(H_Newt - HNS.H_Newt) != 0: error(\"H_Newt\")\n", "if sp.simplify(H_NS_1PN - HNS.H_NS_1PN) != 0: error(\"H_NS_1PN\")\n", "if sp.simplify(H_NS_2PN - HNS.H_NS_2PN) != 0: error(\"H_NS_2PN\")\n", "if sp.simplify(H_NS_3PN - HNS.H_NS_3PN) != 0: error(\"H_NS_3PN\")\n", "\n", "print(\"ALL TESTS PASS\")" ] }, { "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", "[PN-Hamiltonian-Nonspinning.pdf](PN-Hamiltonian-Nonspinning.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": 6, "metadata": { "execution": { "iopub.execute_input": "2021-03-07T17:35:22.292627Z", "iopub.status.busy": "2021-03-07T17:35:22.291934Z", "iopub.status.idle": "2021-03-07T17:35:25.582036Z", "shell.execute_reply": "2021-03-07T17:35:25.582631Z" } }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "Created PN-Hamiltonian-Nonspinning.tex, and compiled LaTeX file to PDF file\n", " PN-Hamiltonian-Nonspinning.pdf\n" ] } ], "source": [ "import os,sys # Standard Python modules for multiplatform OS-level functions\n", "import cmdline_helperNRPyPN as cmd # NRPy+: Multi-platform Python command-line interface\n", "cmd.output_Jupyter_notebook_to_LaTeXed_PDF(\"PN-Hamiltonian-Nonspinning\")" ] } ], "metadata": { "kernelspec": { "display_name": "Python 3 (ipykernel)", "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.10.11" } }, "nbformat": 4, "nbformat_minor": 4 }