Merge branch 'the_space_dimension' into debug

MintonGroup · Dec 4, 2022 · 71390fb · 71390fb
2 parents 5bed910 + eadc6a6
commit 71390fb
Show file tree

Hide file tree

Showing 65 changed files with 1,983 additions and 4,294 deletions.
diff --git a/examples/Basic_Simulation/initial_conditions.ipynb b/examples/Basic_Simulation/initial_conditions.ipynb
diff --git a/examples/Basic_Simulation/initial_conditions.py b/examples/Basic_Simulation/initial_conditions.py
@@ -23,7 +23,7 @@
 from numpy.random import default_rng
 
 # Initialize the simulation object as a variable
-sim = swiftest.Simulation(tstart=0.0, tstop=1.0e3, dt=0.01, tstep_out=1.0e0, dump_cadence=0, fragmentation=True, minimum_fragment_mass = 2.5e-11, mtiny=2.5e-8)
+sim = swiftest.Simulation(fragmentation=True, minimum_fragment_mass = 2.5e-11, mtiny=2.5e-8)
 
 # Add the modern planets and the Sun using the JPL Horizons Database
 sim.add_solar_system_body(["Sun","Mercury","Venus","Earth","Mars","Jupiter","Saturn","Uranus","Neptune","Pluto"])
@@ -42,14 +42,10 @@
 GM_pl       = (np.array([6e23, 8e23, 1e24, 3e24, 5e24]) / sim.param['MU2KG']) * sim.GU
 R_pl        = np.full(npl, (3 * (GM_pl / sim.GU) / (4 * np.pi * density_pl)) ** (1.0 / 3.0))
 Rh_pl       = a_pl * ((GM_pl) / (3 * sim.GU)) ** (1.0 / 3.0)
-Ip1_pl      = [0.4, 0.4, 0.4, 0.4, 0.4]
-Ip2_pl      = [0.4, 0.4, 0.4, 0.4, 0.4]
-Ip3_pl      = [0.4, 0.4, 0.4, 0.4, 0.4]
-rotx_pl     = [0.0, 0.0, 0.0, 0.0, 0.0]
-roty_pl     = [0.0, 0.0, 0.0, 0.0, 0.0]
-rotz_pl     = [0.0, 0.0, 0.0, 0.0, 0.0]
+Ip_pl      = np.full((npl,3),0.4,)
+rot_pl     = np.zeros((npl,3))
 
-sim.add_body(name=name_pl, v1=a_pl, v2=e_pl, v3=inc_pl, v4=capom_pl, v5=omega_pl, v6=capm_pl, Gmass=GM_pl, radius=R_pl, rhill=Rh_pl, Ip1=Ip1_pl, Ip2=Ip2_pl, Ip3=Ip3_pl, rotx=rotx_pl, roty=roty_pl, rotz=rotz_pl)
+sim.add_body(name=name_pl, a=a_pl, e=e_pl, inc=inc_pl, capom=capom_pl, omega=omega_pl, capm=capm_pl, Gmass=GM_pl, radius=R_pl, rhill=Rh_pl, Ip=Ip_pl, rot=rot_pl)
 
 # Add 10 user-defined test particles
 ntp = 10
@@ -62,9 +58,9 @@
 omega_tp    = default_rng().uniform(0.0, 360.0, ntp)
 capm_tp     = default_rng().uniform(0.0, 360.0, ntp)
 
-sim.add_body(name=name_tp, v1=a_tp, v2=e_tp, v3=inc_tp, v4=capom_tp, v5=omega_tp, v6=capm_tp)
+sim.add_body(name=name_tp, a=a_tp, e=e_tp, inc=inc_tp, capom=capom_tp, omega=omega_tp, capm=capm_tp)
 # Display the run configuration parameters
 sim.get_parameter()
 
 # Run the simulation
-sim.run()
+sim.run(tstart=0.0, tstop=1.0e3, dt=0.01, tstep_out=1.0e0, dump_cadence=0)
diff --git a/examples/Fragmentation/Fragmentation_Movie.py b/examples/Fragmentation/Fragmentation_Movie.py
@@ -188,13 +188,13 @@ def data_stream(self, frame=0):
        if run_new:
            sim = swiftest.Simulation(param_file=param_file, rotation=True, init_cond_format = "XV", compute_conservation_values=True)
            sim.add_solar_system_body("Sun")
-           sim.add_body(Gmass=body_Gmass[style], radius=body_radius[style], xh=pos_vectors[style], vh=vel_vectors[style], rot=rot_vectors[style])
+           sim.add_body(Gmass=body_Gmass[style], radius=body_radius[style], rh=pos_vectors[style], vh=vel_vectors[style], rot=rot_vectors[style])
 
            # Set fragmentation parameters
            minimum_fragment_gmass = 0.2 * body_Gmass[style][1] # Make the minimum fragment mass a fraction of the smallest body
            gmtiny = 0.99 * body_Gmass[style][1] # Make GMTINY just smaller than the smallest original body. This will prevent runaway collisional cascades
            sim.set_parameter(fragmentation = True, gmtiny=gmtiny, minimum_fragment_gmass=minimum_fragment_gmass, verbose=False)
-           sim.run(dt=1e-8, tstop=2.e-5)
+           sim.run(dt=2e-5, tstop=2.e-5)
        else:
            sim = swiftest.Simulation(param_file=param_file, read_old_output_file=True)
 

diff --git a/examples/helio_gr_test/grsim/param.gr.in b/examples/helio_gr_test/grsim/param.gr.in
diff --git a/...s/helio_gr_test/swiftest_relativity.ipynb → examples/helio_gr_test/helio_gr_test.ipynb b/...s/helio_gr_test/swiftest_relativity.ipynb → examples/helio_gr_test/helio_gr_test.ipynb
@@ -19,7 +19,7 @@
    "metadata": {},
    "outputs": [],
    "source": [
-    "sim_gr = swiftest.Simulation(param_file=\"grsim/param.gr.in\", output_file_name=\"bin.gr.nc\")\n",
+    "sim_gr = swiftest.Simulation(simdir=\"gr\")\n",
     "sim_gr.add_solar_system_body([\"Sun\",\"Mercury\",\"Venus\",\"Earth\",\"Mars\",\"Jupiter\",\"Saturn\",\"Uranus\",\"Neptune\"])"
    ]
   },
@@ -29,7 +29,7 @@
    "metadata": {},
    "outputs": [],
    "source": [
-    "sim_nogr = swiftest.Simulation(param_file=\"nogrsim/param.nogr.in\", output_file_name=\"bin.nogr.nc\")\n",
+    "sim_nogr = swiftest.Simulation(simdir=\"nogr\")\n",
     "sim_nogr.add_solar_system_body([\"Sun\",\"Mercury\",\"Venus\",\"Earth\",\"Mars\",\"Jupiter\",\"Saturn\",\"Uranus\",\"Neptune\"])"
    ]
   },
@@ -39,8 +39,7 @@
    "metadata": {},
    "outputs": [],
    "source": [
-    "tstep_out = 10.0\n",
-    "sim_gr.run(tstop=1000.0, dt=0.005, tstep_out=tstep_out, integrator=\"helio\",general_relativity=True)"
+    "run_args = {\"tstop\":1000.0, \"dt\":0.005, \"tstep_out\":10.0, \"dump_cadence\": 0,\"integrator\":\"helio\"}"
    ]
   },
   {
@@ -49,7 +48,16 @@
    "metadata": {},
    "outputs": [],
    "source": [
-    "sim_nogr.run(tstop=1000.0, dt=0.005, tstep_out=tstep_out, integrator=\"helio\",general_relativity=False)"
+    "sim_gr.run(**run_args,general_relativity=True)"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "sim_nogr.run(**run_args,general_relativity=False)"
    ]
   },
   {
@@ -80,17 +88,6 @@
     "varpi_obs = el['w'] + el['Omega']"
    ]
   },
-  {
-   "cell_type": "code",
-   "execution_count": null,
-   "metadata": {},
-   "outputs": [],
-   "source": [
-    "# Compute the longitude of the periapsis\n",
-    "sim_gr.data['varpi'] = np.mod(sim_gr.data['omega'] + sim_gr.data['capom'],360)\n",
-    "sim_nogr.data['varpi'] = np.mod(sim_nogr.data['omega'] + sim_nogr.data['capom'],360)"
-   ]
-  },
   {
    "cell_type": "code",
    "execution_count": null,
@@ -108,8 +105,8 @@
    "metadata": {},
    "outputs": [],
    "source": [
-    "dvarpi_gr = np.diff(varpisim_gr)  * 3600 * 100 / tstep_out\n",
-    "dvarpi_nogr = np.diff(varpisim_nogr)  * 3600 * 100 / tstep_out\n",
+    "dvarpi_gr = np.diff(varpisim_gr)  * 3600 * 100 / run_args['tstep_out']\n",
+    "dvarpi_nogr = np.diff(varpisim_nogr)  * 3600 * 100 / run_args['tstep_out']\n",
     "dvarpi_obs = np.diff(varpi_obs) / np.diff(t) * 3600 * 100"
    ]
   },
@@ -135,20 +132,13 @@
     "print(f'Obs - Swiftest GR    : {np.mean(dvarpi_obs - dvarpi_gr)}')\n",
     "print(f'Obs - Swiftest No GR : {np.mean(dvarpi_obs - dvarpi_nogr)}')"
    ]
-  },
-  {
-   "cell_type": "code",
-   "execution_count": null,
-   "metadata": {},
-   "outputs": [],
-   "source": []
   }
  ],
  "metadata": {
   "kernelspec": {
-   "display_name": "swiftest",
+   "display_name": "Python (My debug_env Kernel)",
    "language": "python",
-   "name": "swiftest"
+   "name": "debug_env"
   },
   "language_info": {
    "codemirror_mode": {

diff --git a/...ples/helio_gr_test/swiftest_relativity.py → examples/helio_gr_test/helio_gr_test.py b/...ples/helio_gr_test/swiftest_relativity.py → examples/helio_gr_test/helio_gr_test.py
@@ -5,15 +5,16 @@
 import numpy as np
 import matplotlib.pyplot as plt
 
-sim_gr = swiftest.Simulation(param_file="param.gr.in", output_file_name="bin.gr.nc")
+sim_gr = swiftest.Simulation(simdir="gr")
 sim_gr.add_solar_system_body(["Sun","Mercury","Venus","Earth","Mars","Jupiter","Saturn","Uranus","Neptune"])
 
-sim_nogr = swiftest.Simulation(param_file="param.nogr.in", output_file_name="bin.nogr.nc")
+sim_nogr = swiftest.Simulation(simdir="nogr")
 sim_nogr.add_solar_system_body(["Sun","Mercury","Venus","Earth","Mars","Jupiter","Saturn","Uranus","Neptune"])
 
-tstep_out = 10.0
-sim_gr.run(tstop=1000.0, dt=0.005, tstep_out=tstep_out, integrator="helio",general_relativity=True)
-sim_nogr.run(tstop=1000.0, dt=0.005, tstep_out=tstep_out, integrator="helio",general_relativity=False)
+run_args = {"tstop":1000.0, "dt":0.005, "tstep_out":10.0, "dump_cadence": 0,"integrator":"helio"}
+
+sim_gr.run(**run_args,general_relativity=True)
+sim_nogr.run(**run_args,general_relativity=False)
 
 # Get the start and end date of the simulation so we can compare with the real solar system
 start_date = sim_gr.ephemeris_date
@@ -29,19 +30,14 @@
 t = (el['datetime_jd']-el['datetime_jd'][0]) / 365.25
 varpi_obs = el['w'] + el['Omega']
 
-# Compute the longitude of the periapsis
-sim_gr.data['varpi'] = np.mod(sim_gr.data['omega'] + sim_gr.data['capom'],360)
-sim_nogr.data['varpi'] = np.mod(sim_nogr.data['omega'] + sim_nogr.data['capom'],360)
-
 varpisim_gr= sim_gr.data['varpi'].sel(name="Mercury")
 varpisim_nogr= sim_nogr.data['varpi'].sel(name="Mercury")
 tsim = sim_gr.data['time']
 
-dvarpi_gr = np.diff(varpisim_gr)  * 3600 * 100 / tstep_out
-dvarpi_nogr = np.diff(varpisim_nogr)  * 3600 * 100 / tstep_out
+dvarpi_gr = np.diff(varpisim_gr)  * 3600 * 100 / run_args['tstep_out']
+dvarpi_nogr = np.diff(varpisim_nogr)  * 3600 * 100 / run_args['tstep_out']
 dvarpi_obs = np.diff(varpi_obs) / np.diff(t) * 3600 * 100
 
-
 fig, ax = plt.subplots()
 
 ax.plot(t, varpi_obs, label="JPL Horizons",linewidth=2.5)

diff --git a/examples/helio_gr_test/nogrsim/param.nogr.in b/examples/helio_gr_test/nogrsim/param.nogr.in