From 6703c5e28230363e9d6a4c7c45e4d9c3fe1a5a1f Mon Sep 17 00:00:00 2001 From: David A Minton Date: Tue, 27 Jul 2021 14:03:50 -0400 Subject: [PATCH] Added Rhill to initial conditions generator --- .../9pl_18tp_encounters/cb.in | 2 +- .../9pl_18tp_encounters/cb.swiftest.in | 2 +- .../9pl_18tp_encounters/pl.in | 16 +++--- .../9pl_18tp_encounters/pl.swifter.in | 18 +++--- .../9pl_18tp_encounters/pl.swiftest.in | 16 +++--- python/swiftest/swiftest/init_cond.py | 4 +- python/swiftest/swiftest/io.py | 56 +++++++++++++------ 7 files changed, 68 insertions(+), 46 deletions(-) diff --git a/examples/symba_swifter_comparison/9pl_18tp_encounters/cb.in b/examples/symba_swifter_comparison/9pl_18tp_encounters/cb.in index 673a79459..81c636655 100644 --- a/examples/symba_swifter_comparison/9pl_18tp_encounters/cb.in +++ b/examples/symba_swifter_comparison/9pl_18tp_encounters/cb.in @@ -1,5 +1,5 @@ 0 -0.0002959122081920778 +0.00029591220819207774 0.004650467260962157 4.7535806948127355e-12 -2.2473967953572827e-18 diff --git a/examples/symba_swifter_comparison/9pl_18tp_encounters/cb.swiftest.in b/examples/symba_swifter_comparison/9pl_18tp_encounters/cb.swiftest.in index 673a79459..81c636655 100644 --- a/examples/symba_swifter_comparison/9pl_18tp_encounters/cb.swiftest.in +++ b/examples/symba_swifter_comparison/9pl_18tp_encounters/cb.swiftest.in @@ -1,5 +1,5 @@ 0 -0.0002959122081920778 +0.00029591220819207774 0.004650467260962157 4.7535806948127355e-12 -2.2473967953572827e-18 diff --git a/examples/symba_swifter_comparison/9pl_18tp_encounters/pl.in b/examples/symba_swifter_comparison/9pl_18tp_encounters/pl.in index d85d737f5..cd3c71538 100644 --- a/examples/symba_swifter_comparison/9pl_18tp_encounters/pl.in +++ b/examples/symba_swifter_comparison/9pl_18tp_encounters/pl.in @@ -1,33 +1,33 @@ 8 -1 4.9125474498983625056e-11 +1 4.9125474498983623693e-11 0.0014751237493860230134 1.6306381826061645943e-05 -0.032433320146471017464 0.30732647407569840814 0.0280888997405028297 -0.033622812072399158034 -0.0019305604712619159943 0.0029264451427202888868 -2 7.243452483873647106e-10 +2 7.243452483873646905e-10 0.006759082196678506012 4.0453784346544178454e-05 -0.6608991468450423623 -0.28805695486041710263 0.034183953683804932377 0.007943018642097033136 -0.018635382188272479886 -0.00071410720992500279457 -3 8.997011382166018993e-10 +3 8.9970113821660187435e-10 0.010044863223462002622 4.25875607065040958e-05 0.5665449483756358484 -0.84285201543201082597 3.8152874628327130158e-05 0.0139986033055793102076 0.009533392738922031109 -5.008237574040859916e-07 -4 9.549535102761465872e-11 +4 9.549535102761465607e-11 0.0072467110395904559343 2.265740805092889601e-05 -1.5854600237231359916 0.50600057977052448344 0.049495356229978339224 -0.0037325822023031099417 -0.0121364162752466003825 -0.00016278089870573419053 -5 2.8253459086313549713e-07 +5 2.825345908631354893e-07 0.35527078496549785303 0.00046732617030490929307 4.1105798235203270252 -2.9003636368897538489 -0.07992066204197022239 0.0042645403767569648595 0.006527961423420942065 -0.00012252307659855749943 -6 8.45971518300641563e-08 +6 8.459715183006415395e-08 0.43765573308845887078 0.00038925687730393611812 6.3549393159832749944 -7.6568459312514027815 -0.11978932080537739446 0.0039872926987931916337 0.0035567518157804990653 -0.00022047226166396519348 -7 1.2920249163736673984e-08 +7 1.2920249163736673626e-08 0.46957395507687206725 0.00016953449859497231466 14.81940372833062014 13.046490834898889943 -0.14356031024960910769 -0.002623943559850705834 0.002775224845039696818 4.4157032104701469965e-05 -8 1.5243589003230834746e-08 +8 1.5243589003230834323e-08 0.7813323455417420909 0.000164587904124493665 29.563994989459040141 -4.5855881090096284325 -0.5869609072731380994 0.00046517035659338968386 0.0031282283842968541462 -7.504927375628088796e-05 diff --git a/examples/symba_swifter_comparison/9pl_18tp_encounters/pl.swifter.in b/examples/symba_swifter_comparison/9pl_18tp_encounters/pl.swifter.in index ba4f47b86..e42b53b68 100644 --- a/examples/symba_swifter_comparison/9pl_18tp_encounters/pl.swifter.in +++ b/examples/symba_swifter_comparison/9pl_18tp_encounters/pl.swifter.in @@ -1,36 +1,36 @@ 9 -0 0.00029591220819207776388 +0 0.00029591220819207775568 0.0 0.0 0.0 0.0 0.0 0.0 -1 4.9125474498983625056e-11 0.3870979888395936208 +1 4.9125474498983623693e-11 0.0014751237493860230134 1.6306381826061645943e-05 -0.032433320146471017464 0.30732647407569840814 0.0280888997405028297 -0.033622812072399158034 -0.0019305604712619159943 0.0029264451427202888868 -2 7.243452483873647106e-10 0.7233257098327222634 +2 7.243452483873646905e-10 0.006759082196678506012 4.0453784346544178454e-05 -0.6608991468450423623 -0.28805695486041710263 0.034183953683804932377 0.007943018642097033136 -0.018635382188272479886 -0.00071410720992500279457 -3 8.997011382166018993e-10 1.000013777403204962 +3 8.9970113821660187435e-10 0.010044863223462002622 4.25875607065040958e-05 0.5665449483756358484 -0.84285201543201082597 3.8152874628327130158e-05 0.0139986033055793102076 0.009533392738922031109 -5.008237574040859916e-07 -4 9.549535102761465872e-11 1.523730802569037035 +4 9.549535102761465607e-11 0.0072467110395904559343 2.265740805092889601e-05 -1.5854600237231359916 0.50600057977052448344 0.049495356229978339224 -0.0037325822023031099417 -0.0121364162752466003825 -0.00016278089870573419053 -5 2.8253459086313549713e-07 5.2035172582389037643 +5 2.825345908631354893e-07 0.35527078496549785303 0.00046732617030490929307 4.1105798235203270252 -2.9003636368897538489 -0.07992066204197022239 0.0042645403767569648595 0.006527961423420942065 -0.00012252307659855749943 -6 8.45971518300641563e-08 9.581689367400855417 +6 8.459715183006415395e-08 0.43765573308845887078 0.00038925687730393611812 6.3549393159832749944 -7.6568459312514027815 -0.11978932080537739446 0.0039872926987931916337 0.0035567518157804990653 -0.00022047226166396519348 -7 1.2920249163736673984e-08 19.232766114773518495 +7 1.2920249163736673626e-08 0.46957395507687206725 0.00016953449859497231466 14.81940372833062014 13.046490834898889943 -0.14356031024960910769 -0.002623943559850705834 0.002775224845039696818 4.4157032104701469965e-05 -8 1.5243589003230834746e-08 30.285501702709421323 +8 1.5243589003230834323e-08 0.7813323455417420909 0.000164587904124493665 29.563994989459040141 -4.5855881090096284325 -0.5869609072731380994 0.00046517035659338968386 0.0031282283842968541462 -7.504927375628088796e-05 diff --git a/examples/symba_swifter_comparison/9pl_18tp_encounters/pl.swiftest.in b/examples/symba_swifter_comparison/9pl_18tp_encounters/pl.swiftest.in index d85d737f5..cd3c71538 100644 --- a/examples/symba_swifter_comparison/9pl_18tp_encounters/pl.swiftest.in +++ b/examples/symba_swifter_comparison/9pl_18tp_encounters/pl.swiftest.in @@ -1,33 +1,33 @@ 8 -1 4.9125474498983625056e-11 +1 4.9125474498983623693e-11 0.0014751237493860230134 1.6306381826061645943e-05 -0.032433320146471017464 0.30732647407569840814 0.0280888997405028297 -0.033622812072399158034 -0.0019305604712619159943 0.0029264451427202888868 -2 7.243452483873647106e-10 +2 7.243452483873646905e-10 0.006759082196678506012 4.0453784346544178454e-05 -0.6608991468450423623 -0.28805695486041710263 0.034183953683804932377 0.007943018642097033136 -0.018635382188272479886 -0.00071410720992500279457 -3 8.997011382166018993e-10 +3 8.9970113821660187435e-10 0.010044863223462002622 4.25875607065040958e-05 0.5665449483756358484 -0.84285201543201082597 3.8152874628327130158e-05 0.0139986033055793102076 0.009533392738922031109 -5.008237574040859916e-07 -4 9.549535102761465872e-11 +4 9.549535102761465607e-11 0.0072467110395904559343 2.265740805092889601e-05 -1.5854600237231359916 0.50600057977052448344 0.049495356229978339224 -0.0037325822023031099417 -0.0121364162752466003825 -0.00016278089870573419053 -5 2.8253459086313549713e-07 +5 2.825345908631354893e-07 0.35527078496549785303 0.00046732617030490929307 4.1105798235203270252 -2.9003636368897538489 -0.07992066204197022239 0.0042645403767569648595 0.006527961423420942065 -0.00012252307659855749943 -6 8.45971518300641563e-08 +6 8.459715183006415395e-08 0.43765573308845887078 0.00038925687730393611812 6.3549393159832749944 -7.6568459312514027815 -0.11978932080537739446 0.0039872926987931916337 0.0035567518157804990653 -0.00022047226166396519348 -7 1.2920249163736673984e-08 +7 1.2920249163736673626e-08 0.46957395507687206725 0.00016953449859497231466 14.81940372833062014 13.046490834898889943 -0.14356031024960910769 -0.002623943559850705834 0.002775224845039696818 4.4157032104701469965e-05 -8 1.5243589003230834746e-08 +8 1.5243589003230834323e-08 0.7813323455417420909 0.000164587904124493665 29.563994989459040141 -4.5855881090096284325 -0.5869609072731380994 0.00046517035659338968386 0.0031282283842968541462 -7.504927375628088796e-05 diff --git a/python/swiftest/swiftest/init_cond.py b/python/swiftest/swiftest/init_cond.py index 995b30194..f9a7378c0 100644 --- a/python/swiftest/swiftest/init_cond.py +++ b/python/swiftest/swiftest/init_cond.py @@ -196,8 +196,8 @@ def solar_system_horizons(plname, idval, param, ephemerides_start_date, ds): GMpl.append(GMcb[0] / MSun_over_Mpl[key]) # Generate planet value vectors if (param['RHILL_PRESENT'] == 'YES'): - Rhill.append(pldata[key].elements()['a'][0] * (3 * MSun_over_Mpl[key]) ** (-THIRDLONG)) - pvec = np.vstack([p1, p2, p3, p4, p5, p6, GMpl, Rpl, p7, p8, p9, p10, p11, p12, Rhill]) + Rhill.append(pldata[key].elements()['a'][0] * DCONV * (3 * MSun_over_Mpl[key]) ** (-THIRDLONG)) + pvec = np.vstack([p1, p2, p3, p4, p5, p6, GMpl, Rpl, Rhill, p7, p8, p9, p10, p11, p12]) else: pvec = np.vstack([p1, p2, p3, p4, p5, p6, GMpl, Rpl, p7, p8, p9, p10, p11, p12]) else: diff --git a/python/swiftest/swiftest/io.py b/python/swiftest/swiftest/io.py index aeaf742bb..2dd4ef7b3 100644 --- a/python/swiftest/swiftest/io.py +++ b/python/swiftest/swiftest/io.py @@ -678,7 +678,10 @@ def swiftest_xr2infile(ds, param, framenum=-1): print(pl.id.count().values, file=plfile) for i in pl.id: pli = pl.sel(id=i) - print(i.values, pli['Mass'].values, file=plfile) + if param['RHILL_PRESENT'] == 'YES': + print(i.values, pli['Mass'].values, pli['Rhill'].values, file=plfile) + else: + print(i.values, pli['Mass'].values, file=plfile) print(pli['Radius'].values, file=plfile) print(pli['px'].values, pli['py'].values, pli['pz'].values, file=plfile) print(pli['vx'].values, pli['vy'].values, pli['vz'].values, file=plfile) @@ -697,36 +700,55 @@ def swiftest_xr2infile(ds, param, framenum=-1): # Now make Swiftest files cbfile = FortranFile(param['CB_IN'], 'w') cbfile.write_record(cbid) - MSun = np.double(1.0) cbfile.write_record(np.double(GMSun)) - cbfile.write_record(np.double(rmin)) + cbfile.write_record(np.double(RSun)) cbfile.write_record(np.double(J2)) cbfile.write_record(np.double(J4)) cbfile.close() plfile = FortranFile(param['PL_IN'], 'w') - plfile.write_record(npl) + npl = pl.id.count().values + plid = pl.id.values + px = pl['px'].values + py = pl['py'].values + pz = pl['pz'].values + vx = pl['vx'].values + vy = pl['vy'].values + vz = pl['vz'].values + mass = pl['Mass'].values + radius = pl['Radius'].values + plfile.write_record(npl) plfile.write_record(plid) - plfile.write_record(p_pl[0]) - plfile.write_record(p_pl[1]) - plfile.write_record(p_pl[2]) - plfile.write_record(v_pl[0]) - plfile.write_record(v_pl[1]) - plfile.write_record(v_pl[2]) + plfile.write_record(px) + plfile.write_record(py) + plfile.write_record(pz) + plfile.write_record(vx) + plfile.write_record(vy) + plfile.write_record(vz) plfile.write_record(mass) + if param['RHILL_PRESENT'] == 'YES': + rhill = pl['Rhill'].values + plfile.write_record(rhill) plfile.write_record(radius) plfile.close() tpfile = FortranFile(param['TP_IN'], 'w') - ntp = 1 + ntp = tp.id.count().values + tpid = tp.id.values + px = tp['px'].values + py = tp['py'].values + pz = tp['pz'].values + vx = tp['vx'].values + vy = tp['vy'].values + vz = tp['vz'].values tpfile.write_record(ntp) tpfile.write_record(tpid) - tpfile.write_record(p_tp[0]) - tpfile.write_record(p_tp[1]) - tpfile.write_record(p_tp[2]) - tpfile.write_record(v_tp[0]) - tpfile.write_record(v_tp[1]) - tpfile.write_record(v_tp[2]) + tpfile.write_record(px) + tpfile.write_record(py) + tpfile.write_record(pz) + tpfile.write_record(vx) + tpfile.write_record(vy) + tpfile.write_record(vz) else: print(f"{param['IN_TYPE']} is an unknown file type")