RT Sub-cycling

doc/RTD/source/RadiativeTransfer/GEAR_RT.rst

@@ -22,7 +22,7 @@ Compiling for GEAR RT
 -   You need to choose a Riemann solver for the RT equations. You can choose
     between the ``GLF`` and ``HLL`` solver. For the time being, I recommend 
     sticking to the ``GLF`` solver as the ``HLL`` solver is more expensive,
-    but seemingly offers no advantage, although this remains to be comfirmed
+    but seemingly offers no advantage, although this remains to be confirmed
     in further testing.
 
 -   GEAR RT is only compatible with the Meshless Finite Volume scheme. You'll
@@ -59,6 +59,10 @@ You need to provide the following runtime parameters in the yaml file:
 
        stellar_spectrum_type: 0                           # Which radiation spectrum to use. 0: constant. 1: blackbody spectrum.
 
+   TimeIntegration:
+       max_nr_rt_subcycles: 128         # maximal number of RT subcycles per hydro step
+
+
 The ``photon_groups_Hz`` need to be ``N`` frequency edges (floats) to separate 
 the spectrum into ``N`` groups, where ``N`` is the same number you configured
 with using ``--with_rt=GEAR_N``. The edges are **lower** edges of the bins, and
@@ -111,6 +115,9 @@ to select between:
    e.g. equations 9 - 11 in `Rosdahl et al. 2013. 
    <https://ui.adsabs.harvard.edu/abs/2013MNRAS.436.2188R/abstract>`_
 
+Finally, you will also need to provide an upper threshold for the number of 
+RT-subcycles w.r.t. a single hydro step via ``TimeIntegration:max_nr_rt_subcycles``.
+For more details, refer to :ref:`the subcycling documentation <rt_subcycling>`.
 
 
 
@@ -393,7 +400,7 @@ useful:
    both methods can run on the exact same ICs.
 
 
