Skip to content
GitLab
Commit 54257804 authored by Matthieu Schaller's avatar Matthieu Schaller
Browse files

Code formatting

parent 4cba6cc3
Hide whitespace changes
Inline Side-by-side
examples/main.c
View file @ 54257804
......@@ -88,8 +88,8 @@ void print_help_message() {
int main(int argc, char *argv[]) {
struct clocks_time tic, toc;
int nr_nodes = 1, myrank = 0;
int nr_nodes = 1, myrank = 0;
#ifdef WITH_MPI
/* Start by initializing MPI. */
int res = 0, prov = 0;
......@@ -258,14 +258,28 @@ int main(int argc, char *argv[]) {
MPI_Bcast(params, sizeof(struct swift_params), MPI_BYTE, 0, MPI_COMM_WORLD);
#endif
/* Initialize unit system */
/* Prepare the domain decomposition scheme */
#ifdef WITH_MPI
struct partition initial_partition;
enum repartition_type reparttype;
partition_init(&initial_partition, &reparttype, params, nr_nodes);
/* Let's report what we did */
if (myrank == 0) {
message("Using initial partition %s",
initial_partition_name[initial_partition.type]);
if (initial_partition.type == INITPART_GRID)
message("grid set to [ %i %i %i ].", initial_partition.grid[0],
initial_partition.grid[1], initial_partition.grid[2]);
message("Using %s repartitioning", repartition_name[reparttype]);
}
#endif
/* Initialize unit system and constants */
struct UnitSystem us;
struct phys_const prog_const;
struct external_potential potential;
units_init(&us, params);
initPhysicalConstants(&us, &prog_const);
if( with_external_gravity)
initPotentialProperties(params, &us, &potential);
if (myrank == 0) {
message("Unit system: U_M = %e g.", us.UnitMass_in_cgs);
message("Unit system: U_L = %e cm.", us.UnitLength_in_cgs);
......@@ -274,31 +288,18 @@ int main(int argc, char *argv[]) {
message("Unit system: U_T = %e K.", us.UnitTemperature_in_cgs);
message("Density units: %e a^%f h^%f.",
units_conversion_factor(&us, UNIT_CONV_DENSITY),
units_a_factor(&us, UNIT_CONV_DENSITY),
units_a_factor(&us, UNIT_CONV_DENSITY),
units_h_factor(&us, UNIT_CONV_DENSITY));
message("Entropy units: %e a^%f h^%f.",
units_conversion_factor(&us, UNIT_CONV_ENTROPY),
units_a_factor(&us, UNIT_CONV_ENTROPY),
units_a_factor(&us, UNIT_CONV_ENTROPY),
units_h_factor(&us, UNIT_CONV_ENTROPY));
message("Gravity constant = %e", prog_const.newton_gravity);
}
/* Prepare the domain decomposition scheme */
#ifdef WITH_MPI
struct partition initial_partition;
enum repartition_type reparttype;
partition_init(&initial_partition, &reparttype, params, nr_nodes);
/* Let's report what we did */
if (myrank == 0) {
message("Using initial partition %s",
initial_partition_name[initial_partition.type]);
if (initial_partition.type == INITPART_GRID)
message("grid set to [ %i %i %i ].", initial_partition.grid[0],
initial_partition.grid[1], initial_partition.grid[2]);
message("Using %s repartitioning", repartition_name[reparttype]);
}
#endif
/* Initialise the external potential properties */
struct external_potential potential;
if (with_external_gravity) initPotentialProperties(params, &us, &potential);
/* Read particles and space information from (GADGET) ICs */
char ICfileName[200] = "";
......@@ -475,16 +476,15 @@ int main(int argc, char *argv[]) {
int icount = 0;
space_map_cells_pre(&s, 0, &map_cellcheck, &icount);
message("map_cellcheck picked up %i parts.", icount);
}
if (myrank == 0) {
int data[2] = {s.maxdepth ,0};
int data[2] = {s.maxdepth, 0};
space_map_cells_pre(&s, 0, &map_maxdepth, data);
message("nr of cells at depth %i is %i.", data[0], data[1]);
}
//exit(-99);
// exit(-99);
/* Legend */
if (myrank == 0)
printf("# %6s %14s %14s %10s %10s %16s [%s]\n", "Step", "Time", "Time-step",
......
src/engine.c
View file @ 54257804
......@@ -1895,7 +1895,7 @@ void engine_init_particles(struct engine *e) {
/* We always have sort tasks */
mask |= 1 << task_type_sort;
mask |= 1 << task_type_init;
/* Add the tasks corresponding to hydro operations to the masks */
if ((e->policy & engine_policy_hydro) == engine_policy_hydro) {
......@@ -2343,7 +2343,7 @@ static bool hyperthreads_present(void) {
void engine_init(struct engine *e, struct space *s,
const struct swift_params *params, int nr_nodes, int nodeID,
int policy, int verbose,
int policy, int verbose,
const struct phys_const *physical_constants,
const struct external_potential *potential) {
......
src/physical_constants_cgs.h
View file @ 54257804
......@@ -23,7 +23,7 @@
/* physical constants in cgs */
#define NEWTON_GRAVITY_CGS 6.672e-8
#define SOLAR_MASS_IN_CGS 1.989e33
#define PARSEC_IN_CGS 3.0856776e18
#define PARSEC_IN_CGS 3.0856776e18
#define PROTON_MASS_IN_CGS 1.6726231e24
#define YEAR_IN_CGS 3.154e+7
......
src/potentials.c
View file @ 54257804
......@@ -31,21 +31,31 @@
* @param us The current internal system of units
* @param potential The external potential properties to initialize
*/
void initPotentialProperties(const struct swift_params * parameter_file,
struct UnitSystem* us,
void initPotentialProperties(const struct swift_params* parameter_file,
struct UnitSystem* us,
struct external_potential* potential) {
potential->point_mass.x = parser_get_param_double(parameter_file, "PointMass:position_x");
potential->point_mass.y = parser_get_param_double(parameter_file, "PointMass:position_y");
potential->point_mass.z = parser_get_param_double(parameter_file, "PointMass:position_z");
potential->point_mass.mass = parser_get_param_double(parameter_file, "PointMass:mass");
message("Point mass properties are (x,y,z) = (%e, %e, %e), M = %e", potential->point_mass.x, potential->point_mass.y, potential->point_mass.z, potential->point_mass.mass);
potential->point_mass.x =
parser_get_param_double(parameter_file, "PointMass:position_x");
potential->point_mass.y =
parser_get_param_double(parameter_file, "PointMass:position_y");
potential->point_mass.z =
parser_get_param_double(parameter_file, "PointMass:position_z");
potential->point_mass.mass =
parser_get_param_double(parameter_file, "PointMass:mass");
message("Point mass properties are (x,y,z) = (%e, %e, %e), M = %e",
potential->point_mass.x, potential->point_mass.y,
potential->point_mass.z, potential->point_mass.mass);
/* potential->point_mass.x = */
/* 50000 * PARSEC_IN_CGS / units_conversion_factor(us, UNIT_CONV_LENGTH); */
/* 50000 * PARSEC_IN_CGS / units_conversion_factor(us, UNIT_CONV_LENGTH);
*/
/* potential->point_mass.y = */
/* 50000 * PARSEC_IN_CGS / units_conversion_factor(us, UNIT_CONV_LENGTH); */
/* 50000 * PARSEC_IN_CGS / units_conversion_factor(us, UNIT_CONV_LENGTH);
*/
/* potential->point_mass.z = */
/* 50000 * PARSEC_IN_CGS / units_conversion_factor(us, UNIT_CONV_LENGTH); */
/* 50000 * PARSEC_IN_CGS / units_conversion_factor(us, UNIT_CONV_LENGTH);
*/
/* potential->point_mass.mass = */
/* 1e10 * SOLAR_MASS_IN_CGS / units_conversion_factor(us, UNIT_CONV_MASS); */
/* 1e10 * SOLAR_MASS_IN_CGS / units_conversion_factor(us, UNIT_CONV_MASS);
*/
}
src/potentials.h
View file @ 54257804
......@@ -42,9 +42,12 @@ struct external_potential {
#ifdef EXTERNAL_POTENTIAL_POINTMASS
/* #define const_unit_length_in_cgs (1e3 * PARSEC_IN_CGS) */
/* #define const_unit_mass_in_cgs (SOLAR_MASS_IN_CGS) */
/* #define External_Potential_X (50000 * PARSEC_IN_CGS / const_unit_length_in_cgs) */
/* #define External_Potential_Y (50000 * PARSEC_IN_CGS / const_unit_length_in_cgs) */
/* #define External_Potential_Z (50000 * PARSEC_IN_CGS / const_unit_length_in_cgs) */
/* #define External_Potential_X (50000 * PARSEC_IN_CGS /
* const_unit_length_in_cgs) */
/* #define External_Potential_Y (50000 * PARSEC_IN_CGS /
* const_unit_length_in_cgs) */
/* #define External_Potential_Z (50000 * PARSEC_IN_CGS /
* const_unit_length_in_cgs) */
/* #define External_Potential_Mass \ */
/* (1e10 * SOLAR_MASS_IN_CGS / const_unit_mass_in_cgs) */
......@@ -56,7 +59,7 @@ struct external_potential {
*/
__attribute__((always_inline))
INLINE static float external_gravity_pointmass_timestep(
const struct external_potential* potential,
const struct external_potential* potential,
const struct phys_const* const phys_const,
const struct gpart* const g) {
......@@ -68,12 +71,12 @@ __attribute__((always_inline))
const float drdv = (g->x[0] - potential->point_mass.x) * (g->v_full[0]) +
(g->x[1] - potential->point_mass.y) * (g->v_full[1]) +
(g->x[2] - potential->point_mass.z) * (g->v_full[2]);
const float dota_x = G_newton * potential->point_mass.mass * rinv * rinv * rinv *
(-g->v_full[0] + 3.f * rinv * rinv * drdv * dx);
const float dota_y = G_newton * potential->point_mass.mass * rinv * rinv * rinv *
(-g->v_full[1] + 3.f * rinv * rinv * drdv * dy);
const float dota_z = G_newton * potential->point_mass.mass * rinv * rinv * rinv *
(-g->v_full[2] + 3.f * rinv * rinv * drdv * dz);
const float dota_x = G_newton * potential->point_mass.mass * rinv * rinv *
rinv * (-g->v_full[0] + 3.f * rinv * rinv * drdv * dx);
const float dota_y = G_newton * potential->point_mass.mass * rinv * rinv *
rinv * (-g->v_full[1] + 3.f * rinv * rinv * drdv * dy);
const float dota_z = G_newton * potential->point_mass.mass * rinv * rinv *
rinv * (-g->v_full[2] + 3.f * rinv * rinv * drdv * dz);
const float dota_2 = dota_x * dota_x + dota_y * dota_y + dota_z * dota_z;
const float a_2 = g->a_grav[0] * g->a_grav[0] + g->a_grav[1] * g->a_grav[1] +
g->a_grav[2] * g->a_grav[2];
......@@ -88,8 +91,9 @@ __attribute__((always_inline))
* @param phys_cont The physical constants in internal units.
* @param gp Pointer to the g-particle data.
*/
__attribute__((always_inline)) INLINE static void external_gravity_pointmass(const struct external_potential* potential, const struct phys_const* const phys_const, struct gpart* g) {
__attribute__((always_inline)) INLINE static void external_gravity_pointmass(
const struct external_potential* potential,
const struct phys_const* const phys_const, struct gpart* g) {
/* message(" (x,y,z) = (%e, %e, %e), M= %e", potential->point_mass.x, */
/* potential->point_mass.y, potential->point_mass.z, */
......@@ -102,9 +106,12 @@ __attribute__((always_inline)) INLINE static void external_gravity_pointmass(con
const float dz = g->x[2] - potential->point_mass.z;
const float rinv = 1.f / sqrtf(dx * dx + dy * dy + dz * dz);
g->a_grav[0] += -G_newton * potential->point_mass.mass * dx * rinv * rinv * rinv;
g->a_grav[1] += -G_newton * potential->point_mass.mass * dy * rinv * rinv * rinv;
g->a_grav[2] += -G_newton * potential->point_mass.mass * dz * rinv * rinv * rinv;
g->a_grav[0] +=
-G_newton * potential->point_mass.mass * dx * rinv * rinv * rinv;
g->a_grav[1] +=
-G_newton * potential->point_mass.mass * dy * rinv * rinv * rinv;
g->a_grav[2] +=
-G_newton * potential->point_mass.mass * dz * rinv * rinv * rinv;
}
#endif /* EXTERNAL_POTENTIAL_POINTMASS */
......@@ -115,8 +122,8 @@ __attribute__((always_inline)) INLINE static void external_gravity_pointmass(con
* @param us The current internal system of units
* @param potential The external potential properties to initialize
*/
void initPotentialProperties(const struct swift_params * parameter_file,
struct UnitSystem* us,
void initPotentialProperties(const struct swift_params* parameter_file,
struct UnitSystem* us,
struct external_potential* potential);
#endif /* SWIFT_POTENTIALS_H */
src/runner.c
View file @ 54257804
......@@ -120,7 +120,6 @@ void runner_dograv_external(struct runner *r, struct cell *c) {
/* CentreOfPotential[1] = 0.5 * s->dim[1]; */
/* CentreOfPotential[2] = 0.5 * s->dim[2]; */
/* Recurse? */
if (c->split) {
for (k = 0; k < 8; k++)
......@@ -151,7 +150,6 @@ void runner_dograv_external(struct runner *r, struct cell *c) {
const float dr = sqrtf((dx * dx) + (dy * dy) + (dz * dz));
const float rinv = 1.f / sqrtf(dx * dx + dy * dy + dz * dz);
const int current_dti = g->ti_end - g->ti_begin;
const float dt = 0.5f * current_dti * r->e->timeBase;
const float vx = g->v_full[0] + dt * g->a_grav[0];
......@@ -160,13 +158,15 @@ void runner_dograv_external(struct runner *r, struct cell *c) {
/* E/L */
E = 0.5 * ((vx * vx) + (vy * vy) + (vz * vz)) -
r->e->physical_constants->newton_gravity * r->e->potential->point_mass.mass * rinv;
r->e->physical_constants->newton_gravity *
r->e->potential->point_mass.mass * rinv;
L[0] = dy * vz - dz * vy;
L[1] = dz * vx - dx * vz;
L[2] = dx * vy - dy * vx;
if (abs(g->id) == 1) {
float v2 = vx * vx + vy * vy + vz * vz;
float fg = r->e->physical_constants->newton_gravity * r->e->potential->point_mass.mass * rinv;
float fg = r->e->physical_constants->newton_gravity *
r->e->potential->point_mass.mass * rinv;
float fga = sqrtf((g->a_grav[0] * g->a_grav[0]) +
(g->a_grav[1] * g->a_grav[1]) +
(g->a_grav[2] * g->a_grav[2])) *
......@@ -174,10 +174,12 @@ void runner_dograv_external(struct runner *r, struct cell *c) {
// message("grav_external time= %f\t V_c^2= %f GM/r= %f E= %f L[2]= %f
// x= %f y= %f vx= %f vy= %f\n", r->e->time, v2, fg, E, L[2], g->x[0],
// g->x[1], vx, vy);
message("%f\t %f %f %f %f %f %f %f %f %f %f %f %f %f\n", r->e->time, g->tx,
g->tv, dt, v2, fg, fga, dr, E, L[2], g->x[0], g->x[1], vx, vy);
/* message(" G=%e M=%e\n", r->e->physical_constants->newton_gravity, r->e->potential->point_mass.mass); */
/* exit(-1); */
message("%f\t %f %f %f %f %f %f %f %f %f %f %f %f %f\n", r->e->time,
g->tx, g->tv, dt, v2, fg, fga, dr, E, L[2], g->x[0], g->x[1],
vx, vy);
/* message(" G=%e M=%e\n", r->e->physical_constants->newton_gravity,
* r->e->potential->point_mass.mass); */
/* exit(-1); */
}
}
}
......@@ -940,7 +942,8 @@ void runner_dokick(struct runner *r, struct cell *c, int timer) {
} else {
/* Compute the next timestep (gravity condition) */
float new_dt = gravity_compute_timestep(r->e->potential, r->e->physical_constants, gp);
float new_dt = gravity_compute_timestep(r->e->potential,
r->e->physical_constants, gp);
/* Limit timestep within the allowed range */
new_dt = fminf(new_dt, global_dt_max);
......@@ -1025,8 +1028,8 @@ void runner_dokick(struct runner *r, struct cell *c, int timer) {
/* Compute the next timestep (gravity condition) */
float new_dt_grav = FLT_MAX;
if (p->gpart != NULL)
new_dt_grav =
gravity_compute_timestep(r->e->potential, r->e->physical_constants, p->gpart);
new_dt_grav = gravity_compute_timestep(
r->e->potential, r->e->physical_constants, p->gpart);
float new_dt = fminf(new_dt_hydro, new_dt_grav);
......
src/space.c
View file @ 54257804
......@@ -1348,9 +1348,9 @@ void space_init(struct space *s, const struct swift_params *params,
space_maxsize = parser_get_param_int(params, "Scheduler:cell_max_size");
space_subsize = parser_get_param_int(params, "Scheduler:cell_sub_size");
space_splitsize = parser_get_param_int(params, "Scheduler:cell_split_size");
if(verbose)
if (verbose)
message("max_size set to %d, sub_size set to %d, split_size set to %d",
space_maxsize, space_subsize, space_splitsize);
space_maxsize, space_subsize, space_splitsize);
/* Check that we have enough cells */
if (s->cell_min * 3 > dim[0] || s->cell_min * 3 > dim[1] ||
......
src/task.c
View file @ 54257804
......@@ -47,9 +47,9 @@
/* Task type names. */
const char *taskID_names[task_type_count] = {
"none", "sort", "self", "pair", "sub",
"init", "ghost", "drift", "kick", "send",
"recv", "grav_pp", "grav_mm", "grav_up", "grav_down",
"none", "sort", "self", "pair", "sub",
"init", "ghost", "drift", "kick", "send",
"recv", "grav_pp", "grav_mm", "grav_up", "grav_down",
"grav_external", "part_sort", "gpart_sort", "split_cell", "rewait"};
const char *subtaskID_names[task_type_count] = {"none", "density",
......
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