From 5e096a7f9fad6e2bc6bcc740d6c5c739b003fd49 Mon Sep 17 00:00:00 2001 From: Axel Kohlmeyer Date: Sun, 11 Jul 2021 16:03:42 -0400 Subject: [PATCH] correct ev_tally4() for pair_modify nofdotr --- src/OPENMP/thr_omp.cpp | 31 ++++++++++++++++--------- src/pair.cpp | 51 ++++++++++++++++++++++++++++-------------- 2 files changed, 54 insertions(+), 28 deletions(-) diff --git a/src/OPENMP/thr_omp.cpp b/src/OPENMP/thr_omp.cpp index 397b70a3d0..545e3bbe88 100644 --- a/src/OPENMP/thr_omp.cpp +++ b/src/OPENMP/thr_omp.cpp @@ -725,18 +725,27 @@ void ThrOMP::ev_tally4_thr(Pair * const pair, const int i, const int j, } } - if (pair->vflag_atom) { - v[0] = 0.25 * (drim[0]*fi[0] + drjm[0]*fj[0] + drkm[0]*fk[0]); - v[1] = 0.25 * (drim[1]*fi[1] + drjm[1]*fj[1] + drkm[1]*fk[1]); - v[2] = 0.25 * (drim[2]*fi[2] + drjm[2]*fj[2] + drkm[2]*fk[2]); - v[3] = 0.25 * (drim[0]*fi[1] + drjm[0]*fj[1] + drkm[0]*fk[1]); - v[4] = 0.25 * (drim[0]*fi[2] + drjm[0]*fj[2] + drkm[0]*fk[2]); - v[5] = 0.25 * (drim[1]*fi[2] + drjm[1]*fj[2] + drkm[1]*fk[2]); + if (pair->vflag_either) { + v[0] = (drim[0]*fi[0] + drjm[0]*fj[0] + drkm[0]*fk[0]); + v[1] = (drim[1]*fi[1] + drjm[1]*fj[1] + drkm[1]*fk[1]); + v[2] = (drim[2]*fi[2] + drjm[2]*fj[2] + drkm[2]*fk[2]); + v[3] = (drim[0]*fi[1] + drjm[0]*fj[1] + drkm[0]*fk[1]); + v[4] = (drim[0]*fi[2] + drjm[0]*fj[2] + drkm[0]*fk[2]); + v[5] = (drim[1]*fi[2] + drjm[1]*fj[2] + drkm[1]*fk[2]); + if (pair->vflag_global) v_tally(thr->virial_pair,v); - v_tally(thr->vatom_pair[i],v); - v_tally(thr->vatom_pair[j],v); - v_tally(thr->vatom_pair[k],v); - v_tally(thr->vatom_pair[m],v); + if (pair->vflag_atom) { + v[0] *= 0.25; + v[1] *= 0.25; + v[2] *= 0.25; + v[3] *= 0.25; + v[4] *= 0.25; + v[5] *= 0.25; + v_tally(thr->vatom_pair[i],v); + v_tally(thr->vatom_pair[j],v); + v_tally(thr->vatom_pair[k],v); + v_tally(thr->vatom_pair[m],v); + } } } diff --git a/src/pair.cpp b/src/pair.cpp index 912b071490..5e8c666743 100644 --- a/src/pair.cpp +++ b/src/pair.cpp @@ -1363,22 +1363,40 @@ void Pair::ev_tally4(int i, int j, int k, int m, double evdwl, } } - if (vflag_atom) { - v[0] = 0.25 * (drim[0]*fi[0] + drjm[0]*fj[0] + drkm[0]*fk[0]); - v[1] = 0.25 * (drim[1]*fi[1] + drjm[1]*fj[1] + drkm[1]*fk[1]); - v[2] = 0.25 * (drim[2]*fi[2] + drjm[2]*fj[2] + drkm[2]*fk[2]); - v[3] = 0.25 * (drim[0]*fi[1] + drjm[0]*fj[1] + drkm[0]*fk[1]); - v[4] = 0.25 * (drim[0]*fi[2] + drjm[0]*fj[2] + drkm[0]*fk[2]); - v[5] = 0.25 * (drim[1]*fi[2] + drjm[1]*fj[2] + drkm[1]*fk[2]); + if (vflag_either) { + v[0] = (drim[0]*fi[0] + drjm[0]*fj[0] + drkm[0]*fk[0]); + v[1] = (drim[1]*fi[1] + drjm[1]*fj[1] + drkm[1]*fk[1]); + v[2] = (drim[2]*fi[2] + drjm[2]*fj[2] + drkm[2]*fk[2]); + v[3] = (drim[0]*fi[1] + drjm[0]*fj[1] + drkm[0]*fk[1]); + v[4] = (drim[0]*fi[2] + drjm[0]*fj[2] + drkm[0]*fk[2]); + v[5] = (drim[1]*fi[2] + drjm[1]*fj[2] + drkm[1]*fk[2]); - vatom[i][0] += v[0]; vatom[i][1] += v[1]; vatom[i][2] += v[2]; - vatom[i][3] += v[3]; vatom[i][4] += v[4]; vatom[i][5] += v[5]; - vatom[j][0] += v[0]; vatom[j][1] += v[1]; vatom[j][2] += v[2]; - vatom[j][3] += v[3]; vatom[j][4] += v[4]; vatom[j][5] += v[5]; - vatom[k][0] += v[0]; vatom[k][1] += v[1]; vatom[k][2] += v[2]; - vatom[k][3] += v[3]; vatom[k][4] += v[4]; vatom[k][5] += v[5]; - vatom[m][0] += v[0]; vatom[m][1] += v[1]; vatom[m][2] += v[2]; - vatom[m][3] += v[3]; vatom[m][4] += v[4]; vatom[m][5] += v[5]; + if (vflag_global) { + virial[0] += v[0]; + virial[1] += v[1]; + virial[2] += v[2]; + virial[3] += v[3]; + virial[4] += v[4]; + virial[5] += v[5]; + } + + if (vflag_atom) { + v[0] *= 0.25; + v[1] *= 0.25; + v[2] *= 0.25; + v[3] *= 0.25; + v[4] *= 0.25; + v[5] *= 0.25; + + vatom[i][0] += v[0]; vatom[i][1] += v[1]; vatom[i][2] += v[2]; + vatom[i][3] += v[3]; vatom[i][4] += v[4]; vatom[i][5] += v[5]; + vatom[j][0] += v[0]; vatom[j][1] += v[1]; vatom[j][2] += v[2]; + vatom[j][3] += v[3]; vatom[j][4] += v[4]; vatom[j][5] += v[5]; + vatom[k][0] += v[0]; vatom[k][1] += v[1]; vatom[k][2] += v[2]; + vatom[k][3] += v[3]; vatom[k][4] += v[4]; vatom[k][5] += v[5]; + vatom[m][0] += v[0]; vatom[m][1] += v[1]; vatom[m][2] += v[2]; + vatom[m][3] += v[3]; vatom[m][4] += v[4]; vatom[m][5] += v[5]; + } } } @@ -1534,8 +1552,7 @@ void Pair::v_tally2(int i, int j, double fpair, double *drij) called by AIREBO and Tersoff potential, newton_pair is always on ------------------------------------------------------------------------- */ -void Pair::v_tally3(int i, int j, int k, - double *fi, double *fj, double *drik, double *drjk) +void Pair::v_tally3(int i, int j, int k, double *fi, double *fj, double *drik, double *drjk) { double v[6];