Commit 810d3499 authored by Matthieu Schaller's avatar Matthieu Schaller
Browse files

a --> a_hydro everywhere

parent 51362dbb
......@@ -25,6 +25,7 @@ __attribute__((always_inline))
"h=%.3e, "
"wcount=%d, m=%.3e, dh_drho=%.3e, rho=%.3e, t_begin=%.3e, t_end=%.3e\n",
p->x[0], p->x[1], p->x[2], p->v[0], p->v[1], p->v[2], xp->v_full[0],
xp->v_full[1], xp->v_full[2], p->a[0], p->a[1], p->a[2], 2. * p->h,
(int)p->density.wcount, p->mass, p->rho_dh, p->rho, p->t_begin, p->t_end);
xp->v_full[1], xp->v_full[2], p->a_hydro[0], p->a_hydro[1], p->a_hydro[2],
p->h, (int)p->density.wcount, p->mass, p->rho_dh, p->rho, p->t_begin,
p->t_end);
}
......@@ -45,7 +45,7 @@ __attribute__((always_inline)) INLINE static void hydro_read_particles(
COMPULSORY);
readArray(h_grp, "ParticleIDs", ULONGLONG, N, 1, parts, N_total, offset, id,
COMPULSORY);
readArray(h_grp, "Acceleration", FLOAT, N, 3, parts, N_total, offset, a,
readArray(h_grp, "Acceleration", FLOAT, N, 3, parts, N_total, offset, a_hydro,
OPTIONAL);
readArray(h_grp, "Density", FLOAT, N, 1, parts, N_total, offset, rho,
OPTIONAL);
......@@ -84,7 +84,7 @@ __attribute__((always_inline)) INLINE static void hydro_write_particles(
writeArray(h_grp, fileName, xmfFile, "ParticleIDs", ULONGLONG, N, 1, parts,
N_total, mpi_rank, offset, id, us, UNIT_CONV_NO_UNITS);
writeArray(h_grp, fileName, xmfFile, "Acceleration", FLOAT, N, 3, parts,
N_total, mpi_rank, offset, a, us, UNIT_CONV_ACCELERATION);
N_total, mpi_rank, offset, a_hydro, us, UNIT_CONV_ACCELERATION);
writeArray(h_grp, fileName, xmfFile, "Density", FLOAT, N, 1, parts, N_total,
mpi_rank, offset, rho, us, UNIT_CONV_DENSITY);
}
......@@ -44,7 +44,7 @@ struct part {
float v[3];
/* Particle acceleration. */
float a[3];
float a_hydro[3];
/* Particle cutoff radius. */
float h;
......
......@@ -28,7 +28,9 @@ __attribute__((always_inline)) INLINE static float hydro_compute_timestep(
struct part* p, struct xpart* xp) {
/* Acceleration */
float ac = sqrtf(p->a[0] * p->a[0] + p->a[1] * p->a[1] + p->a[2] * p->a[2]);
float ac = sqrtf(p->a_hydro[0] * p->a_hydro[0] +
p->a_hydro[1] * p->a_hydro[1] +
p->a_hydro[2] * p->a_hydro[2]);
ac = fmaxf(ac, 1e-30);
const float dt_accel = sqrtf(2.f); // MATTHIEU
......@@ -164,9 +166,9 @@ __attribute__((always_inline))
INLINE static void hydro_reset_acceleration(struct part* p) {
/* Reset the acceleration. */
p->a[0] = 0.0f;
p->a[1] = 0.0f;
p->a[2] = 0.0f;
p->a_hydro[0] = 0.0f;
p->a_hydro[1] = 0.0f;
p->a_hydro[2] = 0.0f;
p->force.h_dt = 0.0f;
......
......@@ -28,8 +28,8 @@ __attribute__((always_inline))
"divV=%.3e, curlV=%.3e, rotV=[%.3e,%.3e,%.3e] \n "
"v_sig=%e t_begin=%.3e, t_end=%.3e\n",
p->x[0], p->x[1], p->x[2], p->v[0], p->v[1], p->v[2], xp->v_full[0],
xp->v_full[1], xp->v_full[2], p->a[0], p->a[1], p->a[2], 2. * p->h,
(int)p->density.wcount, p->mass, p->rho_dh, p->rho, p->pressure,
xp->v_full[1], xp->v_full[2], p->a_hydro[0], p->a_hydro[1], p->a_hydro[2],
p->h, (int)p->density.wcount, p->mass, p->rho_dh, p->rho, p->pressure,
p->entropy, p->entropy_dt, p->div_v, p->curl_v, p->rot_v[0], p->rot_v[1],
p->rot_v[2], p->v_sig, p->t_begin, p->t_end);
}
......@@ -272,13 +272,13 @@ __attribute__((always_inline)) INLINE static void runner_iact_force(
/* message("oO"); */
/* Use the force Luke ! */
pi->a[0] -= acc * dx[0];
pi->a[1] -= acc * dx[1];
pi->a[2] -= acc * dx[2];
pi->a_hydro[0] -= acc * dx[0];
pi->a_hydro[1] -= acc * dx[1];
pi->a_hydro[2] -= acc * dx[2];
pj->a[0] += acc * dx[0];
pj->a[1] += acc * dx[1];
pj->a[2] += acc * dx[2];
pj->a_hydro[0] += acc * dx[0];
pj->a_hydro[1] += acc * dx[1];
pj->a_hydro[2] += acc * dx[2];
/* Get the time derivative for h. */
pi->force.h_dt -= mj * dvdr * r_inv / rhoj * wi_dr;
......@@ -366,9 +366,9 @@ __attribute__((always_inline)) INLINE static void runner_iact_nonsym_force(
const float acc = visc_term + sph_term;
/* Use the force Luke ! */
pi->a[0] -= acc * dx[0];
pi->a[1] -= acc * dx[1];
pi->a[2] -= acc * dx[2];
pi->a_hydro[0] -= acc * dx[0];
pi->a_hydro[1] -= acc * dx[1];
pi->a_hydro[2] -= acc * dx[2];
/* Get the time derivative for h. */
pi->force.h_dt -= mj * dvdr * r_inv / rhoj * wi_dr;
......
......@@ -45,7 +45,7 @@ __attribute__((always_inline)) INLINE static void hydro_read_particles(
entropy, COMPULSORY);
readArray(h_grp, "ParticleIDs", ULONGLONG, N, 1, parts, N_total, offset, id,
COMPULSORY);
readArray(h_grp, "Acceleration", FLOAT, N, 3, parts, N_total, offset, a,
readArray(h_grp, "Acceleration", FLOAT, N, 3, parts, N_total, offset, a_hydro,
OPTIONAL);
readArray(h_grp, "Density", FLOAT, N, 1, parts, N_total, offset, rho,
OPTIONAL);
......@@ -85,7 +85,7 @@ __attribute__((always_inline)) INLINE static void hydro_write_particles(
writeArray(h_grp, fileName, xmfFile, "ParticleIDs", ULONGLONG, N, 1, parts,
N_total, mpi_rank, offset, id, us, UNIT_CONV_NO_UNITS);
writeArray(h_grp, fileName, xmfFile, "Acceleration", FLOAT, N, 3, parts,
N_total, mpi_rank, offset, a, us, UNIT_CONV_ACCELERATION);
N_total, mpi_rank, offset, a_hydro, us, UNIT_CONV_ACCELERATION);
writeArray(h_grp, fileName, xmfFile, "Density", FLOAT, N, 1, parts, N_total,
mpi_rank, offset, rho, us, UNIT_CONV_DENSITY);
}
......@@ -41,7 +41,7 @@ struct part {
float v[3];
/* Particle acceleration. */
float a[3];
float a_hydro[3];
/* Particle cutoff radius. */
float h;
......
......@@ -729,9 +729,9 @@ void runner_dodrift(struct runner *r, struct cell *c, int timer) {
p->x[2] += xp->v_full[2] * dt;
/* Predict velocities (for hydro terms) */
p->v[0] += p->a[0] * dt;
p->v[1] += p->a[1] * dt;
p->v[2] += p->a[2] * dt;
p->v[0] += p->a_hydro[0] * dt;
p->v[1] += p->a_hydro[1] * dt;
p->v[2] += p->a_hydro[2] * dt;
/* Predict smoothing length */
w = p->force.h_dt * h_inv * dt;
......@@ -901,13 +901,13 @@ void runner_dokick(struct runner *r, struct cell *c, int timer) {
p->t_end = p->t_begin + new_dt;
/* Kick particles in momentum space */
xp->v_full[0] += p->a[0] * dt;
xp->v_full[1] += p->a[1] * dt;
xp->v_full[2] += p->a[2] * dt;
xp->v_full[0] += p->a_hydro[0] * dt;
xp->v_full[1] += p->a_hydro[1] * dt;
xp->v_full[2] += p->a_hydro[2] * dt;
p->v[0] = xp->v_full[0] - half_dt * p->a[0];
p->v[1] = xp->v_full[1] - half_dt * p->a[1];
p->v[2] = xp->v_full[2] - half_dt * p->a[2];
p->v[0] = xp->v_full[0] - half_dt * p->a_hydro[0];
p->v[1] = xp->v_full[1] - half_dt * p->a_hydro[1];
p->v[2] = xp->v_full[2] - half_dt * p->a_hydro[2];
/* if(p->id == 1000 || p->id == 515050 || p->id == 504849) */
/* message("%lld: current_t=%f t_beg=%f t_end=%f half_dt=%f v=[%.3e
......
......@@ -420,7 +420,7 @@ void engine_single_force(double *dim, long long int pid,
}
/* Dump the result. */
message("part %lli (h=%e) has a=[%.3e,%.3e,%.3e]", p.id, p.h, p.a[0], p.a[1],
p.a[2]);
message("part %lli (h=%e) has a=[%.3e,%.3e,%.3e]", p.id, p.h, p.a_hydro[0], p.a_hydro[1],
p.a_hydro[2]);
fflush(stdout);
}
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