Commit 3e07bc5c authored by Matthieu Schaller's avatar Matthieu Schaller
Browse files

Write the gas density surrounding the stars to the snapshots.

parent d398bbe8
......@@ -37,7 +37,7 @@ Statistics:
InitialConditions:
file_name: lowres8.hdf5 # The file to read
periodic: 0 # Are we running with periodic ICs?
stars_smoothing_length: 0.3
stars_smoothing_length: 0.5
# Parameters for the hydrodynamics scheme
SPH:
......
......@@ -223,7 +223,9 @@ void runner_do_stars_ghost(struct runner *r, struct cell *c, int timer) {
#ifdef SWIFT_DEBUG_CHECKS
if ((f > 0.f && h_new > h_old) || (f < 0.f && h_new < h_old))
error(
"Smoothing length correction not going in the right direction");
"Smoothing length correction not going in the right direction"
"sp->id=%lld sp->h=%e",
sp->id, sp->h);
#endif
/* Safety check: truncate to the range [ h_old/2 , 2h_old ]. */
......
......@@ -33,21 +33,6 @@ __attribute__((always_inline)) INLINE static float stars_compute_timestep(
return FLT_MAX;
}
/**
* @brief Initialises the s-particles for the first time
*
* This function is called only once just after the ICs have been
* read in to do some conversions.
*
* @param sp The particle to act upon
*/
__attribute__((always_inline)) INLINE static void stars_first_init_spart(
struct spart* sp) {
sp->time_bin = 0;
sp->birth_density = -1.f;
}
/**
* @brief Prepares a s-particle for its interactions
*
......@@ -64,6 +49,24 @@ __attribute__((always_inline)) INLINE static void stars_init_spart(
sp->density.wcount = 0.f;
sp->density.wcount_dh = 0.f;
sp->rho_gas = 0.f;
}
/**
* @brief Initialises the s-particles for the first time
*
* This function is called only once just after the ICs have been
* read in to do some conversions.
*
* @param sp The particle to act upon
*/
__attribute__((always_inline)) INLINE static void stars_first_init_spart(
struct spart* sp) {
sp->time_bin = 0;
sp->birth_density = -1.f;
stars_init_spart(sp);
}
/**
......@@ -132,6 +135,7 @@ __attribute__((always_inline)) INLINE static void stars_end_density(
const float h_inv_dim_plus_one = h_inv_dim * h_inv; /* 1/h^(d+1) */
/* Finish the calculation by inserting the missing h-factors */
sp->rho_gas *= h_inv_dim;
sp->density.wcount *= h_inv_dim;
sp->density.wcount_dh *= h_inv_dim_plus_one;
}
......@@ -146,14 +150,10 @@ __attribute__((always_inline)) INLINE static void stars_end_density(
__attribute__((always_inline)) INLINE static void stars_spart_has_no_neighbours(
struct spart* restrict sp, const struct cosmology* cosmo) {
/* Some smoothing length multiples. */
const float h = sp->h;
const float h_inv = 1.0f / h; /* 1/h */
const float h_inv_dim = pow_dimension(h_inv); /* 1/h^d */
/* Re-set problematic values */
sp->density.wcount = kernel_root * h_inv_dim;
sp->density.wcount = 0.f;
sp->density.wcount_dh = 0.f;
sp->rho_gas = 0.f;
}
/**
......
......@@ -16,6 +16,9 @@ runner_iact_nonsym_stars_density(float r2, const float *dx, float hi, float hj,
const struct part *restrict pj, float a,
float H) {
/* Get the gas mass. */
const float mj = pj->mass;
float wi, wi_dx;
/* Get r and 1/r. */
......@@ -31,6 +34,9 @@ runner_iact_nonsym_stars_density(float r2, const float *dx, float hi, float hj,
si->density.wcount += wi;
si->density.wcount_dh -= (hydro_dimension * wi + ui * wi_dx);
/* Compute contribution to the density */
si->rho_gas += mj * wi;
#ifdef DEBUG_INTERACTIONS_STARS
/* Update ngb counters */
if (si->num_ngb_density < MAX_NUM_OF_NEIGHBOURS_STARS)
......
......@@ -62,7 +62,7 @@ INLINE static void stars_write_particles(const struct spart *sparts,
int *num_fields) {
/* Say how much we want to write */
*num_fields = 8;
*num_fields = 9;
/* List what we want to write */
list[0] = io_make_output_field("Coordinates", DOUBLE, 3, UNIT_CONV_LENGTH,
......@@ -77,10 +77,12 @@ INLINE static void stars_write_particles(const struct spart *sparts,
sparts, h);
list[5] = io_make_output_field("BirthDensity", FLOAT, 1, UNIT_CONV_DENSITY,
sparts, birth_density);
list[6] = io_make_output_field("Initial_Masses", FLOAT, 1, UNIT_CONV_MASS,
list[6] = io_make_output_field("InitialMasses", FLOAT, 1, UNIT_CONV_MASS,
sparts, mass_init);
list[7] = io_make_output_field("Birth_time", FLOAT, 1, UNIT_CONV_TIME, sparts,
list[7] = io_make_output_field("Birthtime", FLOAT, 1, UNIT_CONV_TIME, sparts,
birth_time);
list[8] = io_make_output_field("GasDensity", FLOAT, 1, UNIT_CONV_DENSITY,
sparts, rho_gas);
}
/**
......
......@@ -58,9 +58,12 @@ struct spart {
/*! Initial star mass */
float mass_init;
/* Particle cutoff radius. */
/*! Particle smoothing length. */
float h;
/*! Density of the gas surrounding the star. */
float rho_gas;
/*! Particle time bin */
timebin_t time_bin;
......
Supports Markdown
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment