Commit 2e33483d authored by Matthieu Schaller's avatar Matthieu Schaller
Browse files

Also update the other hydro schemes with the new i/o API and the particle wrapping.

parent 1d64aa96
......@@ -53,6 +53,20 @@ void hydro_read_particles(struct part* parts, struct io_props* list,
UNIT_CONV_DENSITY, parts, rho);
}
void convert_part_pos(const struct engine* e, const struct part* p,
double* ret) {
if (e->s->periodic) {
ret[0] = box_wrap(p->x[0], 0.0, e->s->dim[0]);
ret[1] = box_wrap(p->x[1], 0.0, e->s->dim[1]);
ret[2] = box_wrap(p->x[2], 0.0, e->s->dim[2]);
} else {
ret[0] = p->x[0];
ret[1] = p->x[1];
ret[2] = p->x[2];
}
}
/**
* @brief Specifies which particle fields to write to a dataset
*
......@@ -66,8 +80,8 @@ void hydro_write_particles(struct part* parts, struct io_props* list,
*num_fields = 8;
/* List what we want to write */
list[0] = io_make_output_field("Coordinates", DOUBLE, 3, UNIT_CONV_LENGTH,
parts, x);
list[0] = io_make_output_field_convert_part(
"Coordinates", DOUBLE, 3, UNIT_CONV_LENGTH, parts, convert_part_pos);
list[1] =
io_make_output_field("Velocities", FLOAT, 3, UNIT_CONV_SPEED, parts, v);
list[2] =
......
......@@ -125,6 +125,16 @@ struct part {
/* Particle time-bin */
timebin_t time_bin;
#ifdef SWIFT_DEBUG_CHECKS
/* Time of the last drift */
integertime_t ti_drift;
/* Time of the last kick */
integertime_t ti_kick;
#endif
} SWIFT_STRUCT_ALIGN;
#endif /* SWIFT_DEFAULT_HYDRO_PART_H */
......@@ -69,10 +69,11 @@ void hydro_read_particles(struct part* parts, struct io_props* list,
*
* @param e #engine.
* @param p Particle.
* @return Internal energy of the particle
* @param ret (return) Internal energy of the particle
*/
float convert_u(struct engine* e, struct part* p) {
return hydro_get_internal_energy(p);
void convert_u(const struct engine* e, const struct part* p, float* ret) {
ret[0] = hydro_get_internal_energy(p);
}
/**
......@@ -80,10 +81,10 @@ float convert_u(struct engine* e, struct part* p) {
*
* @param e #engine.
* @param p Particle.
* @return Entropic function of the particle
* @param ret (return) Entropic function of the particle
*/
float convert_A(struct engine* e, struct part* p) {
return hydro_get_entropy(p);
void convert_A(const struct engine* e, const struct part* p, float* ret) {
ret[0] = hydro_get_entropy(p);
}
/**
......@@ -93,9 +94,9 @@ float convert_A(struct engine* e, struct part* p) {
* @param p Particle.
* @return Total energy of the particle
*/
float convert_Etot(struct engine* e, struct part* p) {
void convert_Etot(const struct engine* e, const struct part* p, float* ret) {
#ifdef GIZMO_TOTAL_ENERGY
return p->conserved.energy;
ret[0] = p->conserved.energy;
#else
float momentum2;
......@@ -103,10 +104,24 @@ float convert_Etot(struct engine* e, struct part* p) {
p->conserved.momentum[1] * p->conserved.momentum[1] +
p->conserved.momentum[2] * p->conserved.momentum[2];
return p->conserved.energy + 0.5f * momentum2 / p->conserved.mass;
ret[0] = p->conserved.energy + 0.5f * momentum2 / p->conserved.mass;
#endif
}
void convert_part_pos(const struct engine* e, const struct part* p,
double* ret) {
if (e->s->periodic) {
ret[0] = box_wrap(p->x[0], 0.0, e->s->dim[0]);
ret[1] = box_wrap(p->x[1], 0.0, e->s->dim[1]);
ret[2] = box_wrap(p->x[2], 0.0, e->s->dim[2]);
} else {
ret[0] = p->x[0];
ret[1] = p->x[1];
ret[2] = p->x[2];
}
}
/**
* @brief Specifies which particle fields to write to a dataset
*
......@@ -120,8 +135,8 @@ void hydro_write_particles(struct part* parts, struct io_props* list,
*num_fields = 10;
/* List what we want to write */
list[0] = io_make_output_field("Coordinates", DOUBLE, 3, UNIT_CONV_LENGTH,
parts, x);
list[0] = io_make_output_field_convert_part(
"Coordinates", DOUBLE, 3, UNIT_CONV_LENGTH, parts, convert_part_pos);
list[1] = io_make_output_field("Velocities", FLOAT, 3, UNIT_CONV_SPEED, parts,
primitives.v);
list[2] = io_make_output_field("Masses", FLOAT, 1, UNIT_CONV_MASS, parts,
......@@ -130,18 +145,17 @@ void hydro_write_particles(struct part* parts, struct io_props* list,
parts, h);
list[4] = io_make_output_field_convert_part("InternalEnergy", FLOAT, 1,
UNIT_CONV_ENERGY_PER_UNIT_MASS,
parts, primitives.P, convert_u);
parts, convert_u);
list[5] = io_make_output_field("ParticleIDs", ULONGLONG, 1,
UNIT_CONV_NO_UNITS, parts, id);
list[6] = io_make_output_field("Density", FLOAT, 1, UNIT_CONV_DENSITY, parts,
primitives.rho);
list[7] = io_make_output_field_convert_part(
"Entropy", FLOAT, 1, UNIT_CONV_ENTROPY, parts, primitives.P, convert_A);
"Entropy", FLOAT, 1, UNIT_CONV_ENTROPY, parts, convert_A);
list[8] = io_make_output_field("Pressure", FLOAT, 1, UNIT_CONV_PRESSURE,
parts, primitives.P);
list[9] =
io_make_output_field_convert_part("TotEnergy", FLOAT, 1, UNIT_CONV_ENERGY,
parts, conserved.energy, convert_Etot);
list[9] = io_make_output_field_convert_part(
"TotEnergy", FLOAT, 1, UNIT_CONV_ENERGY, parts, convert_Etot);
}
/**
......
......@@ -69,14 +69,28 @@ void hydro_read_particles(struct part* parts, struct io_props* list,
UNIT_CONV_DENSITY, parts, rho);
}
float convert_S(struct engine* e, struct part* p) {
void convert_S(const struct engine* e, const struct part* p, float* ret) {
return hydro_get_entropy(p);
ret[0] = hydro_get_entropy(p);
}
float convert_P(struct engine* e, struct part* p) {
void convert_P(const struct engine* e, const struct part* p, float* ret) {
return hydro_get_pressure(p);
ret[0] = hydro_get_pressure(p);
}
void convert_part_pos(const struct engine* e, const struct part* p,
double* ret) {
if (e->s->periodic) {
ret[0] = box_wrap(p->x[0], 0.0, e->s->dim[0]);
ret[1] = box_wrap(p->x[1], 0.0, e->s->dim[1]);
ret[2] = box_wrap(p->x[2], 0.0, e->s->dim[2]);
} else {
ret[0] = p->x[0];
ret[1] = p->x[1];
ret[2] = p->x[2];
}
}
/**
......@@ -92,8 +106,8 @@ void hydro_write_particles(struct part* parts, struct io_props* list,
*num_fields = 10;
/* List what we want to write */
list[0] = io_make_output_field("Coordinates", DOUBLE, 3, UNIT_CONV_LENGTH,
parts, x);
list[0] = io_make_output_field_convert_part(
"Coordinates", DOUBLE, 3, UNIT_CONV_LENGTH, parts, convert_part_pos);
list[1] =
io_make_output_field("Velocities", FLOAT, 3, UNIT_CONV_SPEED, parts, v);
list[2] =
......@@ -108,11 +122,10 @@ void hydro_write_particles(struct part* parts, struct io_props* list,
UNIT_CONV_ACCELERATION, parts, a_hydro);
list[7] =
io_make_output_field("Density", FLOAT, 1, UNIT_CONV_DENSITY, parts, rho);
list[8] = io_make_output_field_convert_part("Entropy", FLOAT, 1,
UNIT_CONV_ENTROPY_PER_UNIT_MASS,
parts, rho, convert_S);
list[8] = io_make_output_field_convert_part(
"Entropy", FLOAT, 1, UNIT_CONV_ENTROPY_PER_UNIT_MASS, parts, convert_S);
list[9] = io_make_output_field_convert_part(
"Pressure", FLOAT, 1, UNIT_CONV_PRESSURE, parts, rho, convert_P);
"Pressure", FLOAT, 1, UNIT_CONV_PRESSURE, parts, convert_P);
}
/**
......
......@@ -67,14 +67,28 @@ void hydro_read_particles(struct part* parts, struct io_props* list,
UNIT_CONV_DENSITY, parts, rho);
}
float convert_u(struct engine* e, struct part* p) {
void convert_u(const struct engine* e, const struct part* p, float* ret) {
return hydro_get_internal_energy(p);
ret[0] = hydro_get_internal_energy(p);
}
float convert_P(struct engine* e, struct part* p) {
void convert_P(const struct engine* e, const struct part* p, float* ret) {
return hydro_get_pressure(p);
ret[0] = hydro_get_pressure(p);
}
void convert_part_pos(const struct engine* e, const struct part* p,
double* ret) {
if (e->s->periodic) {
ret[0] = box_wrap(p->x[0], 0.0, e->s->dim[0]);
ret[1] = box_wrap(p->x[1], 0.0, e->s->dim[1]);
ret[2] = box_wrap(p->x[2], 0.0, e->s->dim[2]);
} else {
ret[0] = p->x[0];
ret[1] = p->x[1];
ret[2] = p->x[2];
}
}
/**
......@@ -90,27 +104,27 @@ void hydro_write_particles(struct part* parts, struct io_props* list,
*num_fields = 11;
/* List what we want to write */
list[0] = io_make_output_field("Coordinates", DOUBLE, 3, UNIT_CONV_LENGTH,
parts, x);
list[0] = io_make_output_field_convert_part(
"Coordinates", DOUBLE, 3, UNIT_CONV_LENGTH, parts, convert_part_pos);
list[1] =
io_make_output_field("Velocities", FLOAT, 3, UNIT_CONV_SPEED, parts, v);
list[2] =
io_make_output_field("Masses", FLOAT, 1, UNIT_CONV_MASS, parts, mass);
list[3] = io_make_output_field("SmoothingLength", FLOAT, 1, UNIT_CONV_LENGTH,
parts, h);
list[4] = io_make_output_field_convert_part("InternalEnergy", FLOAT, 1,
UNIT_CONV_ENERGY_PER_UNIT_MASS,
parts, entropy, convert_u);
list[4] = io_make_output_field(
"Entropy", FLOAT, 1, UNIT_CONV_ENTROPY_PER_UNIT_MASS, parts, entropy);
list[5] = io_make_output_field("ParticleIDs", ULONGLONG, 1,
UNIT_CONV_NO_UNITS, parts, id);
list[6] = io_make_output_field("Acceleration", FLOAT, 3,
UNIT_CONV_ACCELERATION, parts, a_hydro);
list[7] =
io_make_output_field("Density", FLOAT, 1, UNIT_CONV_DENSITY, parts, rho);
list[8] = io_make_output_field(
"Entropy", FLOAT, 1, UNIT_CONV_ENTROPY_PER_UNIT_MASS, parts, entropy);
list[8] = io_make_output_field_convert_part("InternalEnergy", FLOAT, 1,
UNIT_CONV_ENERGY_PER_UNIT_MASS,
parts, convert_u);
list[9] = io_make_output_field_convert_part(
"Pressure", FLOAT, 1, UNIT_CONV_PRESSURE, parts, rho, convert_P);
"Pressure", FLOAT, 1, UNIT_CONV_PRESSURE, parts, convert_P);
list[10] = io_make_output_field("WeightedDensity", FLOAT, 1,
UNIT_CONV_DENSITY, parts, rho_bar);
}
......
......@@ -22,6 +22,7 @@
#include "approx_math.h"
#include "equation_of_state.h"
#include "hydro_gradients.h"
#include "hydro_properties.h"
#include "hydro_space.h"
#include "voronoi_algorithm.h"
......
......@@ -19,6 +19,7 @@
#include "adiabatic_index.h"
#include "equation_of_state.h"
#include "hydro.h"
#include "hydro_gradients.h"
#include "hydro_slope_limiters.h"
#include "io_properties.h"
......@@ -63,13 +64,8 @@ void hydro_read_particles(struct part* parts, struct io_props* list,
* @param p Particle.
* @return Internal energy of the particle
*/
float convert_u(struct engine* e, struct part* p) {
if (p->primitives.rho > 0.) {
return gas_internal_energy_from_pressure(p->primitives.rho,
p->primitives.P);
} else {
return 0.;
}
void convert_u(const struct engine* e, const struct part* p, float* ret) {
ret[0] = hydro_get_internal_energy(p);
}
/**
......@@ -79,12 +75,8 @@ float convert_u(struct engine* e, struct part* p) {
* @param p Particle.
* @return Entropic function of the particle
*/
float convert_A(struct engine* e, struct part* p) {
if (p->primitives.rho > 0.) {
return gas_entropy_from_pressure(p->primitives.rho, p->primitives.P);
} else {
return 0.;
}
void convert_A(const struct engine* e, const struct part* p, float* ret) {
ret[0] = hydro_get_entropy(p);
}
/**
......@@ -94,7 +86,7 @@ float convert_A(struct engine* e, struct part* p) {
* @param p Particle.
* @return Total energy of the particle
*/
float convert_Etot(struct engine* e, struct part* p) {
void convert_Etot(const struct engine* e, const struct part* p, float* ret) {
#ifdef SHADOWFAX_TOTAL_ENERGY
return p->conserved.energy;
#else
......@@ -105,13 +97,27 @@ float convert_Etot(struct engine* e, struct part* p) {
p->conserved.momentum[1] * p->conserved.momentum[1] +
p->conserved.momentum[2] * p->conserved.momentum[2];
return p->conserved.energy + 0.5f * momentum2 / p->conserved.mass;
ret[0] = p->conserved.energy + 0.5f * momentum2 / p->conserved.mass;
} else {
return 0.;
ret[0] = 0.;
}
#endif
}
void convert_part_pos(const struct engine* e, const struct part* p,
double* ret) {
if (e->s->periodic) {
ret[0] = box_wrap(p->x[0], 0.0, e->s->dim[0]);
ret[1] = box_wrap(p->x[1], 0.0, e->s->dim[1]);
ret[2] = box_wrap(p->x[2], 0.0, e->s->dim[2]);
} else {
ret[0] = p->x[0];
ret[1] = p->x[1];
ret[2] = p->x[2];
}
}
/**
* @brief Specifies which particle fields to write to a dataset
*
......@@ -125,8 +131,8 @@ void hydro_write_particles(struct part* parts, struct io_props* list,
*num_fields = 13;
/* List what we want to write */
list[0] = io_make_output_field("Coordinates", DOUBLE, 3, UNIT_CONV_LENGTH,
parts, x);
list[0] = io_make_output_field_convert_part(
"Coordinates", DOUBLE, 3, UNIT_CONV_LENGTH, parts, convert_part_pos);
list[1] = io_make_output_field("Velocities", FLOAT, 3, UNIT_CONV_SPEED, parts,
primitives.v);
list[2] = io_make_output_field("Masses", FLOAT, 1, UNIT_CONV_MASS, parts,
......@@ -135,7 +141,7 @@ void hydro_write_particles(struct part* parts, struct io_props* list,
parts, h);
list[4] = io_make_output_field_convert_part("InternalEnergy", FLOAT, 1,
UNIT_CONV_ENERGY_PER_UNIT_MASS,
parts, primitives.P, convert_u);
parts, convert_u);
list[5] = io_make_output_field("ParticleIDs", ULONGLONG, 1,
UNIT_CONV_NO_UNITS, parts, id);
list[6] = io_make_output_field("Acceleration", FLOAT, 3,
......@@ -147,12 +153,11 @@ void hydro_write_particles(struct part* parts, struct io_props* list,
list[9] = io_make_output_field("GradDensity", FLOAT, 3, UNIT_CONV_DENSITY,
parts, primitives.gradients.rho);
list[10] = io_make_output_field_convert_part(
"Entropy", FLOAT, 1, UNIT_CONV_ENTROPY, parts, primitives.P, convert_A);
"Entropy", FLOAT, 1, UNIT_CONV_ENTROPY, parts, convert_A);
list[11] = io_make_output_field("Pressure", FLOAT, 1, UNIT_CONV_PRESSURE,
parts, primitives.P);
list[12] =
io_make_output_field_convert_part("TotEnergy", FLOAT, 1, UNIT_CONV_ENERGY,
parts, conserved.energy, convert_Etot);
list[12] = io_make_output_field_convert_part(
"TotEnergy", FLOAT, 1, UNIT_CONV_ENERGY, parts, convert_Etot);
}
/**
......
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