Commit bca94340 authored by Matthieu Schaller's avatar Matthieu Schaller
Browse files

Merge branch 'master' into FOF

parents f8a0b838 080b8536
...@@ -49,9 +49,12 @@ examples/*/*.mpg ...@@ -49,9 +49,12 @@ examples/*/*.mpg
examples/*/*/gravity_checks_*.dat examples/*/*/gravity_checks_*.dat
examples/*/*/coolingtables.tar.gz examples/*/*/coolingtables.tar.gz
examples/*/*/coolingtables examples/*/*/coolingtables
examples/*/*/yieldtables.tar.gz
examples/*/*/yieldtables
examples/Cooling/CoolingRates/cooling_rates examples/Cooling/CoolingRates/cooling_rates
examples/Cooling/CoolingRates/cooling_element_*.dat examples/Cooling/CoolingRates/cooling_element_*.dat
examples/Cooling/CoolingRates/cooling_output.dat examples/Cooling/CoolingRates/cooling_output.dat
examples/SubgridTests/StellarEvolution/StellarEvolutionSolution*
tests/testActivePair tests/testActivePair
tests/testActivePair.sh tests/testActivePair.sh
...@@ -110,6 +113,7 @@ tests/testMaths ...@@ -110,6 +113,7 @@ tests/testMaths
tests/testRandom tests/testRandom
tests/testThreadpool tests/testThreadpool
tests/testParser tests/testParser
tests/testFeedback
tests/parser_output.yml tests/parser_output.yml
tests/testPeriodicBC.sh tests/testPeriodicBC.sh
tests/testPeriodicBCPerturbed.sh tests/testPeriodicBCPerturbed.sh
......
...@@ -35,6 +35,7 @@ Parameters: ...@@ -35,6 +35,7 @@ Parameters:
-M, --multipole-reconstruction Reconstruct the multipoles every time-step. -M, --multipole-reconstruction Reconstruct the multipoles every time-step.
-s, --hydro Run with hydrodynamics. -s, --hydro Run with hydrodynamics.
-S, --stars Run with stars. -S, --stars Run with stars.
-B, --black-holes Run with black holes.
-x, --velociraptor Run with structure finding. -x, --velociraptor Run with structure finding.
--limiter Run with time-step limiter. --limiter Run with time-step limiter.
......
...@@ -83,6 +83,7 @@ Parameters: ...@@ -83,6 +83,7 @@ Parameters:
-M, --multipole-reconstruction Reconstruct the multipoles every time-step. -M, --multipole-reconstruction Reconstruct the multipoles every time-step.
-s, --hydro Run with hydrodynamics. -s, --hydro Run with hydrodynamics.
-S, --stars Run with stars. -S, --stars Run with stars.
-B, --black-holes Run with black holes.
-x, --velociraptor Run with structure finding. -x, --velociraptor Run with structure finding.
--limiter Run with time-step limiter. --limiter Run with time-step limiter.
......
...@@ -1329,7 +1329,8 @@ case "$with_subgrid" in ...@@ -1329,7 +1329,8 @@ case "$with_subgrid" in
with_subgrid_entropy_floor=none with_subgrid_entropy_floor=none
with_subgrid_stars=GEAR with_subgrid_stars=GEAR
with_subgrid_star_formation=GEAR with_subgrid_star_formation=GEAR
with_subgrid_feedback=thermal with_subgrid_feedback=none
with_subgrid_black_holes=none
;; ;;
EAGLE) EAGLE)
with_subgrid_cooling=EAGLE with_subgrid_cooling=EAGLE
...@@ -1338,7 +1339,8 @@ case "$with_subgrid" in ...@@ -1338,7 +1339,8 @@ case "$with_subgrid" in
with_subgrid_entropy_floor=EAGLE with_subgrid_entropy_floor=EAGLE
with_subgrid_stars=EAGLE with_subgrid_stars=EAGLE
with_subgrid_star_formation=EAGLE with_subgrid_star_formation=EAGLE
with_subgrid_feedback=none with_subgrid_feedback=EAGLE
with_subgrid_black_holes=EAGLE
;; ;;
*) *)
AC_MSG_ERROR([Unknown subgrid choice: $with_subgrid]) AC_MSG_ERROR([Unknown subgrid choice: $with_subgrid])
...@@ -1369,7 +1371,7 @@ esac ...@@ -1369,7 +1371,7 @@ esac
# Hydro scheme. # Hydro scheme.
AC_ARG_WITH([hydro], AC_ARG_WITH([hydro],
[AS_HELP_STRING([--with-hydro=<scheme>], [AS_HELP_STRING([--with-hydro=<scheme>],
[Hydro dynamics to use @<:@gadget2, minimal, pressure-entropy, pressure-energy, pressure-energy-monaghan, default, gizmo-mfv, gizmo-mfm, shadowfax, planetary, anarchy-pu debug default: gadget2@:>@] [Hydro dynamics to use @<:@gadget2, minimal, pressure-entropy, pressure-energy, pressure-energy-monaghan, default, gizmo-mfv, gizmo-mfm, shadowfax, planetary, anarchy-pu default: gadget2@:>@]
)], )],
[with_hydro="$withval"], [with_hydro="$withval"],
[with_hydro="gadget2"] [with_hydro="gadget2"]
...@@ -1728,7 +1730,7 @@ case "$with_stars" in ...@@ -1728,7 +1730,7 @@ case "$with_stars" in
AC_DEFINE([STARS_GEAR], [1], [GEAR stellar model]) AC_DEFINE([STARS_GEAR], [1], [GEAR stellar model])
;; ;;
none) none)
AC_DEFINE([STARS_NONE], [1], [None stellar model]) AC_DEFINE([STARS_NONE], [1], [Basic stellar model])
;; ;;
*) *)
...@@ -1739,7 +1741,7 @@ esac ...@@ -1739,7 +1741,7 @@ esac
# Feedback model # Feedback model
AC_ARG_WITH([feedback], AC_ARG_WITH([feedback],
[AS_HELP_STRING([--with-feedback=<model>], [AS_HELP_STRING([--with-feedback=<model>],
[Feedback model to use @<:@none, thermal, debug default: none@:>@] [Feedback model to use @<:@none, EAGLE, debug default: none@:>@]
)], )],
[with_feedback="$withval"], [with_feedback="$withval"],
[with_feedback="none"] [with_feedback="none"]
...@@ -1754,10 +1756,11 @@ if test "$with_subgrid" != "none"; then ...@@ -1754,10 +1756,11 @@ if test "$with_subgrid" != "none"; then
fi fi
case "$with_feedback" in case "$with_feedback" in
thermal) EAGLE)
AC_DEFINE([FEEDBACK_THERMAL], [1], [Thermal Blastwave]) AC_DEFINE([FEEDBACK_EAGLE], [1], [EAGLE stellar feedback and evolution model])
;; ;;
none) none)
AC_DEFINE([FEEDBACK_NONE], [1], [No feedback])
;; ;;
*) *)
...@@ -1765,10 +1768,39 @@ case "$with_feedback" in ...@@ -1765,10 +1768,39 @@ case "$with_feedback" in
;; ;;
esac esac
# Black hole model.
AC_ARG_WITH([black-holes],
[AS_HELP_STRING([--with-black-holes=<model>],
[Black holes model to use @<:@none, default: none@:>@]
)],
[with_black_holes="$withval"],
[with_black_holes="none"]
)
if test "$with_subgrid" != "none"; then
if test "$with_black_holes" != "none"; then
AC_MSG_ERROR([Cannot provide with-subgrid and with-black-holes together])
else
with_black_holes="$with_subgrid_black_holes"
fi
fi
case "$with_black_holes" in
none)
AC_DEFINE([BLACK_HOLES_NONE], [1], [No black hole model])
;;
EAGLE)
AC_DEFINE([BLACK_HOLES_EAGLE], [1], [EAGLE black hole model])
;;
*)
AC_MSG_ERROR([Unknown black-hole model: $with_black_holes])
;;
esac
# External potential # External potential
AC_ARG_WITH([ext-potential], AC_ARG_WITH([ext-potential],
[AS_HELP_STRING([--with-ext-potential=<pot>], [AS_HELP_STRING([--with-ext-potential=<pot>],
[external potential @<:@none, point-mass, point-mass-ring, point-mass-softened, isothermal, softened-isothermal, nfw, hernquist, disc-patch, sine-wave, default: none@:>@] [external potential @<:@none, point-mass, point-mass-ring, point-mass-softened, isothermal, nfw, hernquist, disc-patch, sine-wave, default: none@:>@]
)], )],
[with_potential="$withval"], [with_potential="$withval"],
[with_potential="none"] [with_potential="none"]
...@@ -1886,6 +1918,9 @@ AM_CONDITIONAL([HAVE_DOXYGEN], [test "$ac_cv_path_ac_pt_DX_DOXYGEN" != ""]) ...@@ -1886,6 +1918,9 @@ AM_CONDITIONAL([HAVE_DOXYGEN], [test "$ac_cv_path_ac_pt_DX_DOXYGEN" != ""])
# Check if using EAGLE cooling # Check if using EAGLE cooling
AM_CONDITIONAL([HAVEEAGLECOOLING], [test $with_cooling = "EAGLE"]) AM_CONDITIONAL([HAVEEAGLECOOLING], [test $with_cooling = "EAGLE"])
# Check if using EAGLE feedback
AM_CONDITIONAL([HAVEEAGLEFEEDBACK], [test $with_feedback = "EAGLE"])
# Handle .in files. # Handle .in files.
AC_CONFIG_FILES([Makefile src/Makefile examples/Makefile examples/Cooling/CoolingRates/Makefile doc/Makefile doc/Doxyfile tests/Makefile]) AC_CONFIG_FILES([Makefile src/Makefile examples/Makefile examples/Cooling/CoolingRates/Makefile doc/Makefile doc/Doxyfile tests/Makefile])
AC_CONFIG_FILES([argparse/Makefile tools/Makefile]) AC_CONFIG_FILES([argparse/Makefile tools/Makefile])
...@@ -1962,7 +1997,8 @@ AC_MSG_RESULT([ ...@@ -1962,7 +1997,8 @@ AC_MSG_RESULT([
Tracers : $with_tracers Tracers : $with_tracers
Stellar model : $with_stars Stellar model : $with_stars
Star formation model : $with_star_formation Star formation model : $with_star_formation
Feedback model : $with_feedback Star feedback model : $with_feedback
Black holes model : $with_black_holes
Individual timers : $enable_timers Individual timers : $enable_timers
Task debugging : $enable_task_debugging Task debugging : $enable_task_debugging
......
...@@ -763,7 +763,6 @@ INPUT = @top_srcdir@ @top_srcdir@/src @top_srcdir@/tests @top_ ...@@ -763,7 +763,6 @@ INPUT = @top_srcdir@ @top_srcdir@/src @top_srcdir@/tests @top_
INPUT += @top_srcdir@/src/hydro/Minimal INPUT += @top_srcdir@/src/hydro/Minimal
INPUT += @top_srcdir@/src/hydro/Gadget2 INPUT += @top_srcdir@/src/hydro/Gadget2
INPUT += @top_srcdir@/src/gravity/Default INPUT += @top_srcdir@/src/gravity/Default
INPUT += @top_srcdir@/src/stars/Default
INPUT += @top_srcdir@/src/riemann INPUT += @top_srcdir@/src/riemann
INPUT += @top_srcdir@/src/potential/point_mass INPUT += @top_srcdir@/src/potential/point_mass
INPUT += @top_srcdir@/src/equation_of_state/ideal_gas INPUT += @top_srcdir@/src/equation_of_state/ideal_gas
...@@ -776,6 +775,8 @@ INPUT += @top_srcdir@/src/entropy_floor/EAGLE ...@@ -776,6 +775,8 @@ INPUT += @top_srcdir@/src/entropy_floor/EAGLE
INPUT += @top_srcdir@/src/star_formation/EAGLE INPUT += @top_srcdir@/src/star_formation/EAGLE
INPUT += @top_srcdir@/src/tracers/EAGLE INPUT += @top_srcdir@/src/tracers/EAGLE
INPUT += @top_srcdir@/src/stars/EAGLE INPUT += @top_srcdir@/src/stars/EAGLE
INPUT += @top_srcdir@/src/feedback/EAGLE
INPUT += @top_srcdir@/src/black_holes/EAGLE
# This tag can be used to specify the character encoding of the source files # This tag can be used to specify the character encoding of the source files
# that doxygen parses. Internally doxygen uses the UTF-8 encoding. Doxygen uses # that doxygen parses. Internally doxygen uses the UTF-8 encoding. Doxygen uses
...@@ -810,7 +811,7 @@ RECURSIVE = NO ...@@ -810,7 +811,7 @@ RECURSIVE = NO
# Note that relative paths are relative to the directory from which doxygen is # Note that relative paths are relative to the directory from which doxygen is
# run. # run.
EXCLUDE = EXCLUDE = @top_srcdir@/src/cooling/EAGLE/newton_cooling.c
# The EXCLUDE_SYMLINKS tag can be used to select whether or not files or # The EXCLUDE_SYMLINKS tag can be used to select whether or not files or
# directories that are symbolic links (a Unix file system feature) are excluded # directories that are symbolic links (a Unix file system feature) are excluded
......
...@@ -30,6 +30,7 @@ can be found by typing ``./swift -h``:: ...@@ -30,6 +30,7 @@ can be found by typing ``./swift -h``::
-M, --multipole-reconstruction Reconstruct the multipoles every time-step. -M, --multipole-reconstruction Reconstruct the multipoles every time-step.
-s, --hydro Run with hydrodynamics. -s, --hydro Run with hydrodynamics.
-S, --stars Run with stars. -S, --stars Run with stars.
-B, --black-holes Run with black holes.
-x, --velociraptor Run with structure finding. -x, --velociraptor Run with structure finding.
--limiter Run with time-step limiter. --limiter Run with time-step limiter.
......
...@@ -656,6 +656,17 @@ which stops these from being done at the scale of the leaf cells, of which ...@@ -656,6 +656,17 @@ which stops these from being done at the scale of the leaf cells, of which
there can be a large number. In this case cells with gravity tasks must be at there can be a large number. In this case cells with gravity tasks must be at
least 4 levels above the leaf cells (when possible). least 4 levels above the leaf cells (when possible).
To control the depth at which the ghost tasks are placed, there are
two parameters (one for the gas, one for the stars). These specify the
maximum number of particles allowed in such a task before splitting
into finer ones. These parameters are:
.. code:: YAML
engine_max_parts_per_ghost: 1000
engine_max_sparts_per_ghost: 1000
Extra space is required when particles are created in the system (to the time Extra space is required when particles are created in the system (to the time
of the next rebuild). These are controlled by: of the next rebuild). These are controlled by:
......
...@@ -20,7 +20,10 @@ particles. Two floors are used in conjonction. Both are implemented as ...@@ -20,7 +20,10 @@ particles. Two floors are used in conjonction. Both are implemented as
polytropic "equations of states":math:`P = P_c polytropic "equations of states":math:`P = P_c
\left(\rho/\rho_c\right)^\gamma` (all done in physical coordinates), with \left(\rho/\rho_c\right)^\gamma` (all done in physical coordinates), with
the constants derived from the user input given in terms of temperature and the constants derived from the user input given in terms of temperature and
Hydrogen number density. Hydrogen number density. The code computing the entropy floor
is located in the directory ``src/entropy_floor/EAGLE/`` and the floor
is applied in the drift and kick operations of the hydro scheme. It is
also used in some of the other subgrid schemes.
The first limit, labelled as ``Cool``, is typically used to prevent The first limit, labelled as ``Cool``, is typically used to prevent
low-density high-metallicity particles to cool below the warm phase because low-density high-metallicity particles to cool below the warm phase because
...@@ -57,13 +60,19 @@ critical density at redshift zero [#f1]_, and :math:`\rho_{\rm com}` the ...@@ -57,13 +60,19 @@ critical density at redshift zero [#f1]_, and :math:`\rho_{\rm com}` the
gas co-moving density. Typical values for :math:`\Delta_{\rm floor}` are of gas co-moving density. Typical values for :math:`\Delta_{\rm floor}` are of
order 10. order 10.
The model is governed by 4 parameters for each of the two The model is governed by 4 parameters for each of the two limits. These are
limits. These are given in the ``EAGLEEntropyFloor`` section of the given in the ``EAGLEEntropyFloor`` section of the YAML file. The parameters
YAML file. The parameters are the Hydrogen number density (in are the Hydrogen number density (in :math:`cm^{-3}`) and temperature (in
:math:`cm^{-3}`) and temperature (in :math:`K`) of the anchor point of :math:`K`) of the anchor point of each floor as well as the power-law slope
each floor as well as the power-law slope of each floor and the of each floor and the minimal over-density required to apply the
minimal over-density required to apply the limit. For a normal limit. Note that, even though the anchor points are given in terms of
EAGLE run, that section of the parameter file reads: temperatures, the slopes are expressed using a power-law in terms of
entropy and *not* in terms of temperature. For a slope of :math:`\gamma` in
the parameter file, the temperature as a function of density will be
limited to be above a power-law with slope :math:`\gamma - 1` (as shown on
the figure above).
For a normal EAGLE run, that section of the parameter file reads:
.. code:: YAML .. code:: YAML
...@@ -87,6 +96,10 @@ floor by a factor :math:`\frac{\mu_{\rm neutral}}{\mu_{ionised}} ...@@ -87,6 +96,10 @@ floor by a factor :math:`\frac{\mu_{\rm neutral}}{\mu_{ionised}}
\approx \frac{1.22}{0.59} \approx 2` due to the different ionisation \approx \frac{1.22}{0.59} \approx 2` due to the different ionisation
states of the gas. states of the gas.
Recall that we additionally impose an absolute minium temperature at all
densities with a value provided in the :ref:`Parameters_SPH` section of the parameter
file. This minimal temperature is typically set to 100 Kelvin.
Note that the model only makes sense if the ``Cool`` threshold is at a lower Note that the model only makes sense if the ``Cool`` threshold is at a lower
density than the ``Jeans`` threshold. density than the ``Jeans`` threshold.
...@@ -457,8 +470,46 @@ Stellar enrichment: Wiersma+2009b ...@@ -457,8 +470,46 @@ Stellar enrichment: Wiersma+2009b
.. _EAGLE_feedback: .. _EAGLE_feedback:
Supernova feedback: Dalla Vecchia+2012 Supernova feedback: Dalla Vecchia+2012 & Schaye+2015
~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
.. code:: YAML
# EAGLE stellar enrichment and feedback model
EAGLEFeedback:
use_SNII_feedback: 1 # Global switch for SNII thermal (stochastic) feedback.
use_SNIa_feedback: 1 # Global switch for SNIa thermal (continuous) feedback.
use_AGB_enrichment: 1 # Global switch for enrichement from AGB stars.
use_SNII_enrichment: 1 # Global switch for enrichement from SNII stars.
use_SNIa_enrichment: 1 # Global switch for enrichement from SNIa stars.
filename: ./yieldtables/ # Path to the directory containing the EAGLE yield tables.
IMF_min_mass_Msun: 0.1 # Minimal stellar mass considered for the Chabrier IMF in solar masses.
IMF_max_mass_Msun: 100.0 # Maximal stellar mass considered for the Chabrier IMF in solar masses.
SNII_min_mass_Msun: 6.0 # Minimal mass considered for SNII feedback (not SNII enrichment!) in solar masses.
SNII_max_mass_Msun: 100.0 # Maximal mass considered for SNII feedback (not SNII enrichment!) in solar masses.
SNII_wind_delay_Gyr: 0.03 # Time in Gyr between a star's birth and the SNII thermal feedback event.
SNII_delta_T_K: 3.16228e7 # Change in temperature to apply to the gas particle in a SNII thermal feedback event in Kelvin.
SNII_energy_erg: 1.0e51 # Energy of one SNII explosion in ergs.
SNII_energy_fraction_min: 3.0 # Maximal fraction of energy applied in a SNII feedback event.
SNII_energy_fraction_max: 0.3 # Minimal fraction of energy applied in a SNII feedback event.
SNII_energy_fraction_Z_0: 0.0012663729 # Pivot point for the metallicity dependance of the SNII energy fraction (metal mass fraction).
SNII_energy_fraction_n_0_H_p_cm3: 0.67 # Pivot point for the birth density dependance of the SNII energy fraction in cm^-3.
SNII_energy_fraction_n_Z: 0.8686 # Power-law for the metallicity dependance of the SNII energy fraction.
SNII_energy_fraction_n_n: 0.8686 # Power-law for the birth density dependance of the SNII energy fraction.
SNIa_max_mass_Msun: 8.0 # Maximal mass considered for SNIa feedback and enrichment in solar masses.
SNIa_timescale_Gyr: 2.0 # Time-scale of the exponential decay of the SNIa rates in Gyr.
SNIa_efficiency_p_Msun: 0.002 # Normalisation of the SNIa rates in inverse solar masses.
SNIa_energy_erg: 1.0e51 # Energy of one SNIa explosion in ergs.
AGB_ejecta_velocity_km_p_s: 10.0 # Velocity of the AGB ejectas in km/s.
SNII_yield_factor_Hydrogen: 1.0 # (Optional) Correction factor to apply to the Hydrogen yield from the SNII channel.
SNII_yield_factor_Helium: 1.0 # (Optional) Correction factor to apply to the Helium yield from the SNII channel.
SNII_yield_factor_Carbon: 0.5 # (Optional) Correction factor to apply to the Carbon yield from the SNII channel.
SNII_yield_factor_Nitrogen: 1.0 # (Optional) Correction factor to apply to the Nitrogen yield from the SNII channel.
SNII_yield_factor_Oxygen: 1.0 # (Optional) Correction factor to apply to the Oxygen yield from the SNII channel.
SNII_yield_factor_Neon: 1.0 # (Optional) Correction factor to apply to the Neon yield from the SNII channel.
SNII_yield_factor_Magnesium: 2.0 # (Optional) Correction factor to apply to the Magnesium yield from the SNII channel.
SNII_yield_factor_Silicon: 1.0 # (Optional) Correction factor to apply to the Silicon yield from the SNII channel.
SNII_yield_factor_Iron: 0.5 # (Optional) Correction factor to apply to the Iron yield from the SNII channel.
.. _EAGLE_black_hole_seeding: .. _EAGLE_black_hole_seeding:
......
...@@ -52,21 +52,21 @@ plot([1e-1, 1e-1], [20, 4000], "k:", lw=0.6) ...@@ -52,21 +52,21 @@ plot([1e-1, 1e-1], [20, 4000], "k:", lw=0.6)
plot([3e-6, 3e-4], [28000, 28000], "k--", lw=0.6) plot([3e-6, 3e-4], [28000, 28000], "k--", lw=0.6)
text( text(
3e-6, 1e-6,
22500, 22500,
"$n_{\\rm H}$^Cool_gamma_effective", "$n_{\\rm H}$^Cool_gamma_effective - 1",
va="top", va="top",
fontsize=6, fontsize=5.5,
family="monospace", family="monospace",
) )
plot([3e-1, 3e1], [15000.0, 15000.0 * 10.0 ** (2.0 / 3.0)], "k--", lw=0.6) plot([3e-1, 3e1], [15000.0, 15000.0 * 10.0 ** (2.0 / 3.0)], "k--", lw=0.6)
text( text(
3e-1, 1.2e-1,
200000, 190000,
"$n_{\\rm H}$^Jeans_gamma_effective", "$n_{\\rm H}$^Jeans_gamma_effective - 1",
va="top", va="top",
rotation=43, rotation=43,
fontsize=6, fontsize=5.5,
family="monospace", family="monospace",
) )
text( text(
...@@ -76,7 +76,7 @@ text( ...@@ -76,7 +76,7 @@ text(
rotation=90, rotation=90,
va="bottom", va="bottom",
ha="right", ha="right",
fontsize=6, fontsize=5.5,
family="monospace", family="monospace",
) )
text( text(
...@@ -89,9 +89,9 @@ text( ...@@ -89,9 +89,9 @@ text(
fontsize=5.5, fontsize=5.5,
family="monospace", family="monospace",
) )
text(5e-8, 8800, "Cool_temperature_norm_K", va="bottom", fontsize=6, family="monospace") text(5e-8, 8800, "Cool_temperature_norm_K", va="bottom", fontsize=5.5, family="monospace")
text( text(
5e-8, 4400, "Jeans_temperature_norm_K", va="bottom", fontsize=6, family="monospace" 5e-8, 4400, "Jeans_temperature_norm_K", va="bottom", fontsize=5.5, family="monospace"
) )
fill_between([1e-5, 1e5], [10, 10], [8000, 8000], color="0.9") fill_between([1e-5, 1e5], [10, 10], [8000, 8000], color="0.9")
fill_between([1e-1, 1e5], [4000, 400000], color="0.9") fill_between([1e-1, 1e5], [4000, 400000], color="0.9")
......
...@@ -10,7 +10,7 @@ This section of the documentation includes information on the task system ...@@ -10,7 +10,7 @@ This section of the documentation includes information on the task system
available in SWIFT, as well as how to implement your own task. available in SWIFT, as well as how to implement your own task.
SWIFT can produce a graph containing all the dependencies. SWIFT can produce a graph containing all the dependencies.
Everything is described in :ref:`_analysistools`. Everything is described in :ref:`Analysis_Tools`.
.. toctree:: .. toctree::
......
...@@ -104,9 +104,9 @@ INLINE static double eagle_print_metal_cooling_rate( ...@@ -104,9 +104,9 @@ INLINE static double eagle_print_metal_cooling_rate(
/* calculate cooling rates */ /* calculate cooling rates */
for (int j = 0; j < eagle_cooling_N_metal + 2; j++) element_lambda[j] = 0.0; for (int j = 0; j < eagle_cooling_N_metal + 2; j++) element_lambda[j] = 0.0;
lambda_net = eagle_metal_cooling_rate( lambda_net =
log10(u), cosmo->z, n_h, abundance_ratio, n_h_i, d_n_h, He_i, d_He, eagle_metal_cooling_rate(log10(u), cosmo->z, n_h, abundance_ratio, n_h_i,
cooling, /*dLambdaNet_du=*/NULL, element_lambda); d_n_h, He_i, d_He, cooling, element_lambda);
/* write cooling rate contributions to their own files. */ /* write cooling rate contributions to their own files. */
for (int j = 0; j < eagle_cooling_N_metal + 2; j++) { for (int j = 0; j < eagle_cooling_N_metal + 2; j++) {
...@@ -164,6 +164,7 @@ int main(int argc, char **argv) { ...@@ -164,6 +164,7 @@ int main(int argc, char **argv) {
struct part p; struct part p;
struct xpart xp; struct xpart xp;
struct phys_const internal_const; struct phys_const internal_const;
struct hydro_props hydro_properties;
struct cooling_function_data cooling; struct cooling_function_data cooling;
struct cosmology cosmo; struct cosmology cosmo;
struct space s; struct space s;
...@@ -210,8 +211,13 @@ int main(int argc, char **argv) { ...@@ -210,8 +211,13 @@ int main(int argc, char **argv) {
// Init units // Init units
units_init_from_params(&us, params, "InternalUnitSystem"); units_init_from_params(&us, params, "InternalUnitSystem");
// Init physical constants
phys_const_init(&us, params, &internal_const); phys_const_init(&us, params, &internal_const);
// Init porperties of hydro
hydro_props_init(&hydro_properties, &internal_const, &us, params);
// Init chemistry // Init chemistry
chemistry_init(params, &us, &internal_const, &chem_data); chemistry_init(params, &us, &internal_const, &chem_data);
chemistry_first_init_part(&internal_const, &us, &cosmo, &chem_data, &p, &xp); chemistry_first_init_part(&internal_const, &us, &cosmo, &chem_data, &p, &xp);
...@@ -228,19 +234,27 @@ int main(int argc, char **argv) { ...@@ -228,19 +234,27 @@ int main(int argc, char **argv) {
message("Redshift is %f", cosmo.z); message("Redshift is %f", cosmo.z);
// Init cooling // Init cooling
cooling_init(params, &us, &internal_const, &cooling); cooling_init(params, &us, &internal_const, &hydro_properties, &cooling);
cooling.H_reion_done = 1;
cooling_print(&cooling); cooling_print(&cooling);
cooling_update(&cosmo, &cooling, &s); cooling_update(&cosmo, &cooling, &s);
// Copy over the raw metals into the smoothed metals
memcpy(&p.chemistry_data.smoothed_metal_mass_fraction,
&p.chemistry_data.metal_mass_fraction,
chemistry_element_count * sizeof(float));
p.chemistry_data.smoothed_metal_mass_fraction_total =
p.chemistry_data.metal_mass_fraction_total;
// Calculate abundance ratios // Calculate abundance ratios
float abundance_ratio[(chemistry_element_count + 2)]; float abundance_ratio[(chemistry_element_count + 2)];
abundance_ratio_to_solar(&p, &cooling, abundance_ratio); abundance_ratio_to_solar(&p, &cooling, abundance_ratio);
// extract mass fractions, calculate table indices and offsets // extract mass fractions, calculate table indices and offsets
float XH = p.chemistry_data.metal_mass_fraction[chemistry_element_H]; float XH = p.chemistry_data.smoothed_metal_mass_fraction[chemistry_element_H];
float HeFrac = float XHe =
p.chemistry_data.metal_mass_fraction[chemistry_element_He] / p.chemistry_data.smoothed_metal_mass_fraction[chemistry_element_He];
(XH + p.chemistry_data.metal_mass_fraction[chemistry_element_He]); float HeFrac = XHe / (XH + XHe);
int He_i, n_h_i; int He_i, n_h_i;
float d_He, d_n_h; float d_He, d_n_h;
get_index_1d(cooling.HeFrac, eagle_cooling_N_He_frac, HeFrac, &He_i, &d_He); get_index_1d(cooling.HeFrac, eagle_cooling_N_He_frac, HeFrac, &He_i, &d_He);
...@@ -280,7 +294,7 @@ int main(int argc, char **argv) { ...@@ -280,7 +294,7 @@ int main(int argc, char **argv) {
// calculate cooling rates // calculate cooling rates
const double temperature = eagle_convert_u_to_temp( const double temperature = eagle_convert_u_to_temp(
log10(u), cosmo.z, 0, NULL, n_h_i, He_i, d_n_h, d_He, &cooling); log10(u), cosmo.z, n_h_i, He_i, d_n_h, d_He, &cooling);
const double cooling_du_dt = eagle_print_metal_cooling_rate(