Skip to content

Commit

Permalink
Major NRPy+ update: Update ADM spacetime initial data notebooks to be…
Browse files Browse the repository at this point in the history
… consistent with their respective, updated Python modules.
  • Loading branch information
zachetienne committed Dec 28, 2022
1 parent cf5a43d commit 7d2cbfd
Show file tree
Hide file tree
Showing 5 changed files with 226 additions and 196 deletions.
3 changes: 2 additions & 1 deletion BSSN/UIUCBlackHole.py
Original file line number Diff line number Diff line change
Expand Up @@ -107,13 +107,14 @@ def UIUCBlackHole():
KDD[1][2] = KDD[2][1] = -((2*a**3*M*rBL*sp.cos(th)*sp.sin(th)**3)/ \
(SIG*sp.sqrt(AA*SIG)))*(r - rp/4)*sp.sqrt((rBL - rm)/r)

alpha = sp.sympify(1) # We generally choose alpha = 1/psi**2 (psi = BSSN conformal factor) for these initial data
betaU = ixp.zerorank1() # We generally choose \beta^i = 0 for these initial data
BU = ixp.zerorank1() # We generally choose B^i = 0 for these initial data

# Validated against original SENR: KDD[0][2], KDD[1][2], gammaDD[2][2], gammaDD[0][0], gammaDD[1][1]
#print(sp.mathematica_code(gammaDD[1][1]))

# Finally set alpha. We generally choose alpha = 1/psi**2 (psi = BSSN conformal factor)
# for these initial data
try:
cf_type = par.parval_from_str("EvolvedConformalFactor_cf")
except:
Expand Down
82 changes: 41 additions & 41 deletions Tutorial-ADM_Initial_Data-Brill-Lindquist.ipynb
Original file line number Diff line number Diff line change
Expand Up @@ -146,14 +146,14 @@
"psi += BH2_mass / ( 2 * sp.sqrt((Cartxyz[0]-BH2_posn_x)**2 + (Cartxyz[1]-BH2_posn_y)**2 + (Cartxyz[2]-BH2_posn_z)**2) )\n",
"\n",
"# Step 2.c: Set all needed ADM variables in Cartesian coordinates\n",
"gammaCartDD = ixp.zerorank2()\n",
"KCartDD = ixp.zerorank2() # K_{ij} = 0 for these initial data\n",
"gammaDD = ixp.zerorank2()\n",
"KDD = ixp.zerorank2() # K_{ij} = 0 for these initial data\n",
"for i in range(DIM):\n",
" gammaCartDD[i][i] = psi**4\n",
" gammaDD[i][i] = psi**4\n",
"\n",
"alphaCart = 1/psi**2\n",
"betaCartU = ixp.zerorank1() # We generally choose \\beta^i = 0 for these initial data\n",
"BCartU = ixp.zerorank1() # We generally choose B^i = 0 for these initial data"
"alpha = 1/psi**2\n",
"betaU = ixp.zerorank1() # We generally choose \\beta^i = 0 for these initial data\n",
"BU = ixp.zerorank1() # We generally choose B^i = 0 for these initial data"
]
},
{
Expand Down Expand Up @@ -188,31 +188,31 @@
"output_type": "stream",
"text": [
"Consistency check between Brill-Lindquist tutorial and NRPy+ BSSN.BrillLindquist module. ALL SHOULD BE ZERO.\n",
"alphaCart - bl.alphaCart = 0\n",
"betaCartU[0] - bl.betaCartU[0] = 0\n",
"BCartU[0] - bl.BaCartU[0] = 0\n",
"gammaCartDD[0][0] - bl.gammaCartDD[0][0] = 0\n",
"KCartDD[0][0] - bl.KCartDD[0][0] = 0\n",
"gammaCartDD[0][1] - bl.gammaCartDD[0][1] = 0\n",
"KCartDD[0][1] - bl.KCartDD[0][1] = 0\n",
"gammaCartDD[0][2] - bl.gammaCartDD[0][2] = 0\n",
"KCartDD[0][2] - bl.KCartDD[0][2] = 0\n",
"betaCartU[1] - bl.betaCartU[1] = 0\n",
"BCartU[1] - bl.BaCartU[1] = 0\n",
"gammaCartDD[1][0] - bl.gammaCartDD[1][0] = 0\n",
"KCartDD[1][0] - bl.KCartDD[1][0] = 0\n",
"gammaCartDD[1][1] - bl.gammaCartDD[1][1] = 0\n",
"KCartDD[1][1] - bl.KCartDD[1][1] = 0\n",
"gammaCartDD[1][2] - bl.gammaCartDD[1][2] = 0\n",
"KCartDD[1][2] - bl.KCartDD[1][2] = 0\n",
"betaCartU[2] - bl.betaCartU[2] = 0\n",
"BCartU[2] - bl.BaCartU[2] = 0\n",
"gammaCartDD[2][0] - bl.gammaCartDD[2][0] = 0\n",
"KCartDD[2][0] - bl.KCartDD[2][0] = 0\n",
"gammaCartDD[2][1] - bl.gammaCartDD[2][1] = 0\n",
"KCartDD[2][1] - bl.KCartDD[2][1] = 0\n",
"gammaCartDD[2][2] - bl.gammaCartDD[2][2] = 0\n",
"KCartDD[2][2] - bl.KCartDD[2][2] = 0\n"
"alpha - bl.alpha = 0\n",
"betaU[0] - bl.betaU[0] = 0\n",
"BU[0] - bl.BaCartU[0] = 0\n",
"gammaDD[0][0] - bl.gammaDD[0][0] = 0\n",
"KDD[0][0] - bl.KDD[0][0] = 0\n",
"gammaDD[0][1] - bl.gammaDD[0][1] = 0\n",
"KDD[0][1] - bl.KDD[0][1] = 0\n",
"gammaDD[0][2] - bl.gammaDD[0][2] = 0\n",
"KDD[0][2] - bl.KDD[0][2] = 0\n",
"betaU[1] - bl.betaU[1] = 0\n",
"BU[1] - bl.BaCartU[1] = 0\n",
"gammaDD[1][0] - bl.gammaDD[1][0] = 0\n",
"KDD[1][0] - bl.KDD[1][0] = 0\n",
"gammaDD[1][1] - bl.gammaDD[1][1] = 0\n",
"KDD[1][1] - bl.KDD[1][1] = 0\n",
"gammaDD[1][2] - bl.gammaDD[1][2] = 0\n",
"KDD[1][2] - bl.KDD[1][2] = 0\n",
"betaU[2] - bl.betaU[2] = 0\n",
"BU[2] - bl.BaCartU[2] = 0\n",
"gammaDD[2][0] - bl.gammaDD[2][0] = 0\n",
"KDD[2][0] - bl.KDD[2][0] = 0\n",
"gammaDD[2][1] - bl.gammaDD[2][1] = 0\n",
"KDD[2][1] - bl.KDD[2][1] = 0\n",
"gammaDD[2][2] - bl.gammaDD[2][2] = 0\n",
"KDD[2][2] - bl.KDD[2][2] = 0\n"
]
}
],
Expand All @@ -230,16 +230,16 @@
"bl.BrillLindquist()\n",
"\n",
"print(\"Consistency check between Brill-Lindquist tutorial and NRPy+ BSSN.BrillLindquist module. ALL SHOULD BE ZERO.\")\n",
"print(\"alphaCart - bl.alphaCart = \"+str(sp.simplify(alphaCart - bl.alphaCart)))\n",
"print(\"alpha - bl.alpha = \"+str(sp.simplify(alpha - bl.alpha)))\n",
"for i in range(DIM):\n",
" print(\"betaCartU[\"+str(i)+\"] - bl.betaCartU[\"+str(i)+\"] = \"+\\\n",
" str(sp.simplify(betaCartU[i] - bl.betaCartU[i])))\n",
" print(\"BCartU[\"+str(i)+\"] - bl.BaCartU[\"+str(i)+\"] = \"+str(sp.simplify(BCartU[i] - bl.BCartU[i])))\n",
" print(\"betaU[\"+str(i)+\"] - bl.betaU[\"+str(i)+\"] = \"+\\\n",
" str(sp.simplify(betaU[i] - bl.betaU[i])))\n",
" print(\"BU[\"+str(i)+\"] - bl.BaCartU[\"+str(i)+\"] = \"+str(sp.simplify(BU[i] - bl.BU[i])))\n",
" for j in range(DIM):\n",
" print(\"gammaCartDD[\"+str(i)+\"][\"+str(j)+\"] - bl.gammaCartDD[\"+str(i)+\"][\"+str(j)+\"] = \"+\\\n",
" str(sp.simplify(gammaCartDD[i][j] - bl.gammaCartDD[i][j])))\n",
" print(\"KCartDD[\"+str(i)+\"][\"+str(j)+\"] - bl.KCartDD[\"+str(i)+\"][\"+str(j)+\"] = \"+\\\n",
" str(sp.simplify(KCartDD[i][j] - bl.KCartDD[i][j])))"
" print(\"gammaDD[\"+str(i)+\"][\"+str(j)+\"] - bl.gammaDD[\"+str(i)+\"][\"+str(j)+\"] = \"+\\\n",
" str(sp.simplify(gammaDD[i][j] - bl.gammaDD[i][j])))\n",
" print(\"KDD[\"+str(i)+\"][\"+str(j)+\"] - bl.KDD[\"+str(i)+\"][\"+str(j)+\"] = \"+\\\n",
" str(sp.simplify(KDD[i][j] - bl.KDD[i][j])))"
]
},
{
Expand Down Expand Up @@ -284,7 +284,7 @@
],
"metadata": {
"kernelspec": {
"display_name": "Python 3",
"display_name": "Python 3 (ipykernel)",
"language": "python",
"name": "python3"
},
Expand All @@ -298,7 +298,7 @@
"name": "python",
"nbconvert_exporter": "python",
"pygments_lexer": "ipython3",
"version": "3.8.2"
"version": "3.11.1"
}
},
"nbformat": 4,
Expand Down
104 changes: 52 additions & 52 deletions Tutorial-ADM_Initial_Data-ShiftedKerrSchild.ipynb
Original file line number Diff line number Diff line change
Expand Up @@ -185,15 +185,15 @@
"rho2 = rKS*rKS + a*a*sp.cos(th)**2\n",
"\n",
"# alpha = 1/sqrt{1 + M*rKS/rho^2}\n",
"alphaSph = 1/(sp.sqrt(1 + 2*M*rKS/rho2))\n",
"alpha = 1/(sp.sqrt(1 + 2*M*rKS/rho2))\n",
"\n",
"# Initialize the shift vector, \\beta^i, to zero.\n",
"betaSphU = ixp.zerorank1()\n",
"betaU = ixp.zerorank1()\n",
"# beta^r = alpha^2*2Mr/rho^2\n",
"betaSphU[0] = alphaSph*alphaSph*2*M*rKS/rho2\n",
"betaU[0] = alpha*alpha*2*M*rKS/rho2\n",
"\n",
"# Time derivative of shift vector beta^i, B^i, is zero.\n",
"BSphU = ixp.zerorank1()"
"BU = ixp.zerorank1()"
]
},
{
Expand Down Expand Up @@ -241,19 +241,19 @@
"# gamma_{thetatheta}, gamma_{phiphi}\n",
"\n",
"# Initialize \\gamma_{ij} to zero.\n",
"gammaSphDD = ixp.zerorank2()\n",
"gammaDD = ixp.zerorank2()\n",
"\n",
"# gammaDD{rKS rKS} = 1 +2M*rKS/rho^2\n",
"gammaSphDD[0][0] = 1 + 2*M*rKS/rho2\n",
"gammaDD[0][0] = 1 + 2*M*rKS/rho2\n",
"\n",
"# gammaDD{rKS phi} = -a*gammaDD{r r}*sin^2(theta)\n",
"gammaSphDD[0][2] = gammaSphDD[2][0] = -a*gammaSphDD[0][0]*sp.sin(th)**2\n",
"gammaDD[0][2] = gammaDD[2][0] = -a*gammaDD[0][0]*sp.sin(th)**2\n",
"\n",
"# gammaDD{theta theta} = rho^2\n",
"gammaSphDD[1][1] = rho2\n",
"gammaDD[1][1] = rho2\n",
"\n",
"# gammaDD{phi phi} = (rKS^2 + a^2 + 2Mr/rho^2*a^2*sin^2(theta))*sin^2(theta)\n",
"gammaSphDD[2][2] = (rKS*rKS + a*a + 2*M*rKS*a*a*sp.sin(th)**2/rho2)*sp.sin(th)**2"
"gammaDD[2][2] = (rKS*rKS + a*a + 2*M*rKS*a*a*sp.sin(th)**2/rho2)*sp.sin(th)**2"
]
},
{
Expand Down Expand Up @@ -354,18 +354,18 @@
"source": [
"# Step 4: Define the extrinsic curvature in spherical polar coordinates\n",
"# Establish the 3x3 zero-matrix\n",
"KSphDD = ixp.zerorank2()\n",
"KDD = ixp.zerorank2()\n",
"\n",
"# *** Fill in the nonzero components ***\n",
"# *** This will create an upper-triangular matrix ***\n",
"# K_{r r} = D(A+2Mr)/(A^2*B)[4M(a^2*cos(2theta) + a^2 - 2r^2)]\n",
"KSphDD[0][0] = D*(A+2*M*rKS)/(A*A*B)*(4*M*(a*a*sp.cos(2*th)+a*a-2*rKS*rKS))\n",
"KDD[0][0] = D*(A+2*M*rKS)/(A*A*B)*(4*M*(a*a*sp.cos(2*th)+a*a-2*rKS*rKS))\n",
"\n",
"# K_{r theta} = D/(AB)[8a^2*Mr*sin(theta)cos(theta)]\n",
"KSphDD[0][1] = KSphDD[1][0] = D/(A*B)*(8*a*a*M*rKS*sp.sin(th)*sp.cos(th))\n",
"KDD[0][1] = KDD[1][0] = D/(A*B)*(8*a*a*M*rKS*sp.sin(th)*sp.cos(th))\n",
"\n",
"# K_{r phi} = D/A^2[-2aMsin^2(theta)(a^2cos(2theta)+a^2-2r^2)]\n",
"KSphDD[0][2] = KSphDD[2][0] = D/(A*A)*(-2*a*M*sp.sin(th)**2*(a*a*sp.cos(2*th)+a*a-2*rKS*rKS))"
"KDD[0][2] = KDD[2][0] = D/(A*A)*(-2*a*M*sp.sin(th)**2*(a*a*sp.cos(2*th)+a*a-2*rKS*rKS))"
]
},
{
Expand Down Expand Up @@ -399,14 +399,14 @@
"outputs": [],
"source": [
"# K_{theta theta} = D/B[4Mr^2]\n",
"KSphDD[1][1] = D/B*(4*M*rKS*rKS)\n",
"KDD[1][1] = D/B*(4*M*rKS*rKS)\n",
"\n",
"# K_{theta phi} = D/(AB)*(-8*a^3*Mr*sin^3(theta)cos(theta))\n",
"KSphDD[1][2] = KSphDD[2][1] = D/(A*B)*(-8*a**3*M*rKS*sp.sin(th)**3*sp.cos(th))\n",
"KDD[1][2] = KDD[2][1] = D/(A*B)*(-8*a**3*M*rKS*sp.sin(th)**3*sp.cos(th))\n",
"\n",
"# K_{phi phi} = D/(A^2*B)[2Mr*sin^2(theta)(a^4(M+3r)\n",
"# +4a^2r^2(2r-M)+4a^2r*cos(2theta)(a^2+r(M+2r))+8r^5)]\n",
"KSphDD[2][2] = D/(A*A*B)*(2*M*rKS*sp.sin(th)**2*(a**4*(rKS-M)*sp.cos(4*th)\\\n",
"KDD[2][2] = D/(A*A*B)*(2*M*rKS*sp.sin(th)**2*(a**4*(rKS-M)*sp.cos(4*th)\\\n",
" + a**4*(M+3*rKS)+4*a*a*rKS*rKS*(2*rKS-M)\\\n",
" + 4*a*a*rKS*sp.cos(2*th)*(a*a + rKS*(M + 2*rKS)) + 8*rKS**5))"
]
Expand Down Expand Up @@ -442,32 +442,32 @@
"name": "stdout",
"output_type": "stream",
"text": [
"Consistency check between Brill-Lindquist tutorial and NRPy+ BSSN.BrillLindquist module. ALL SHOULD BE ZERO.\n",
"alphaSph - sks.alphaSph = 0\n",
"betaSphU[0] - sks.betaSphU[0] = 0\n",
"BSphU[0] - sks.BaSphU[0] = 0\n",
"gammaSphDD[0][0] - sks.gammaSphDD[0][0] = 0\n",
"KSphDD[0][0] - sks.KSphDD[0][0] = 0\n",
"gammaSphDD[0][1] - sks.gammaSphDD[0][1] = 0\n",
"KSphDD[0][1] - sks.KSphDD[0][1] = 0\n",
"gammaSphDD[0][2] - sks.gammaSphDD[0][2] = 0\n",
"KSphDD[0][2] - sks.KSphDD[0][2] = 0\n",
"betaSphU[1] - sks.betaSphU[1] = 0\n",
"BSphU[1] - sks.BaSphU[1] = 0\n",
"gammaSphDD[1][0] - sks.gammaSphDD[1][0] = 0\n",
"KSphDD[1][0] - sks.KSphDD[1][0] = 0\n",
"gammaSphDD[1][1] - sks.gammaSphDD[1][1] = 0\n",
"KSphDD[1][1] - sks.KSphDD[1][1] = 0\n",
"gammaSphDD[1][2] - sks.gammaSphDD[1][2] = 0\n",
"KSphDD[1][2] - sks.KSphDD[1][2] = 0\n",
"betaSphU[2] - sks.betaSphU[2] = 0\n",
"BSphU[2] - sks.BaSphU[2] = 0\n",
"gammaSphDD[2][0] - sks.gammaSphDD[2][0] = 0\n",
"KSphDD[2][0] - sks.KSphDD[2][0] = 0\n",
"gammaSphDD[2][1] - sks.gammaSphDD[2][1] = 0\n",
"KSphDD[2][1] - sks.KSphDD[2][1] = 0\n",
"gammaSphDD[2][2] - sks.gammaSphDD[2][2] = 0\n",
"KSphDD[2][2] - sks.KSphDD[2][2] = 0\n"
"Consistency check between ShiftedKerrSchild tutorial and NRPy+ BSSN.ShiftedKerrSchild module. ALL SHOULD BE ZERO.\n",
"alpha - sks.alpha = 0\n",
"betaU[0] - sks.betaU[0] = 0\n",
"BU[0] - sks.BSphU[0] = 0\n",
"gammaDD[0][0] - sks.gammaDD[0][0] = 0\n",
"KDD[0][0] - sks.KDD[0][0] = 0\n",
"gammaDD[0][1] - sks.gammaDD[0][1] = 0\n",
"KDD[0][1] - sks.KDD[0][1] = 0\n",
"gammaDD[0][2] - sks.gammaDD[0][2] = 0\n",
"KDD[0][2] - sks.KDD[0][2] = 0\n",
"betaU[1] - sks.betaU[1] = 0\n",
"BU[1] - sks.BSphU[1] = 0\n",
"gammaDD[1][0] - sks.gammaDD[1][0] = 0\n",
"KDD[1][0] - sks.KDD[1][0] = 0\n",
"gammaDD[1][1] - sks.gammaDD[1][1] = 0\n",
"KDD[1][1] - sks.KDD[1][1] = 0\n",
"gammaDD[1][2] - sks.gammaDD[1][2] = 0\n",
"KDD[1][2] - sks.KDD[1][2] = 0\n",
"betaU[2] - sks.betaU[2] = 0\n",
"BU[2] - sks.BSphU[2] = 0\n",
"gammaDD[2][0] - sks.gammaDD[2][0] = 0\n",
"KDD[2][0] - sks.KDD[2][0] = 0\n",
"gammaDD[2][1] - sks.gammaDD[2][1] = 0\n",
"KDD[2][1] - sks.KDD[2][1] = 0\n",
"gammaDD[2][2] - sks.gammaDD[2][2] = 0\n",
"KDD[2][2] - sks.KDD[2][2] = 0\n"
]
}
],
Expand All @@ -483,17 +483,17 @@
"sks.ShiftedKerrSchild()\n",
"\n",
"# It is SAFE to ignore the warning(s) from re-initializing parameters.\n",
"print(\"Consistency check between Brill-Lindquist tutorial and NRPy+ BSSN.BrillLindquist module. ALL SHOULD BE ZERO.\")\n",
"print(\"alphaSph - sks.alphaSph = \"+str(sp.simplify(alphaSph - sks.alphaSph)))\n",
"print(\"Consistency check between ShiftedKerrSchild tutorial and NRPy+ BSSN.ShiftedKerrSchild module. ALL SHOULD BE ZERO.\")\n",
"print(\"alpha - sks.alpha = \"+str(sp.simplify(alpha - sks.alpha)))\n",
"for i in range(DIM):\n",
" print(\"betaSphU[\"+str(i)+\"] - sks.betaSphU[\"+str(i)+\"] = \"+\\\n",
" str(sp.simplify(betaSphU[i] - sks.betaSphU[i])))\n",
" print(\"BSphU[\"+str(i)+\"] - sks.BaSphU[\"+str(i)+\"] = \"+str(sp.simplify(BSphU[i] - sks.BSphU[i])))\n",
" print(\"betaU[\"+str(i)+\"] - sks.betaU[\"+str(i)+\"] = \"+\\\n",
" str(sp.simplify(betaU[i] - sks.betaU[i])))\n",
" print(\"BU[\"+str(i)+\"] - sks.BSphU[\"+str(i)+\"] = \"+str(sp.simplify(BU[i] - sks.BU[i])))\n",
" for j in range(DIM):\n",
" print(\"gammaSphDD[\"+str(i)+\"][\"+str(j)+\"] - sks.gammaSphDD[\"+str(i)+\"][\"+str(j)+\"] = \"+\\\n",
" str(sp.simplify(gammaSphDD[i][j] - sks.gammaSphDD[i][j])))\n",
" print(\"KSphDD[\"+str(i)+\"][\"+str(j)+\"] - sks.KSphDD[\"+str(i)+\"][\"+str(j)+\"] = \"+\\\n",
" str(sp.simplify(KSphDD[i][j] - sks.KSphDD[i][j])))"
" print(\"gammaDD[\"+str(i)+\"][\"+str(j)+\"] - sks.gammaDD[\"+str(i)+\"][\"+str(j)+\"] = \"+\\\n",
" str(sp.simplify(gammaDD[i][j] - sks.gammaDD[i][j])))\n",
" print(\"KDD[\"+str(i)+\"][\"+str(j)+\"] - sks.KDD[\"+str(i)+\"][\"+str(j)+\"] = \"+\\\n",
" str(sp.simplify(KDD[i][j] - sks.KDD[i][j])))"
]
},
{
Expand Down Expand Up @@ -552,7 +552,7 @@
"name": "python",
"nbconvert_exporter": "python",
"pygments_lexer": "ipython3",
"version": "3.10.4"
"version": "3.11.1"
}
},
"nbformat": 4,
Expand Down

0 comments on commit 7d2cbfd

Please sign in to comment.