Fix API

configure.ac

@@ -1222,7 +1222,7 @@ if test "x$with_grackle" != "xno"; then
       [AC_DEFINE([HAVE_GRACKLE],1,[The GRACKLE library appears to be present.])
         AC_DEFINE([CONFIG_BFLOAT_8],1,[Use doubles in grackle])
       ],
-      [AC_MSG_ERROR(Cannot find grackle library!)],
+      [AC_MSG_ERROR(Cannot find grackle library! You need to have a version >= 3.1.1)],
       [$GRACKLE_LIBS])
 fi
 AC_SUBST([GRACKLE_LIBS])

doc/RTD/source/SubgridModels/GEAR/index.rst

@@ -57,7 +57,7 @@ The code uses an iterative method in order to find the correct initial compositi
 In order to compile SWIFT with Grackle, you need to provide the options ``with-chemistry=GEAR`` and ``with-grackle=$GRACKLE_ROOT``
 where ``$GRACKLE_ROOT`` is the root of the install directory (not the ``lib``).
 
-You will need a Grackle version later than 3.1. To compile it, run
+You will need a Grackle version later than 3.1.1. To compile it, run
 the following commands from the root directory of Grackle:
 ``./configure; cd src/clib``.
 Update the variables ``LOCAL_HDF5_INSTALL`` and ``MACH_INSTALL_PREFIX`` in
@@ -103,6 +103,9 @@ A particle is in the star forming regime if:
 
 A star will be able to form if a randomly drawn number is below :math:`\frac{m_g}{m_\star}\left(1 - \exp\left(-c_\star \Delta t / t_\textrm{ff}\right)\right)` where :math:`t_\textrm{ff}` is the free fall time, :math:`\Delta t` is the time step of the particle and :math:`c_\star` is the star formation coefficient (``GEARStarFormation:star_formation_efficiency``), :math:`m_g` the mass of the gas particle and :math:`m_\star` the mass of the possible future star. The mass of the star is computed from the average gas mass in the initial conditions divided by the number of possible stars formed per gas particle (``GEARStarFormation:n_stars_per_particle``). When we cannot have enough mass to form a second star (defined with the fraction of mass ``GEARStarFormation:min_mass_frac``), we fully convert the gas particle into a stellar particle. Once the star is formed, we move it a bit in a random direction and fraction of the smoothing length in order to avoid any division by 0.
 
+Currently, only the following hydro schemes are compatible: SPHENIX and Gadget2.
+Implementing the other hydro schemes is not complicated but requires some careful thinking about the cosmological terms in the definition of the velocity divergence (comoving vs non comoving coordinates and if the Hubble flow is included or not).
+
 .. code:: YAML
 
   GEARStarFormation:

examples/Cooling/CoolingBox/coolingBox.yml

@@ -39,20 +39,22 @@ LambdaCooling:
   lambda_nH2_cgs:              1e-22 # Cooling rate divided by square Hydrogen number density (in cgs units [erg * s^-1 * cm^3])
   cooling_tstep_mult:          1.0        # Dimensionless pre-factor for the time-step condition
 
-# Cooling with Grackle 2.0
+# Cooling with Grackle 3.0
 GrackleCooling:
-  CloudyTable: CloudyData_UVB=HM2012.h5 # Name of the Cloudy Table (available on the grackle bitbucket repository)
-  WithUVbackground: 0 # Enable or not the UV background
-  Redshift: 0 # Redshift to use (-1 means time based redshift)
-  WithMetalCooling: 1 # Enable or not the metal cooling
-  ProvideVolumetricHeatingRates: 0 # User provide volumetric heating rates
-  ProvideSpecificHeatingRates: 0 # User provide specific heating rates
-  SelfShieldingMethod: 0 # Grackle (<= 3) or Gear self shielding method
-  MaxSteps: 1000
-  ConvergenceLimit: 1e-2
-  
-GearChemistry:
-  InitialMetallicity: 0.01295
+  cloudy_table: CloudyData_UVB=HM2012.h5       # Name of the Cloudy Table (available on the grackle bitbucket repository)
+  with_UV_background: 0                        # Enable or not the UV background
+  redshift: 0                                  # Redshift to use (-1 means time based redshift)
+  with_metal_cooling: 1                        # Enable or not the metal cooling
+  provide_volumetric_heating_rates: 0          # (optional) User provide volumetric heating rates
+  provide_specific_heating_rates: 0            # (optional) User provide specific heating rates
+  max_steps: 10000                             # (optional) Max number of step when computing the initial composition
+  convergence_limit: 1e-2                      # (optional) Convergence threshold (relative) for initial composition
+  thermal_time_myr: 5                          # (optional) Time (in Myr) for adiabatic cooling after a feedback event.
+  self_shielding_method: 0                    # (optional) Grackle (1->3 for Grackle's ones, 0 for none and -1 for GEAR)
+  self_shielding_threshold_atom_per_cm3: 0.007 # Required only with GEAR's self shielding. Density threshold of the self shielding
+
+GEARChemistry:
+  initial_metallicity: 0.01295
 
 EAGLECooling:
   dir_name:                ./coolingtables/
@@ -74,3 +76,5 @@ EAGLEChemistry:             # Solar abundances
   init_abundance_Silicon:   6.825874e-4
   init_abundance_Iron:      1.1032152e-3
 
+GEARPressureFloor:
+  jeans_factor: 0.       # Number of particles required to suppose a resolved clump and avoid the pressure floor.

src/chemistry/none/chemistry.h

@@ -49,6 +49,17 @@ chemistry_get_element_name(enum chemistry_element elem) {
   return chemistry_element_names[elem];
 }
 
+/**
+ * @brief Copies the chemistry properties of the gas particle over to the
+ * star particle.
+ *
+ * @param p the gas particles.
+ * @param xp the additional properties of the gas particles.
+ * @param sp the new created star particle with its properties.
+ */
+INLINE static void chemistry_copy_star_formation_properties(
+    const struct part* p, const struct xpart* xp, struct spart* sp) {}
+
 /**
  * @brief Initialises the chemistry properties.
  *

src/chemistry/none/chemistry_struct.h

@@ -43,6 +43,13 @@ struct chemistry_global_data {};
  */
 struct chemistry_part_data {};
 
+/**
+ * @brief Chemistry properties carried by the #spart.
+ *
+ * Nothing here.
+ */
+struct chemistry_spart_data {};
+
 /**
  * @brief Chemistry properties carried by the #bpart.
  *

src/cooling/grackle/cooling.c

@@ -18,7 +18,7 @@
  ******************************************************************************/
 /**
  * @file src/cooling/grackle/cooling.c
- * @brief Cooling using the GRACKLE 3.0 library.
+ * @brief Cooling using the GRACKLE 3.1.1 library.
  */
 
 #include "../config.h"

src/cooling/grackle/cooling.h

@@ -21,7 +21,7 @@
 
 /**
  * @file src/cooling/grackle/cooling.h
- * @brief Cooling using the GRACKLE 3.0 library.
+ * @brief Cooling using the GRACKLE 3.1.1 library.
  */
 
 /* Some standard headers. */

src/parallel_io.c

