Added the quick Lyman-alpha star formation model

src/Makefile.am

@@ -210,6 +210,8 @@ nobase_noinst_HEADERS = align.h approx_math.h atomic.h barrier.h cycle.h error.h
                  potential/sine_wave/potential.h \
 		 star_formation/none/star_formation.h star_formation/none/star_formation_struct.h \
 		 star_formation/none/star_formation_io.h star_formation/none/star_formation_iact.h \
+		 star_formation/QLA/star_formation.h star_formation/QLA/star_formation_struct.h \
+		 star_formation/QLA/star_formation_io.h star_formation/QLA/star_formation_iact.h \
 		 star_formation/EAGLE/star_formation.h star_formation/EAGLE/star_formation_struct.h \
 		 star_formation/EAGLE/star_formation_io.h star_formation/EAGLE/star_formation_iact.h \
 		 star_formation/GEAR/star_formation.h star_formation/GEAR/star_formation_struct.h \

src/star_formation.h

@@ -33,6 +33,8 @@
 /* Import the right star formation law definition */
 #if defined(STAR_FORMATION_NONE)
 #include "./star_formation/none/star_formation.h"
+#elif defined(STAR_FORMATION_QLA)
+#include "./star_formation/QLA/star_formation.h"
 #elif defined(STAR_FORMATION_EAGLE)
 #include "./star_formation/EAGLE/star_formation.h"
 #elif defined(STAR_FORMATION_GEAR)

src/star_formation/QLA/star_formation.h

