Merge branch 'OOPSymba' into OOPTides

examples/symba_swifter_comparison/8pl_16tp_encounters/cb.in

@@ -0,0 +1,5 @@
+0
+0.00029591220819207774
+0.004650467260962157
+4.7535806948127355e-12
+-2.2473967953572827e-18

examples/symba_swifter_comparison/8pl_16tp_encounters/cb.swiftest.in

@@ -0,0 +1,5 @@
+0
+0.00029591220819207774
+0.004650467260962157
+4.7535806948127355e-12
+-2.2473967953572827e-18

examples/symba_swifter_comparison/8pl_16tp_encounters/init_cond.py

@@ -0,0 +1,134 @@
+#!/usr/bin/env python3
+import numpy as np
+import swiftest
+import swiftest.io as swio
+import astropy.constants as const
+import sys
+import xarray as xr
+
+# Both codes use the same tp input file
+tpin = "tp.in"
+
+swifter_input  = "param.swifter.in"
+swifter_pl     = "pl.swifter.in"
+swifter_bin    = "bin.swifter.dat"
+swifter_enc    = "enc.swifter.dat"
+
+swiftest_input = "param.swiftest.in"
+swiftest_pl    = "pl.swiftest.in"
+swiftest_cb    = "cb.swiftest.in"
+swiftest_bin   = "bin.swiftest.dat"
+swiftest_enc   = "enc.swiftest.dat"
+
+sim = swiftest.Simulation()
+
+sim.param['T0'] = 0.0
+sim.param['DT'] = 1.0 
+sim.param['TSTOP'] = 365.25e1
+sim.param['ISTEP_OUT']  = 11
+sim.param['ISTEP_DUMP'] = 1
+sim.param['CHK_QMIN_COORD'] = "HELIO"
+sim.param['CHK_QMIN'] = swiftest.RSun / swiftest.AU2M
+sim.param['CHK_QMIN_RANGE'] = f"{swiftest.RSun / swiftest.AU2M} 1000.0"
+sim.param['CHK_RMIN'] = swiftest.RSun / swiftest.AU2M
+sim.param['CHK_RMAX'] = 1000.0
+sim.param['CHK_EJECT'] = 1000.0
+sim.param['OUT_FORM'] = "XV"
+sim.param['OUT_STAT'] = "UNKNOWN"
+sim.param['GR'] = 'NO'
+sim.param['CHK_CLOSE'] = 'YES'
+sim.param['RHILL_PRESENT'] = 'YES'
+sim.param['MTINY'] = 1.0e-12
+
+sim.param['MU2KG'] = swiftest.MSun
+sim.param['TU2S'] = swiftest.JD2S
+sim.param['DU2M'] = swiftest.AU2M
+
+bodyid = {
+   "Sun": 0,
+   "Mercury": 1,
+   "Venus": 2,
+   "Earth": 3,
+   "Mars": 4,
+   "Jupiter": 5,
+   "Saturn": 6,
+   "Uranus": 7,
+   "Neptune": 8,
+}
+
+for name, id in bodyid.items():
+   sim.add(name, idval=id)
+
+ds = sim.ds
+cb = ds.sel(id=0)
+pl = ds.where(ds.id > 0, drop=True)
+npl = pl.id.size
+
+ntp = 16
+dims = ['time', 'id', 'vec']
+tp = []
+t = np.array([0.0])
+clab, plab, tlab = swio.make_swiftest_labels(sim.param)
+
+# For each planet, we will initialize a pair of test particles. One on its way in, and one on its way out. We will also initialize two additional particles that don't encounter anything
+tpnames = np.arange(101, 101 + ntp)
+tpxv1 = np.empty((6))
+tpxv2 = np.empty((6))
+
+p1 = []
+p2 = []
+p3 = []
+p4 = []
+p5 = []
+p6 = []
+
+for i in pl.id:
+    pli = pl.sel(id=i)
+    rstart = 2 * np.double(pli['Radius'])  # Start the test particles at a multiple of the planet radius away
+    vstart = 1.5 * np.sqrt(2 * np.double(pli['Mass'])  / rstart)  # Start the test particle velocities at a multiple of the escape speed
+    xvstart = np.array([rstart / np.sqrt(2.0), rstart / np.sqrt(2.0), 0.0, vstart, 0.0, 0.0])
+    # The positions and velocities of each pair of test particles will be in reference to a planet
+    plvec = np.array([np.double(pli['px']),
+                      np.double(pli['py']),
+                      np.double(pli['pz']),
+                      np.double(pli['vx']),
+                      np.double(pli['vy']),
+                      np.double(pli['vz'])])
+    tpxv1 = plvec + xvstart
+    tpxv2 = plvec - xvstart
+    p1.append(tpxv1[0])
+    p1.append(tpxv2[0])
+    p2.append(tpxv1[1])
+    p2.append(tpxv2[1])
+    p3.append(tpxv1[2])
+    p3.append(tpxv2[2])
+    p4.append(tpxv1[3])
+    p4.append(tpxv2[3])
+    p5.append(tpxv1[4])
+    p5.append(tpxv2[4])
+    p6.append(tpxv1[5])
+    p6.append(tpxv2[5])
+
+tvec = np.vstack([p1, p2, p3, p4, p5, p6])
+tpframe = np.expand_dims(tvec.T, axis=0)
+tpxr = xr.DataArray(tpframe, dims = dims, coords = {'time' : t, 'id' : tpnames, 'vec' : tlab})
+
+tp = [tpxr]
+tpda = xr.concat(tp,dim='time')
+tpds = tpda.to_dataset(dim = 'vec')
+
+sim.ds = xr.combine_by_coords([sim.ds, tpds])
+swio.swiftest_xr2infile(sim.ds, sim.param)
+
+sim.param['PL_IN'] = swiftest_pl
+sim.param['TP_IN'] = tpin
+sim.param['CB_IN'] = swiftest_cb
+sim.param['BIN_OUT'] = swiftest_bin
+sim.param['ENC_OUT'] = swiftest_enc
+sim.save(swiftest_input)
+
+sim.param['PL_IN'] = swifter_pl
+sim.param['TP_IN'] = tpin
+sim.param['BIN_OUT'] = swifter_bin
+sim.param['ENC_OUT'] = swifter_enc
+sim.save(swifter_input, codename="Swifter")

examples/symba_swifter_comparison/8pl_16tp_encounters/param.swifter.in

@@ -0,0 +1,26 @@
+! VERSION        Swifter parameter file converted from Swiftest
+T0               0.0
+TSTOP            3652.5
+DT               1.0
+ISTEP_OUT        11
+ISTEP_DUMP       1
+OUT_FORM         XV
+OUT_TYPE         REAL8
+OUT_STAT         UNKNOWN
+IN_TYPE          ASCII
+PL_IN            pl.swifter.in
+TP_IN            tp.in
+BIN_OUT          bin.swifter.dat
+ENC_OUT          enc.swifter.dat
+CHK_QMIN         0.004650467260962157
+CHK_RMIN         0.004650467260962157
+CHK_RMAX         1000.0
+CHK_EJECT        1000.0
+CHK_QMIN_COORD   HELIO
+CHK_QMIN_RANGE   0.004650467260962157 1000.0
+EXTRA_FORCE      NO
+BIG_DISCARD      NO
+CHK_CLOSE        YES
+RHILL_PRESENT    YES
+J2               4.7535806948127355e-12
+J4               -2.2473967953572827e-18

examples/symba_swifter_comparison/8pl_16tp_encounters/param.swiftest.in

@@ -0,0 +1,36 @@
+! VERSION        Swiftest parameter input
+T0               0.0
+TSTOP            3652.5
+DT               1.0
+ISTEP_OUT        11
+ISTEP_DUMP       1
+OUT_FORM         XV
+OUT_TYPE         REAL8
+OUT_STAT         UNKNOWN
+IN_TYPE          ASCII
+PL_IN            pl.swiftest.in
+TP_IN            tp.in
+CB_IN            cb.swiftest.in
+BIN_OUT          bin.swiftest.dat
+ENC_OUT          enc.swiftest.dat
+CHK_QMIN         0.004650467260962157
+CHK_RMIN         0.004650467260962157
+CHK_RMAX         1000.0
+CHK_EJECT        1000.0
+CHK_QMIN_COORD   HELIO
+CHK_QMIN_RANGE   0.004650467260962157 1000.0
+MU2KG            1.988409870698051e+30
+TU2S             86400
+DU2M             149597870700.0
+EXTRA_FORCE      NO
+BIG_DISCARD      NO
+CHK_CLOSE        YES
+RHILL_PRESENT    YES
+FRAGMENTATION    NO
+ROTATION         NO
+TIDES            NO
+ENERGY           NO
+GR               NO
+YARKOVSKY        NO
+YORP             NO
+MTINY            1e-12

examples/symba_swifter_comparison/8pl_16tp_encounters/pl.in

@@ -0,0 +1,33 @@
+8
+1 4.9125474498983623693e-11 0.0014751238438755500459
+1.6306381826061645943e-05
+-0.065841771551149230746 0.30388831943526661838 0.030872485461978960153
+-0.033141166233939436947 -0.0049297226604189817514 0.0026371811668407158825
+2 7.243452483873646905e-10 0.006759080797928606587
+4.0453784346544178454e-05
+-0.65269716062695148917 -0.3065765656441301057 0.033456491497379246824
+0.008459831335658639026 -0.0184014319837384685 -0.0007407193515014080928
+3 8.9970113821660187435e-10 0.010044868190633438806
+4.25875607065040958e-05
+0.58046286084934750615 -0.8332000042504307258 3.7646553415201541957e-05
+0.013836557832279990782 0.009770187318278569788 -5.1179589633921335467e-07
+4 9.549535102761465607e-11 0.0072467082986392815006
+2.265740805092889601e-05
+-1.5891417403740180081 0.4938480736359250889 0.049330990309104823244
+-0.0036308073545784510204 -0.012168467501132099878 -0.00016594932370266260858
+5 2.825345908631354893e-07 0.3552707649709459117
+0.00046732617030490929307
+4.1148395833578952363 -2.8938323061728068453 -0.080043092204059404504
+0.0042549773877191511204 0.006534697671907701254 -0.00012233719535540690457
+6 8.459715183006415395e-08 0.43765596788571493287
+0.00038925687730393611812
+6.3589256477393849565 -7.653288021415167286 -0.12000977499446359442
+0.003985370599203661747 0.0035590677039893160206 -0.00022043610541731448703
+7 1.2920249163736673626e-08 0.46957663585116591335
+0.00016953449859497231466
+14.816779495279050138 13.049265812461410263 -0.14351615042000470668
+-0.0026245225263081049631 0.002774730265364384104 4.416262654344997005e-05
+8 1.5243589003230834323e-08 0.7813355837717117843
+0.000164587904124493665
+29.564459991843019537 -4.5824598513731222837 -0.5870359532621901577
+0.0004648344125208179762 0.0031282868879460171488 -7.5042704502708602616e-05

examples/symba_swifter_comparison/8pl_16tp_encounters/pl.swifter.in

@@ -0,0 +1,36 @@
+9
+0 0.00029591220819207775568
+0.0 0.0 0.0
+0.0 0.0 0.0
+1 4.9125474498983623693e-11 0.0014751238438755500459
+1.6306381826061645943e-05
+-0.065841771551149230746 0.30388831943526661838 0.030872485461978960153
+-0.033141166233939436947 -0.0049297226604189817514 0.0026371811668407158825
+2 7.243452483873646905e-10 0.006759080797928606587
+4.0453784346544178454e-05
+-0.65269716062695148917 -0.3065765656441301057 0.033456491497379246824
+0.008459831335658639026 -0.0184014319837384685 -0.0007407193515014080928
+3 8.9970113821660187435e-10 0.010044868190633438806
+4.25875607065040958e-05
+0.58046286084934750615 -0.8332000042504307258 3.7646553415201541957e-05
+0.013836557832279990782 0.009770187318278569788 -5.1179589633921335467e-07
+4 9.549535102761465607e-11 0.0072467082986392815006
+2.265740805092889601e-05
+-1.5891417403740180081 0.4938480736359250889 0.049330990309104823244
+-0.0036308073545784510204 -0.012168467501132099878 -0.00016594932370266260858
+5 2.825345908631354893e-07 0.3552707649709459117
+0.00046732617030490929307
+4.1148395833578952363 -2.8938323061728068453 -0.080043092204059404504
+0.0042549773877191511204 0.006534697671907701254 -0.00012233719535540690457
+6 8.459715183006415395e-08 0.43765596788571493287
+0.00038925687730393611812
+6.3589256477393849565 -7.653288021415167286 -0.12000977499446359442
+0.003985370599203661747 0.0035590677039893160206 -0.00022043610541731448703
+7 1.2920249163736673626e-08 0.46957663585116591335
+0.00016953449859497231466
+14.816779495279050138 13.049265812461410263 -0.14351615042000470668
+-0.0026245225263081049631 0.002774730265364384104 4.416262654344997005e-05
+8 1.5243589003230834323e-08 0.7813355837717117843
+0.000164587904124493665
+29.564459991843019537 -4.5824598513731222837 -0.5870359532621901577
+0.0004648344125208179762 0.0031282868879460171488 -7.5042704502708602616e-05

examples/symba_swifter_comparison/8pl_16tp_encounters/pl.swiftest.in

@@ -0,0 +1,33 @@
+8
+1 4.9125474498983623693e-11 0.0014751238438755500459
+1.6306381826061645943e-05
+-0.065841771551149230746 0.30388831943526661838 0.030872485461978960153
+-0.033141166233939436947 -0.0049297226604189817514 0.0026371811668407158825
+2 7.243452483873646905e-10 0.006759080797928606587
+4.0453784346544178454e-05
+-0.65269716062695148917 -0.3065765656441301057 0.033456491497379246824
+0.008459831335658639026 -0.0184014319837384685 -0.0007407193515014080928
+3 8.9970113821660187435e-10 0.010044868190633438806
+4.25875607065040958e-05
+0.58046286084934750615 -0.8332000042504307258 3.7646553415201541957e-05
+0.013836557832279990782 0.009770187318278569788 -5.1179589633921335467e-07
+4 9.549535102761465607e-11 0.0072467082986392815006
+2.265740805092889601e-05
+-1.5891417403740180081 0.4938480736359250889 0.049330990309104823244
+-0.0036308073545784510204 -0.012168467501132099878 -0.00016594932370266260858
+5 2.825345908631354893e-07 0.3552707649709459117
+0.00046732617030490929307
+4.1148395833578952363 -2.8938323061728068453 -0.080043092204059404504
+0.0042549773877191511204 0.006534697671907701254 -0.00012233719535540690457
+6 8.459715183006415395e-08 0.43765596788571493287
+0.00038925687730393611812
+6.3589256477393849565 -7.653288021415167286 -0.12000977499446359442
+0.003985370599203661747 0.0035590677039893160206 -0.00022043610541731448703
+7 1.2920249163736673626e-08 0.46957663585116591335
+0.00016953449859497231466
+14.816779495279050138 13.049265812461410263 -0.14351615042000470668
+-0.0026245225263081049631 0.002774730265364384104 4.416262654344997005e-05
+8 1.5243589003230834323e-08 0.7813355837717117843
+0.000164587904124493665
+29.564459991843019537 -4.5824598513731222837 -0.5870359532621901577
+0.0004648344125208179762 0.0031282868879460171488 -7.5042704502708602616e-05