@@ -1010,6 +1010,8 @@ void read_ic_parallel(char* fileName, const struct unit_system* internal_units,
         if (with_stars) {
           Nparticles = *Nstars;
           stars_read_particles(*sparts, list, &num_fields);
+          num_fields +=
+              star_formation_read_particles(*sparts, list + num_fields);
         }
         break;
 

src/serial_io.c

@@ -804,6 +804,8 @@ void read_ic_serial(char* fileName, const struct unit_system* internal_units,
             if (with_stars) {
               Nparticles = *Nstars;
               stars_read_particles(*sparts, list, &num_fields);
+              num_fields +=
+                  star_formation_read_particles(*sparts, list + num_fields);
             }
             break;
 

src/single_io.c

@@ -675,6 +675,8 @@ void read_ic_single(const char* fileName,
         if (with_stars) {
           Nparticles = *Nstars;
           stars_read_particles(*sparts, list, &num_fields);
+          num_fields +=
+              star_formation_read_particles(*sparts, list + num_fields);
         }
         break;
 

src/star_formation/EAGLE/star_formation_io.h

@@ -25,6 +25,19 @@
 /* Local includes */
 #include "io_properties.h"
 
+/**
+ * @brief Specifies which s-particle fields to read from a dataset
+ *
+ * @param sparts The s-particle array.
+ * @param list The list of i/o properties to read.
+ *
+ * @return num_fields The number of i/o fields to read.
+ */
+INLINE static int star_formation_read_particles(struct spart* sparts,
+                                                struct io_props* list) {
+  return 0;
+}
+
 /**
  * @brief Specifies which particle fields to write to a dataset
  *

src/star_formation/GEAR/star_formation.h

@@ -21,6 +21,7 @@
 #define SWIFT_GEAR_STAR_FORMATION_H
 
 /* Local includes */
+#include "chemistry.h"
 #include "cooling.h"
 #include "cosmology.h"
 #include "engine.h"

src/star_formation/GEAR/star_formation_io.h

@@ -26,6 +26,24 @@
 /* Local includes */
 #include "io_properties.h"
 
+/**
+ * @brief Specifies which s-particle fields to read from a dataset
+ *
+ * @param sparts The s-particle array.
+ * @param list The list of i/o properties to read.
+ *
+ * @return num_fields The number of i/o fields to read.
+ */
+INLINE static int star_formation_read_particles(struct spart* sparts,
+                                                struct io_props* list) {
+
+  /* List what we want to read */
+  list[0] = io_make_input_field("BirthMass", FLOAT, 1, OPTIONAL, UNIT_CONV_MASS,
+                                sparts, sf_data.birth_mass);
+
+  return 1;
+}
+
 /**
  * @brief Specifies which particle fields to write to a dataset
  *
@@ -113,7 +131,9 @@ INLINE static void starformation_init_backend(
   starform->min_mass_frac_plus_one += 1.;
 
   /* Get the jeans factor */
-  starform->n_jeans_2_3 = pow(pressure_floor_props.n_jeans, 2. / 3.);
+  starform->n_jeans_2_3 =
+      parser_get_param_float(parameter_file, "GEARPressureFloor:jeans_factor");
+  starform->n_jeans_2_3 = pow(starform->n_jeans_2_3, 2. / 3.);
 
   /* Apply unit change */
   starform->maximal_temperature *=

src/star_formation/QLA/star_formation_io.h

@@ -25,6 +25,19 @@
 /* Local includes */
 #include "io_properties.h"
 
+/**
+ * @brief Specifies which s-particle fields to read from a dataset
+ *
+ * @param sparts The s-particle array.
+ * @param list The list of i/o properties to read.
+ *
+ * @return num_fields The number of i/o fields to read.
+ */
+INLINE static int star_formation_read_particles(struct spart* sparts,
+                                                struct io_props* list) {
+  return 0;
+}
+
 /**
  * @brief Specifies which particle fields to write to a dataset
  *

src/star_formation/none/star_formation_io.h

@@ -25,6 +25,19 @@
 /* Local includes */
 #include "io_properties.h"
 
+/**
+ * @brief Specifies which s-particle fields to read from a dataset
+ *
+ * @param sparts The s-particle array.
+ * @param list The list of i/o properties to read.
+ *
+ * @return num_fields The number of i/o fields to read.
+ */
+INLINE static int star_formation_read_particles(struct spart* sparts,
+                                                struct io_props* list) {
+  return 0;
+}
+
 /**
  * @brief Specifies which particle fields to write to a dataset
  *

src/stars/GEAR/stars.h

@@ -68,8 +68,6 @@ __attribute__((always_inline)) INLINE static void stars_first_init_spart(
     const int with_cosmology, const double scale_factor, const double time) {
 
   sp->time_bin = 0;
-  sp->sf_data.birth_density = 0.f;
-  // sp->birth_time = 0.;
 
   stars_init_spart(sp);
 }

src/stars/GEAR/stars_io.h

@@ -34,7 +34,7 @@ INLINE static void stars_read_particles(struct spart *sparts,
                                         int *num_fields) {
 
   /* Say how much we want to read */
-  *num_fields = 7;
+  *num_fields = 6;
 
   /* List what we want to read */
   list[0] = io_make_input_field("Coordinates", DOUBLE, 3, COMPULSORY,
@@ -47,10 +47,7 @@ INLINE static void stars_read_particles(struct spart *sparts,
                                 UNIT_CONV_NO_UNITS, sparts, id);
   list[4] = io_make_input_field("SmoothingLength", FLOAT, 1, OPTIONAL,
                                 UNIT_CONV_LENGTH, sparts, h);
-  list[5] = io_make_input_field("BirthMass", FLOAT, 1, OPTIONAL, UNIT_CONV_MASS,
-                                sparts, sf_data.birth_mass);
-
-  list[6] = io_make_input_field("BirthTime", FLOAT, 1, OPTIONAL, UNIT_CONV_MASS,
+  list[5] = io_make_input_field("BirthTime", FLOAT, 1, OPTIONAL, UNIT_CONV_MASS,
                                 sparts, birth_time);
 }
 

tools/task_plots/plot_tasks.py

@@ -584,7 +584,7 @@ for rank in ranks:
     ax.yaxis.set_major_locator(loc)
     ax.grid(True, which="major", axis="y", linestyle="-")
 
-    pl.show()
+    # pl.show()
     if mpimode:
         outpng = outbase + str(rank) + ".png"
     else: