Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

improve OpenMP usage #2

Open
wants to merge 2 commits into
base: master
Choose a base branch
from
Open
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
60 changes: 44 additions & 16 deletions src/Hydro.cc
Original file line number Diff line number Diff line change
Expand Up @@ -103,7 +103,9 @@ void Hydro::init() {
cftot = Memory::alloc<double2>(nums);

// initialize hydro vars
#pragma omp parallel for schedule(static)
#pragma omp parallel
{
#pragma omp for schedule(static)
for (int zch = 0; zch < numzch; ++zch) {
int zfirst = mesh->zchzfirst[zch];
int zlast = mesh->zchzlast[zch];
Expand All @@ -116,6 +118,7 @@ void Hydro::init() {
if (!subrgn.empty()) {
const double eps = 1.e-12;
#pragma ivdep
#pragma omp simd
for (int z = zfirst; z < zlast; ++z) {
if (zx[z].x > (subrgn[0] - eps) &&
zx[z].x < (subrgn[1] + eps) &&
Expand All @@ -128,13 +131,14 @@ void Hydro::init() {
}

#pragma ivdep
#pragma omp simd
for (int z = zfirst; z < zlast; ++z) {
zm[z] = zr[z] * zvol[z];
zetot[z] = ze[z] * zm[z];
}
} // for sch

#pragma omp parallel for schedule(static)
#pragma omp for schedule(static)
for (int pch = 0; pch < numpch; ++pch) {
int pfirst = mesh->pchpfirst[pch];
int plast = mesh->pchplast[pch];
Expand All @@ -144,8 +148,11 @@ void Hydro::init() {
fill(&pu[pfirst], &pu[plast], double2(0., 0.));
} // for pch

resetDtHydro();

#pragma omp single
{
resetDtHydro();
}
} // omp parallel
}


Expand All @@ -157,6 +164,7 @@ void Hydro::initRadialVel(
const double eps = 1.e-12;

#pragma ivdep
#pragma omp simd
for (int p = pfirst; p < plast; ++p) {
double pmag = length(px[p]);
if (pmag > eps)
Expand Down Expand Up @@ -194,7 +202,9 @@ void Hydro::doCycle(
double* zdl = mesh->zdl;

// Begin hydro cycle
#pragma omp parallel for schedule(static)
#pragma omp parallel
{
#pragma omp for schedule(static)
for (int pch = 0; pch < numpch; ++pch) {
int pfirst = mesh->pchpfirst[pch];
int plast = mesh->pchplast[pch];
Expand All @@ -208,7 +218,7 @@ void Hydro::doCycle(
advPosHalf(px0, pu0, dt, pxp, pfirst, plast);
} // for pch

#pragma omp parallel for schedule(static)
#pragma omp for schedule(static)
for (int sch = 0; sch < numsch; ++sch) {
int sfirst = mesh->schsfirst[sch];
int slast = mesh->schslast[sch];
Expand Down Expand Up @@ -241,13 +251,17 @@ void Hydro::doCycle(
qcs->calcForce(sfq, sfirst, slast);
sumCrnrForce(sfp, sfq, sft, cftot, sfirst, slast);
} // for sch
mesh->checkBadSides();

// sum corner masses, forces to points
mesh->sumToPoints(cmaswt, pmaswt);
mesh->sumToPoints(cftot, pf);
#pragma omp single
{
mesh->checkBadSides();

#pragma omp parallel for schedule(static)
// sum corner masses, forces to points
mesh->sumToPoints(cmaswt, pmaswt);
mesh->sumToPoints(cftot, pf);
}

#pragma omp for schedule(static)
for (int pch = 0; pch < numpch; ++pch) {
int pfirst = mesh->pchpfirst[pch];
int plast = mesh->pchplast[pch];
Expand All @@ -266,10 +280,12 @@ void Hydro::doCycle(
// 6. advance mesh to end of time step
advPosFull(px0, pu0, pap, dt, px, pu, pfirst, plast);
} // for pch
#pragma omp single
{
resetDtHydro();
}

resetDtHydro();

#pragma omp parallel for schedule(static)
#pragma omp for schedule(static)
for (int sch = 0; sch < numsch; ++sch) {
int sfirst = mesh->schsfirst[sch];
int slast = mesh->schslast[sch];
Expand All @@ -286,9 +302,12 @@ void Hydro::doCycle(
calcWork(sfp, sfq, pu0, pu, pxp, dt, zw, zetot,
sfirst, slast);
} // for sch
mesh->checkBadSides();
#pragma omp single
{
mesh->checkBadSides();
}

#pragma omp parallel for schedule(static)
#pragma omp for schedule(static)
for (int zch = 0; zch < mesh->numzch; ++zch) {
int zfirst = mesh->zchzfirst[zch];
int zlast = mesh->zchzlast[zch];
Expand All @@ -304,6 +323,7 @@ void Hydro::doCycle(
calcDtHydro(zdl, zvol, zvol0, dt, zfirst, zlast);
} // for zch

} // omp parallel
}


Expand All @@ -318,6 +338,7 @@ void Hydro::advPosHalf(
double dth = 0.5 * dt;

#pragma ivdep
#pragma omp simd
for (int p = pfirst; p < plast; ++p) {
pxp[p] = px0[p] + pu0[p] * dth;
}
Expand All @@ -335,6 +356,7 @@ void Hydro::advPosFull(
const int plast) {

#pragma ivdep
#pragma omp simd
for (int p = pfirst; p < plast; ++p) {
pu[p] = pu0[p] + pa[p] * dt;
px[p] = px0[p] + 0.5 * (pu[p] + pu0[p]) * dt;
Expand All @@ -352,6 +374,7 @@ void Hydro::calcCrnrMass(
const int slast) {

#pragma ivdep
#pragma omp simd
for (int s = sfirst; s < slast; ++s) {
int s3 = mesh->mapss3[s];
int z = mesh->mapsz[s];
Expand All @@ -371,6 +394,7 @@ void Hydro::sumCrnrForce(
const int slast) {

#pragma ivdep
#pragma omp simd
for (int s = sfirst; s < slast; ++s) {
int s3 = mesh->mapss3[s];

Expand All @@ -391,6 +415,7 @@ void Hydro::calcAccel(
const double fuzz = 1.e-99;

#pragma ivdep
#pragma omp simd
for (int p = pfirst; p < plast; ++p) {
pa[p] = pf[p] / max(pmass[p], fuzz);
}
Expand All @@ -406,6 +431,7 @@ void Hydro::calcRho(
const int zlast) {

#pragma ivdep
#pragma omp simd
for (int z = zfirst; z < zlast; ++z) {
zr[z] = zm[z] / zvol[z];
}
Expand Down Expand Up @@ -461,6 +487,7 @@ void Hydro::calcWorkRate(
const int zlast) {
double dtinv = 1. / dt;
#pragma ivdep
#pragma omp simd
for (int z = zfirst; z < zlast; ++z) {
double dvol = zvol[z] - zvol0[z];
zwrate[z] = (zw[z] + zp[z] * dvol) * dtinv;
Expand All @@ -478,6 +505,7 @@ void Hydro::calcEnergy(

const double fuzz = 1.e-99;
#pragma ivdep
#pragma omp simd
for (int z = zfirst; z < zlast; ++z) {
ze[z] = zetot[z] / (zm[z] + fuzz);
}
Expand Down