+/*******************************************************************************
+ * This file is part of SWIFT.
+ * Copyright (c) 2020 Matthieu Schaller (schaller@strw.leidenuniv.nl)
+ *
+ * This program is free software: you can redistribute it and/or modify
+ * it under the terms of the GNU Lesser General Public License as published
+ * by the Free Software Foundation, either version 3 of the License, or
+ * (at your option) any later version.
+ *
+ * This program is distributed in the hope that it will be useful,
+ * but WITHOUT ANY WARRANTY; without even the implied warranty of
+ * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
+ * GNU General Public License for more details.
+ *
+ * You should have received a copy of the GNU Lesser General Public License
+ * along with this program.  If not, see <http://www.gnu.org/licenses/>.
+ *
+ *******************************************************************************/
+#ifndef SWIFT_QLA_STAR_FORMATION_H
+#define SWIFT_QLA_STAR_FORMATION_H
+
+/* Local includes */
+#include "cosmology.h"
+#include "engine.h"
+#include "hydro.h"
+#include "parser.h"
+#include "part.h"
+#include "physical_constants.h"
+#include "stars.h"
+#include "units.h"
+
+/**
+ * @file src/star_formation/EAGLE/star_formation.h
+ * @brief Star formation model used in the EAGLE model
+ */
+
+/**
+ * @brief Properties of the EAGLE star formation model.
+ */
+struct star_formation {
+
+  /*! Critical overdensity */
+  double over_density;
+};
+
+/**
+ * @brief Calculate if the gas has the potential of becoming
+ * a star.
+ *
+ * @param starform the star formation law properties to use.
+ * @param p the gas particles.
+ * @param xp the additional properties of the gas particles.
+ * @param phys_const the physical constants in internal units.
+ * @param cosmo the cosmological parameters and properties.
+ * @param hydro_props The properties of the hydro scheme.
+ * @param us The internal system of units.
+ * @param cooling The cooling data struct.
+ * @param entropy_floor_props The entropy floor assumed in this run.
+ */
+INLINE static int star_formation_is_star_forming(
+    const struct part* p, const struct xpart* xp,
+    const struct star_formation* starform, const struct phys_const* phys_const,
+    const struct cosmology* cosmo, const struct hydro_props* hydro_props,
+    const struct unit_system* us, const struct cooling_function_data* cooling,
+    const struct entropy_floor_properties* entropy_floor_props) {
+
+  /* Minimal density (converted from mean baryonic density)
+     for star formation */
+  const double rho_mean_b_times_min_over_den =
+      cosmo->mean_density_Omega_b * starform->over_density;
+
+  /* Physical density of the particle */
+  const double physical_density = hydro_get_physical_density(p, cosmo);
+
+  /* Verify whether we are above the over-density threshold */
+  return (physical_density > rho_mean_b_times_min_over_den);
+}
+
+/**
+ * @brief Compute the star-formation rate of a given particle and store
+ * it into the #xpart.
+ *
+ * Nothing to do here. Particles that pass the SF criterion get automcatically
+ * converted to a star. No need to compute or store a star formation rate.
+ *
+ * @param p #part.
+ * @param xp the #xpart.
+ * @param starform the star formation law properties to use
+ * @param phys_const the physical constants in internal units.
+ * @param hydro_props The properties of the hydro scheme.
+ * @param cosmo the cosmological parameters and properties.
+ * @param dt_star The time-step of this particle.
+ */
+INLINE static void star_formation_compute_SFR(
+    const struct part* restrict p, struct xpart* restrict xp,
+    const struct star_formation* starform, const struct phys_const* phys_const,
+    const struct hydro_props* hydro_props, const struct cosmology* cosmo,
+    const double dt_star) {
+
+  xp->sf_data.convert_to_star = 1;
+}
+
+/**
+ * @brief Decides whether a particle should be converted into a
+ * star or not.
+ *
+ * Equation 21 of Schaye & Dalla Vecchia 2008.
+ *
+ * @param p The #part.
+ * @param xp The #xpart.
+ * @param starform The properties of the star formation model.
+ * @param e The #engine (for random numbers).
+ * @param dt_star The time-step of this particle
+ * @return 1 if a conversion should be done, 0 otherwise.
+ */
+INLINE static int star_formation_should_convert_to_star(
+    const struct part* p, const struct xpart* xp,
+    const struct star_formation* starform, const struct engine* e,
+    const double dt_star) {
+
+  return xp->sf_data.convert_to_star;
+}
+
+/**
+ * @brief Update the SF properties of a particle that is not star forming.
+ *
+ * Nothing to do here in the quick Lyman-alpha model.
+ *
+ * @param p The #part.
+ * @param xp The #xpart.
+ * @param e The #engine.
+ * @param starform The properties of the star formation model.
+ * @param with_cosmology Are we running with cosmology switched on?
+ */
+INLINE static void star_formation_update_part_not_SFR(
+    struct part* p, struct xpart* xp, const struct engine* e,
+    const struct star_formation* starform, const int with_cosmology) {}
+
+/**
+ * @brief Copies the properties of the gas particle over to the
+ * star particle
+ *
+ * @param e The #engine
+ * @param p the gas particles.
+ * @param xp the additional properties of the gas particles.
+ * @param sp the new created star particle with its properties.
+ * @param starform the star formation law properties to use.
+ * @param cosmo the cosmological parameters and properties.
+ * @param with_cosmology if we run with cosmology.
+ * @param phys_const the physical constants in internal units.
+ * @param hydro_props The properties of the hydro scheme.
+ * @param us The internal system of units.
+ * @param cooling The cooling data struct.
+ */
+INLINE static void star_formation_copy_properties(
+    const struct part* p, const struct xpart* xp, struct spart* sp,
+    const struct engine* e, const struct star_formation* starform,
+    const struct cosmology* cosmo, const int with_cosmology,
+    const struct phys_const* phys_const,
+    const struct hydro_props* restrict hydro_props,
+    const struct unit_system* restrict us,
+    const struct cooling_function_data* restrict cooling) {
+
+  /* Store the current mass */
+  sp->mass = hydro_get_mass(p);
+
+  /* Store either the birth_scale_factor or birth_time depending  */
+  if (with_cosmology) {
+    sp->birth_scale_factor = cosmo->a;
+  } else {
+    sp->birth_time = e->time;
+  }
+}
+
+/**
+ * @brief initialization of the star formation law
+ *
+ * @param parameter_file The parsed parameter file
+ * @param phys_const Physical constants in internal units
+ * @param us The current internal system of units.
+ * @param hydro_props The propertis of the hydro model.
+ * @param starform the star formation law properties to initialize
+ */
+INLINE static void starformation_init_backend(
+    struct swift_params* parameter_file, const struct phys_const* phys_const,
+    const struct unit_system* us, const struct hydro_props* hydro_props,
+    struct star_formation* starform) {
+
+  /* Read the critical density contrast from the parameter file*/
+  starform->over_density =
+      parser_get_param_double(parameter_file, "QLAStarFormation:over_density");
+}
+
+/**
+ * @brief Prints the used parameters of the star formation law
+ *
+ * @param starform the star formation law properties.
+ * */
+INLINE static void starformation_print_backend(
+    const struct star_formation* starform) {
+
+  message("Star formation law is Quick Lyman-alpha");
+  message("Over-density for star formation: %f", starform->over_density);
+}
+
+/**
+ * @brief Finishes the density calculation.
+ *
+ * Nothing to do here. We do not need to compute any quantity in the hydro
+ * density loop for the quick Lyman-alpha star formation model.
+ *
+ * @param p The particle to act upon
+ * @param xp The extra particle to act upon
+ * @param cd The global star_formation information.
+ * @param cosmo The current cosmological model.
+ */
+__attribute__((always_inline)) INLINE static void star_formation_end_density(
+    struct part* restrict p, struct xpart* restrict xp,
+    const struct star_formation* cd, const struct cosmology* cosmo) {}
+
+/**
+ * @brief Sets all particle fields to sensible values when the #part has 0 ngbs.
+ *
+ * Nothing to do here. We do not need to compute any quantity in the hydro
+ * density loop for the quick Lyman-alpha star formation model.
+ *
+ * @param p The particle to act upon
+ * @param xp The extended particle data to act upon
+ * @param cd #star_formation containing star_formation informations.
+ * @param cosmo The current cosmological model.
+ */
+__attribute__((always_inline)) INLINE static void
+star_formation_part_has_no_neighbours(struct part* restrict p,
+                                      struct xpart* restrict xp,
+                                      const struct star_formation* cd,
+                                      const struct cosmology* cosmo) {}
+
+/**
+ * @brief Sets the star_formation properties of the (x-)particles to a valid
+ * state to start the density loop.
+ *
+ * Nothing to do here. We do not need to compute any quantity in the hydro
+ * density loop for the quick Lyman-alpha star formation model.
+ *
+ * @param data The global star_formation information used for this run.
+ * @param p Pointer to the particle data.
+ */
+__attribute__((always_inline)) INLINE static void star_formation_init_part(
+    struct part* restrict p, const struct star_formation* data) {}
+
+/**
+ * @brief Sets the star_formation properties of the (x-)particles to a valid
+ * start state at the beginning of the simulation after the ICs have been read.
+ *
+ * Mark the particles as not needing to be converted to stars.
+ *
+ * @param phys_const The physical constant in internal units.
+ * @param us The unit system.
+ * @param cosmo The current cosmological model.
+ * @param data The global star_formation information used for this run.
+ * @param p Pointer to the particle data.
+ * @param xp Pointer to the extended particle data.
+ */
+__attribute__((always_inline)) INLINE static void
+star_formation_first_init_part(const struct phys_const* restrict phys_const,
+                               const struct unit_system* restrict us,
+                               const struct cosmology* restrict cosmo,
+                               const struct star_formation* data,
+                               const struct part* restrict p,
+                               struct xpart* restrict xp) {
+
+  xp->sf_data.convert_to_star = 0;
+}
+
+/**
+ * @brief Split the star formation content of a particle into n pieces
+ *
+ * Nothing to do here.
+ *
+ * @param p The #part.
+ * @param xp The #xpart.
+ * @param n The number of pieces to split into.
+ */
+__attribute__((always_inline)) INLINE static void star_formation_split_part(
+    struct part* p, struct xpart* xp, const double n) {}
+
+#endif /* SWIFT_QLA_STAR_FORMATION_H */

src/star_formation/QLA/star_formation_iact.h

+/*******************************************************************************
+ * This file is part of SWIFT.
+ * Copyright (c) 2020 Matthieu Schaller (schaller@strw.leidenuniv.nl)
+ *
+ * This program is free software: you can redistribute it and/or modify
+ * it under the terms of the GNU Lesser General Public License as published
+ * by the Free Software Foundation, either version 3 of the License, or
+ * (at your option) any later version.
+ *
+ * This program is distributed in the hope that it will be useful,
+ * but WITHOUT ANY WARRANTY; without even the implied warranty of
+ * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
+ * GNU General Public License for more details.
+ *
+ * You should have received a copy of the GNU Lesser General Public License
+ * along with this program.  If not, see <http://www.gnu.org/licenses/>.
+ *
+ ******************************************************************************/
+#ifndef SWIFT_QLA_STAR_FORMATION_IACT_H
+#define SWIFT_QLA_STAR_FORMATION_IACT_H
+
+/**
+ * @file EAGLE/star_formation_iact.h
+ * @brief Density computation
+ */
+
+/**
+ * @brief do star_formation computation after the runner_iact_density (symmetric
+ * version)
+ *
+ * @param r2 Comoving square distance between the two particles.
+ * @param dx Comoving vector separating both particles (pi - pj).
+ * @param hi Comoving smoothing-length of particle i.
+ * @param hj Comoving smoothing-length of particle j.
+ * @param pi First particle.
+ * @param pj Second particle.
+ * @param a Current scale factor.
+ * @param H Current Hubble parameter.
+ */
+__attribute__((always_inline)) INLINE static void runner_iact_star_formation(
+    float r2, const float *dx, float hi, float hj, struct part *restrict pi,
+    struct part *restrict pj, float a, float H) {
+
+  /* Nothing to do here. We do not need to compute any quantity in the hydro
+     density loop for the QLA star formation model. */
+}
+
+/**
+ * @brief do star_formation computation after the runner_iact_density (non
+ * symmetric version)
+ *
+ * @param r2 Comoving square distance between the two particles.
+ * @param dx Comoving vector separating both particles (pi - pj).
+ * @param hi Comoving smoothing-length of particle i.
+ * @param hj Comoving smoothing-length of particle j.
+ * @param pi First particle.
+ * @param pj Second particle (not updated).
+ * @param a Current scale factor.
+ * @param H Current Hubble parameter.
+ */
+__attribute__((always_inline)) INLINE static void
+runner_iact_nonsym_star_formation(float r2, const float *dx, float hi, float hj,
+                                  struct part *restrict pi,
+                                  const struct part *restrict pj, float a,
+                                  float H) {
+
+  /* Nothing to do here. We do not need to compute any quantity in the hydro
+     density loop for the QLA star formation model. */
+}
+
+#endif /* SWIFT_QLA_STAR_FORMATION_IACT_H */

src/star_formation/QLA/star_formation_io.h

+/*******************************************************************************
+ * This file is part of SWIFT.
+ * Copyright (c) 2020 Matthieu Schaller (schaller@strw.leidenuniv.nl)
+ *
+ * This program is free software: you can redistribute it and/or modify
+ * it under the terms of the GNU Lesser General Public License as published
+ * by the Free Software Foundation, either version 3 of the License, or
+ * (at your option) any later version.
+ *
+ * This program is distributed in the hope that it will be useful,
+ * but WITHOUT ANY WARRANTY; without even the implied warranty of
+ * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
+ * GNU General Public License for more details.
+ *
+ * You should have received a copy of the GNU Lesser General Public License
+ * along with this program.  If not, see <http://www.gnu.org/licenses/>.
+ *
+ ******************************************************************************/
+#ifndef SWIFT_STAR_FORMATION_QLA_IO_H
+#define SWIFT_STAR_FORMATION_QLA_IO_H
+
+/* Config parameters. */
+#include "../config.h"
+
+/* Local includes */
+#include "io_properties.h"
+
+/**
+ * @brief Specifies which particle fields to write to a dataset
+ *
+ * @param parts The particle array.
+ * @param xparts The extended data particle array.
+ * @param list The list of i/o properties to write.
+ *
+ * @return Returns the number of fields to write.
+ */
+__attribute__((always_inline)) INLINE static int star_formation_write_particles(
+    const struct part* parts, const struct xpart* xparts,
+    struct io_props* list) {
+
+  /* No SF-related fields in the Quick Lyman-alpha model */
+
+  return 0;
+}
+
+#endif /* SWIFT_STAR_FORMATION_QLA_IO_H */

src/star_formation/QLA/star_formation_struct.h

+/*******************************************************************************
+ * This file is part of SWIFT.
+ * Copyright (c) 2020 Matthieu Schaller (schaller@strw.leidenuniv.nl)
+ *
+ * This program is free software: you can redistribute it and/or modify
+ * it under the terms of the GNU Lesser General Public License as published
+ * by the Free Software Foundation, either version 3 of the License, or
+ * (at your option) any later version.
+ *
+ * This program is distributed in the hope that it will be useful,
+ * but WITHOUT ANY WARRANTY; without even the implied warranty of
+ * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
+ * GNU General Public License for more details.
+ *
+ * You should have received a copy of the GNU Lesser General Public License
+ * along with this program.  If not, see <http://www.gnu.org/licenses/>.
+ *
+ ******************************************************************************/
+#ifndef SWIFT_QLA_STAR_FORMATION_STRUCT_H
+#define SWIFT_QLA_STAR_FORMATION_STRUCT_H
+
+/**
+ * @brief Star-formation-related properties stored in the extended particle
+ * data.
+ */
+struct star_formation_xpart_data {
+
+  char convert_to_star;
+};
+
+#endif /* SWIFT_QLA_STAR_FORMATION_STRUCT_H */

src/star_formation_iact.h

@@ -31,6 +31,8 @@
 /* Import the right star formation law definition */
 #if defined(STAR_FORMATION_NONE)
 #include "./star_formation/none/star_formation_iact.h"
+#elif defined(STAR_FORMATION_QLA)
+#include "./star_formation/QLA/star_formation_iact.h"
 #elif defined(STAR_FORMATION_EAGLE)
 #include "./star_formation/EAGLE/star_formation_iact.h"
 #elif defined(STAR_FORMATION_GEAR)

src/star_formation_io.h

@@ -30,6 +30,8 @@
 /* Import the right cooling definition */
 #if defined(STAR_FORMATION_NONE)
 #include "./star_formation/none/star_formation_io.h"
+#elif defined(STAR_FORMATION_QLA)
+#include "./star_formation/QLA/star_formation_io.h"
 #elif defined(STAR_FORMATION_EAGLE)
 #include "./star_formation/EAGLE/star_formation_io.h"
 #elif defined(STAR_FORMATION_GEAR)

src/star_formation_logger.h

@@ -30,6 +30,8 @@
 /* Import the right SFH logger definition */
 #if defined(STAR_FORMATION_NONE)
 #include "./star_formation/none/star_formation_logger.h"
+#elif defined(STAR_FORMATION_QLA)
+#include "./star_formation/none/star_formation_logger.h"
 #elif defined(STAR_FORMATION_EAGLE)
 #include "./star_formation/EAGLE/star_formation_logger.h"
 #elif defined(STAR_FORMATION_GEAR)

src/star_formation_logger_struct.h

@@ -30,6 +30,8 @@
 /* Import the right SFH logger struct definition */
 #if defined(STAR_FORMATION_NONE)
 #include "./star_formation/none/star_formation_logger_struct.h"
+#elif defined(STAR_FORMATION_QLA)
+#include "./star_formation/none/star_formation_logger_struct.h"
 #elif defined(STAR_FORMATION_EAGLE)
 #include "./star_formation/EAGLE/star_formation_logger_struct.h"
 #elif defined(STAR_FORMATION_GEAR)

src/star_formation_struct.h

@@ -30,6 +30,8 @@
 /* Import the right star formation definition */
 #if defined(STAR_FORMATION_NONE)
 #include "./star_formation/none/star_formation_struct.h"
+#elif defined(STAR_FORMATION_QLA)
+#include "./star_formation/QLA/star_formation_struct.h"
 #elif defined(STAR_FORMATION_EAGLE)
 #include "./star_formation/EAGLE/star_formation_struct.h"
 #elif defined(STAR_FORMATION_GEAR)