diff --git a/examples/symba_gr_test/init_cond.py b/examples/symba_gr_test/init_cond.py index ca9e6726a..e59294190 100755 --- a/examples/symba_gr_test/init_cond.py +++ b/examples/symba_gr_test/init_cond.py @@ -29,7 +29,7 @@ sim.param['OUT_TYPE'] = "NETCDF_DOUBLE" sim.param['RHILL_PRESENT'] = "YES" sim.param['GR'] = 'YES' -sim.param['GMTINY'] = '1e-7' +sim.param['GMTINY'] = '1e-12' bodyid = { "Sun": 0, @@ -61,7 +61,7 @@ GMcb = sim.ds.isel(id=0)['Gmass'].values GU = swiftest.GC / (sim.param['DU2M']**3 / (sim.param['MU2KG'] * sim.param['TU2S']**2)) dens = 3000.0 / (sim.param['MU2KG'] / sim.param['DU2M']**3) # Assume a bulk density of 3 g/cm^3 -GM_pl = 2e-7 +GM_pl = 1e-11 M_pl = GM_pl / GU R_pl = (3 * M_pl / (4 * np.pi * dens))**(1.0 / 3.0) Rh_pl = Me_a * (GM_pl / (3 * GMcb))**(1.0/3.0) diff --git a/examples/symba_gr_test/param.swiftest.in b/examples/symba_gr_test/param.swiftest.in index 5a8aacd25..2a6312ce5 100644 --- a/examples/symba_gr_test/param.swiftest.in +++ b/examples/symba_gr_test/param.swiftest.in @@ -33,4 +33,4 @@ ENERGY NO GR YES INTERACTION_LOOPS ADAPTIVE ENCOUNTER_CHECK ADAPTIVE -GMTINY 1e-7 +GMTINY 1e-12 diff --git a/examples/symba_gr_test/pl.swiftest.in b/examples/symba_gr_test/pl.swiftest.in index 9d4ac811d..8c019ff5e 100644 --- a/examples/symba_gr_test/pl.swiftest.in +++ b/examples/symba_gr_test/pl.swiftest.in @@ -31,7 +31,7 @@ Neptune 0.0020336100526728302882 0.7749718408665500834 0.000164587904124493665 30.038959912561129073 0.008955570159296157365 1.7711193542961420899 131.82211597488270627 284.47484279411258967 308.45137222693909962 -Planetesimal 1.9999999999999999095e-07 0.0004609743061784924139 -6.2096743394181305454e-06 +Planetesimal 9.999999999999999395e-12 1.698243864025464137e-05 +2.2876635862715332271e-07 0.38709858430953941744 0.20562340108970009189 7.0033025080013837638 48.296118373786072198 29.20442403952453958 339.33948746828792764 diff --git a/src/symba/symba_step.f90 b/src/symba/symba_step.f90 index cd48667af..0d077d22d 100644 --- a/src/symba/symba_step.f90 +++ b/src/symba/symba_step.f90 @@ -250,6 +250,9 @@ module subroutine symba_step_reset_system(self, param) select type(tp => system%tp) class is (symba_tp) associate(npl => pl%nbody, ntp => tp%nbody) + nenc_old = system%plplenc_list%nenc + call system%plplenc_list%setup(0) + call system%plplcollision_list%setup(0) if (npl > 0) then pl%lcollision(1:npl) = .false. call pl%reset_kinship([(i, i=1, npl)]) @@ -261,24 +264,21 @@ module subroutine symba_step_reset_system(self, param) pl%lcollision(1:npl) = .false. pl%ldiscard(1:npl) = .false. pl%lmask(1:npl) = .true. - nenc_old = system%plplenc_list%nenc - call system%plplenc_list%setup(0) - call system%plplenc_list%setup(nenc_old) - system%plplenc_list%nenc = 0 - call system%plplcollision_list%setup(0) - call system%pl%set_renc(0) + call pl%set_renc(0) + call system%plplenc_list%setup(nenc_old) ! This resizes the pl-pl encounter list to be the same size as it was the last step, to decrease the number of potential resize operations that have to be one inside the step + system%plplenc_list%nenc = 0 ! Sets the true number of encounters back to 0 after resizing end if + nenc_old = system%pltpenc_list%nenc + call system%pltpenc_list%setup(0) if (ntp > 0) then tp%nplenc(1:ntp) = 0 tp%levelg(1:ntp) = -1 tp%levelm(1:ntp) = -1 tp%lmask(1:ntp) = .true. tp%ldiscard(1:ntp) = .false. - nenc_old = system%pltpenc_list%nenc - call system%pltpenc_list%setup(0) - call system%pltpenc_list%setup(nenc_old) - system%pltpenc_list%nenc = 0 + call system%pltpenc_list%setup(nenc_old)! This resizes the pl-tp encounter list to be the same size as it was the last step, to decrease the number of potential resize operations that have to be one inside the step + system%pltpenc_list%nenc = 0 ! Sets the true number of encounters back to 0 after resizing end if call system%pl_adds%setup(0, param)