tx258/lammps-gran-kokkos
1
0
Fork 0
Code Issues Pull Requests Actions Packages Projects Releases Wiki Activity

fixes for OPENMP versions of dielectric pair styles

Browse Source
This commit is contained in:
Axel Kohlmeyer
2022-08-28 05:05:28 -04:00
parent 2a6bd1aa6b
commit 0e5c758fb8
4 changed files with 12 additions and 33 deletions
Download Patch File Download Diff File Expand all files Collapse all files

3
src/DIELECTRIC/pair_lj_cut_coul_cut_dielectric.cpp
View File

@ -93,15 +93,16 @@ void PairLJCutCoulCutDielectric::compute(int eflag, int vflag)
for (ii = 0; ii < inum; ii++) {
i = ilist[ii];
qtmp = q[i];
etmp = eps[i];
xtmp = x[i][0];
ytmp = x[i][1];
ztmp = x[i][2];
etmp = eps[i];
itype = type[i];
jlist = firstneigh[i];
jnum = numneigh[i];
// self term Eq. (55) for I_{ii} and Eq. (52) and in Barros et al
double curvature_threshold = sqrt(area[i]);
if (curvature[i] < curvature_threshold) {
double sf = curvature[i] / (4.0 * MY_PIS * curvature_threshold) * area[i] * q[i];

21
src/OPENMP/pair_lj_cut_coul_cut_dielectric_omp.cpp
View File

@ -99,7 +99,7 @@ void PairLJCutCoulCutDielectricOMP::eval(int iifrom, int iito, ThrData *const th
{
int i, j, ii, jj, jnum, itype, jtype;
double qtmp, etmp, xtmp, ytmp, ztmp, delx, dely, delz, evdwl, ecoul;
double fpair_i, fpair_j;
double fpair_i;
double rsq, r2inv, r6inv, forcecoul, forcelj, factor_coul, factor_lj;
double efield_i, epot_i;
int *ilist, *jlist, *numneigh, **firstneigh;
@ -167,9 +167,10 @@ void PairLJCutCoulCutDielectricOMP::eval(int iifrom, int iito, ThrData *const th
if (rsq < cutsq[itype][jtype]) {
r2inv = 1.0 / rsq;
const double rinv = sqrt(r2inv);
if (rsq < cut_coulsq[itype][jtype] && rsq > EPSILON) {
efield_i = qqrd2e * q[j] * sqrt(r2inv);
efield_i = qqrd2e * q[j] * rinv;
forcecoul = qtmp * efield_i;
epot_i = efield_i;
} else
@ -182,7 +183,6 @@ void PairLJCutCoulCutDielectricOMP::eval(int iifrom, int iito, ThrData *const th
forcelj = 0.0;
fpair_i = (factor_coul * etmp * forcecoul + factor_lj * forcelj) * r2inv;
fxtmp += delx * fpair_i;
fytmp += dely * fpair_i;
fztmp += delz * fpair_i;
@ -193,29 +193,18 @@ void PairLJCutCoulCutDielectricOMP::eval(int iifrom, int iito, ThrData *const th
eztmp += delz * efield_i;
epot[i] += epot_i;
if (NEWTON_PAIR || j >= nlocal) {
fpair_j = (factor_coul * eps[j] * forcecoul + factor_lj * forcelj) * r2inv;
f[j].x -= delx * fpair_j;
f[j].y -= dely * fpair_j;
f[j].z -= delz * fpair_j;
}
if (EFLAG) {
if (rsq < cut_coulsq[itype][jtype]) {
ecoul = factor_coul * qqrd2e * qtmp * q[j] * (etmp + eps[j]) * sqrt(r2inv);
ecoul = factor_coul * qqrd2e * qtmp * q[j] * 0.5 * (etmp + eps[j]) * rinv;
} else
ecoul = 0.0;
ecoul *= 0.5;
if (rsq < cut_ljsq[itype][jtype]) {
evdwl = r6inv * (lj3[itype][jtype] * r6inv - lj4[itype][jtype]) - offset[itype][jtype];
evdwl *= factor_lj;
} else
evdwl = 0.0;
}
if (EVFLAG)
ev_tally_thr(this, i, j, nlocal, NEWTON_PAIR, evdwl, ecoul, fpair_i, delx, dely, delz,
thr);
if (EVFLAG) ev_tally_full_thr(this, i, evdwl, ecoul, fpair_i, delx, dely, delz, thr);
}
}
f[i].x += fxtmp;

17
src/OPENMP/pair_lj_cut_coul_debye_dielectric_omp.cpp
View File

@ -99,7 +99,7 @@ void PairLJCutCoulDebyeDielectricOMP::eval(int iifrom, int iito, ThrData *const
{
int i, j, ii, jj, jnum, itype, jtype;
double qtmp, etmp, xtmp, ytmp, ztmp, delx, dely, delz, evdwl, ecoul;
double fpair_i, fpair_j;
double fpair_i;
double rsq, r2inv, r6inv, forcecoul, forcelj, factor_coul, factor_lj, efield_i, epot_i;
double r, rinv, screening;
int *ilist, *jlist, *numneigh, **firstneigh;
@ -185,7 +185,6 @@ void PairLJCutCoulDebyeDielectricOMP::eval(int iifrom, int iito, ThrData *const
forcelj = 0.0;
fpair_i = (factor_coul * etmp * forcecoul + factor_lj * forcelj) * r2inv;
fxtmp += delx * fpair_i;
fytmp += dely * fpair_i;
fztmp += delz * fpair_i;
@ -196,19 +195,11 @@ void PairLJCutCoulDebyeDielectricOMP::eval(int iifrom, int iito, ThrData *const
eztmp += delz * efield_i;
epot[i] += epot_i;
if (NEWTON_PAIR || j >= nlocal) {
fpair_j = (factor_coul * eps[j] * forcecoul + factor_lj * forcelj) * r2inv;
f[j].x -= delx * fpair_j;
f[j].y -= dely * fpair_j;
f[j].z -= delz * fpair_j;
}
if (EFLAG) {
if (rsq < cut_coulsq[itype][jtype]) {
ecoul = factor_coul * qqrd2e * qtmp * q[j] * (etmp + eps[j]) * rinv * screening;
ecoul = factor_coul * qqrd2e * qtmp * q[j] * 0.5 * (etmp + eps[j]) * rinv * screening;
} else
ecoul = 0.0;
ecoul *= 0.5;
if (rsq < cut_ljsq[itype][jtype]) {
evdwl = r6inv * (lj3[itype][jtype] * r6inv - lj4[itype][jtype]) - offset[itype][jtype];
evdwl *= factor_lj;
@ -216,9 +207,7 @@ void PairLJCutCoulDebyeDielectricOMP::eval(int iifrom, int iito, ThrData *const
evdwl = 0.0;
}
if (EVFLAG)
ev_tally_thr(this, i, j, nlocal, NEWTON_PAIR, evdwl, ecoul, fpair_i, delx, dely, delz,
thr);
if (EVFLAG) ev_tally_full_thr(this, i, evdwl, ecoul, fpair_i, delx, dely, delz, thr);
}
}
f[i].x += fxtmp;

4
src/OPENMP/pair_lj_cut_coul_long_dielectric_omp.cpp
View File

@ -229,10 +229,10 @@ void PairLJCutCoulLongDielectricOMP::eval(int iifrom, int iito, ThrData *const t
if (EFLAG) {
if (rsq < cut_coulsq) {
if (!ncoultablebits || rsq <= tabinnersq)
ecoul = prefactor * (etmp + eps[j]) * erfc;
ecoul = prefactor * 0.5 * (etmp + eps[j]) * erfc;
else {
table = etable[itable] + fraction * detable[itable];
ecoul = qtmp * q[j] * (etmp + eps[j]) * table;
ecoul = qtmp * q[j] * 0.5 * (etmp + eps[j]) * table;
}
if (factor_coul < 1.0) ecoul -= (1.0 - factor_coul) * prefactor;
} else