-.. [#f2] For example, choosing cgs units as the interal units may lead to
+.. [#f2] For example, choosing cgs units as the internal units may lead to
    trouble with grackle. (Trouble like a gas at 10^6K without any heating
    sources heating up instead of cooling down.) The library is set up to work 
    with units geared towards cosmology. According to Britton Smith (private comm), 

doc/RTD/source/RadiativeTransfer/RT_subcycling.rst

+.. RT Subcycling
+    Mladen Ivkovic 07.2022
+
+.. _rt_subcycling:
+   
+RT Subcycling
+-------------
+
+.. warning::
+    The radiative transfer schemes are still in development and are not useable
+    at this moment. This page is currently a placeholder to document new
+    features and requirements as the code grows.
+
+
+SWIFT allows to sub-cycle the solution of radiative transfer steps (both 
+photon propagation and thermochemistry) with respect to the hydrodynamics
+time steps. Basically you can tell SWIFT to run up to X radiative transfer
+steps during a single hydrodynamics step for all particles in the simulation.
+The aim is to not waste time doing unnecessary hydrodynamics updates, which
+typically allow for much higher time steps compared to radiation due to the
+propagation speed of the respective advected quantity.
+
+You will need to provide an upper limit on how many RT subcycles per hydro
+step you want to allow. That is governed by the
+
+.. code:: yaml
+
+   TimeIntegration:
+       max_nr_rt_subcycles: 128         # maximal number of RT subcycles per hydro step
+
+parameter, which is mandatory for any RT runs. To turn off subcycling and 
+couple the radiative transfer and the hydrodynamics time steps one-to-one,
+set this parameter to either 0 or 1.
+
+Due to the discretization of individual particle time steps in time bins
+with a factor of 2 difference in time step size from a lower to a higher
+time bin, the ``max_nr_rt_subcycles`` parameter itself is required to be
+a power of 2 as well.
+
+Note that this parameter will set an upper limit to the number of subcycles
+per hydro step. If the ratio of hydro-to-RT time step is greater than what
+``max_nr_rt_subcycles`` allows for, then the hydro time step will be reduced
+to fit the maximal threshold. If it is smaller, the particle will simply do 
+fewer subcycles.
+
+.. warning::
+   Contrary to the documentation above, in the current implementation the 
+   ``max_nr_rt_subcycles`` parameters is abused as the fixed number of RT 
+   subcycles per hydro step. This means that RT time steps will be reduced
+   to lower values than necessary to fit the exact ratio should they be too
+   large to begin with.
+   Once the development advances, the behaviour will be set to be as 
+   documented above.
+

doc/RTD/source/RadiativeTransfer/index.rst

@@ -24,4 +24,6 @@ schemes.
    requirements
    GEAR_RT
    SPHM1_RT
+   RT_subcycling
+   RT_notes_for_developers
    additional_tools

doc/RTD/source/Task/current_dependencies.rst

@@ -68,3 +68,6 @@ As the hydrodynamics are described in :ref:`hydro`, we are only showing the grav
     In the second one, the gas particles tagged as "to be swallowed" are effectively swallowed.
     In the third one, the sink particles tagged as "to be swallowed" are effectively swallowed.
     This was done with SWIFT v0.9.0.
+
+For documentation on the radiative transfer tasking system, please refer to its 
+:ref:`own page <rt_task_system>`.

examples/RadiativeTransferTests/Advection_1D/rt_advection1D.yml

@@ -11,6 +11,7 @@ InternalUnitSystem:
 
 # Parameters governing the time integration
 TimeIntegration:
+  max_nr_rt_subcycles: 1
   time_begin: 0.    # The starting time of the simulation (in internal units).
   time_end:   4.e-1  # The end time of the simulation (in internal units).
   dt_min:     1.e-8 # The minimal time-step size of the simulation (in internal units).

examples/RadiativeTransferTests/Advection_1D/run.sh

@@ -18,7 +18,7 @@ fi
     --stars \
     --feedback \
     --external-gravity \
-    -e \
+    --fpe \
     ./rt_advection1D.yml 2>&1 | tee output.log
 
 python3 ./plotSolution.py

examples/RadiativeTransferTests/Advection_2D/rt_advection2D.yml

@@ -11,6 +11,7 @@ InternalUnitSystem:
 
 # Parameters governing the time integration
 TimeIntegration:
+  max_nr_rt_subcycles: 1
   time_begin: 0.     # The starting time of the simulation (in internal units).
   time_end:   3.4e-1 # The end time of the simulation (in internal units).
   dt_min:     1.e-08 # The minimal time-step size of the simulation (in internal units).

examples/RadiativeTransferTests/Advection_2D/run.sh

@@ -25,7 +25,7 @@ fi
     --stars \
     --feedback \
     --external-gravity \
-    -e \
+    --fpe \
     ./rt_advection2D.yml 2>&1 | tee output.log
 
 python3 ./plotSolution.py

examples/RadiativeTransferTests/CoolingTest/rt_cooling_test.yml

@@ -12,6 +12,7 @@ InternalUnitSystem:
 
 # Parameters governing the time integration
 TimeIntegration:
+  max_nr_rt_subcycles: 1
   time_begin: 0.     # The starting time of the simulation (in internal units).
   time_end:   0.100  # The end time of the simulation (in internal units).
   dt_min:     1.e-8  # The minimal time-step size of the simulation (in internal units).

examples/RadiativeTransferTests/HeatingTest/rt_heating_test.yml

@@ -12,6 +12,7 @@ InternalUnitSystem:
 
 # Parameters governing the time integration
 TimeIntegration:
+  max_nr_rt_subcycles: 1
   time_begin: 0.     # The starting time of the simulation (in internal units).
   dt_min:     1.0208453377e-08  # 0.01 yr
   dt_max:     0.0010208453      # 1 kyr

examples/RadiativeTransferTests/IonMassFractionAdvectionTest_2D/advect_ions.yml

@@ -11,6 +11,7 @@ InternalUnitSystem:
 
 # Parameters governing the time integration
 TimeIntegration:
+  max_nr_rt_subcycles: 1
   time_begin: 0.    # The starting time of the simulation (in internal units).
   time_end:   0.01   # end time: radiation reaches edge of box
   dt_min:     1.e-12 # The minimal time-step size of the simulation (in internal units).

examples/RadiativeTransferTests/IonizationEquilibriumICSetupTest/ion_equil.yml

@@ -11,6 +11,7 @@ InternalUnitSystem:
 
 # Parameters governing the time integration
 TimeIntegration:
+  max_nr_rt_subcycles: 1
   time_begin: 0.    # The starting time of the simulation (in internal units).
   time_end:   6.250000e-03  # The end time of the simulation (in internal units).
   dt_min:     1.e-6 # The minimal time-step size of the simulation (in internal units).

examples/RadiativeTransferTests/RandomizedBox_3D/README

@@ -15,3 +15,12 @@ check that everything's right.
 
 Compile swift with:
     --with-rt=GEAR_1 --with-rt-riemann-solver=GLF --with-hydro-dimension=3 --with-hydro=gizmo-mfv --with-riemann-solver=hllc --with-stars=GEAR --with-feedback=none --enable-debugging-checks --with-grackle=$GRACKLE_ROOT
+
+
+
+Some specific notes for debugging:
+
+-   when running on 2 MPI ranks (at least in the first few steps) cells 27 and
+    151 have proxies on the other rank. However, those are made for the gravity
+    tasks, not for the RT tasks. So their foreign counterparts will have no RT
+    tasks at all.

examples/RadiativeTransferTests/RandomizedBox_3D/randomized-rt.yml

@@ -11,9 +11,10 @@ InternalUnitSystem:
 
 # Parameters governing the time integration
 TimeIntegration:
+  max_nr_rt_subcycles: 4
   time_begin: 0.    # The starting time of the simulation (in internal units).
   time_end:   0.1   # The end time of the simulation (in internal units).
-  dt_min:     1.e-8 # The minimal time-step size of the simulation (in internal units).
+  dt_min:     1.e-6 # The minimal time-step size of the simulation (in internal units).
   dt_max:     1.e-3 # The maximal time-step size of the simulation (in internal units).
 
 # Parameters governing the snapshots
@@ -49,7 +50,7 @@ Scheduler:
 Restarts:
   delta_hours:        72        # (Optional) decimal hours between dumps of restart files.
   enable:             1          # (Optional) whether to enable dumping restarts at fixed intervals.
-  stop_steps:         500        # (Optional) how many steps to process before checking if the <subdir>/stop file exists. When present the application will attempt to exit early, dumping restart files first.
+  stop_steps:         128        # (Optional) how many steps to process before checking if the <subdir>/stop file exists. When present the application will attempt to exit early, dumping restart files first.
 
 # Parameters for the self-gravity scheme
 Gravity:

examples/RadiativeTransferTests/RandomizedBox_3D/run.sh

@@ -9,12 +9,26 @@ if [ ! -f 'randomized-sine.hdf5' ]; then
     python3 makeIC.py
 fi
 
+# use cmdline args as shortcut to run with debugger/MPI
+# -g run with gdb
+# -m run with MPI
+# -mg run with MPI and gdb in individual xterm windows
+# -ml run with MPI and write individual output file per MPI rank
 cmd=../../../swift
 if [ $# -gt 0 ]; then
     case "$1" in 
-    g | gdb)
+    -g | g | gdb)
         cmd='gdb --args ../../../swift'
         ;;
+    -m | m | mpi)
+        cmd='mpirun -n 3 ../../../swift_mpi' 
+        ;;
+    -mg | -gm | gm | mg | gmpi | gdbmpi )
+        cmd='mpirun -n 3 xterm -e gdb -ex run --args ../../../swift_mpi'
+        ;;
+    -ml | ml | lm | mpilog | logmpi)
+        cmd='mpirun -n 2 --output-filename individual_rank_output --merge-stderr-to-stdout ../../../swift_mpi'
+        ;;
     *)
         echo unknown cmdline param, running without gdb
         ;;
@@ -24,8 +38,8 @@ fi
 # Run SWIFT with RT
 $cmd \
     --hydro \
-    --threads=9 \
-    --verbose=0  \
+    --threads=3 \
+    --verbose=0 \
     --radiation \
     --self-gravity \
     --stars \

examples/RadiativeTransferTests/StromgrenSphere_2D/plotSolution.py

@@ -308,3 +308,4 @@ if __name__ == "__main__":
 
     for f in snaplist:
         plot_result(f)
+        gc.collect()