ICODE-ADC/lammps
4
0
Fork 0
Code Issues Pull Requests Actions Packages Projects Releases Wiki Activity

enforce only newton pair on with smatb pair styles. add unit tests

Browse Source
This commit is contained in:
Axel Kohlmeyer
2022-04-27 18:08:05 -04:00
parent 2d45e3340f
commit b76594e551
9 changed files with 319 additions and 55 deletions
Download Patch File Download Diff File Expand all files Collapse all files

3
src/MANYBODY/pair_eam.cpp
View File

@ -319,8 +319,7 @@ void PairEAM::compute(int eflag, int vflag)
}
if (eflag) evdwl = scale[itype][jtype]*phi;
if (evflag) ev_tally(i,j,nlocal,newton_pair,
evdwl,0.0,fpair,delx,dely,delz);
if (evflag) ev_tally(i,j,nlocal,newton_pair,evdwl,0.0,fpair,delx,dely,delz);
}
}
}

42
src/SMTBQ/pair_smatb.cpp
View File

@ -110,11 +110,7 @@ void PairSMATB::compute(int eflag, int vflag)
int newton_pair = force->newton_pair;
// zero out on_eb
if (newton_pair) {
memset(on_eb, 0, nall * sizeof(on_eb[0]));
} else {
memset(on_eb, 0, nlocal * sizeof(on_eb[0]));
}
memset(on_eb, 0, nall * sizeof(double));
int inum = list->inum;
int *ilist = list->ilist;
@ -156,14 +152,14 @@ void PairSMATB::compute(int eflag, int vflag)
qsiexpq = qsiexpq * qsiexpq;
}
on_eb[i] += qsiexpq;
on_eb[j] += qsiexpq;
if (newton_pair) on_eb[j] += qsiexpq;
}
}
}
// communicate the squared bonding energy between the various bins
comm->reverse_comm(this);
if (newton_pair) comm->reverse_comm(this);
// Support Loop: take the square root of the bonding energy and
// accumulate it in the energy accumulator if needed the store the
@ -178,7 +174,7 @@ void PairSMATB::compute(int eflag, int vflag)
} else {
on_eb[i] = 0.0;
}
//if needed the bonding energy is accumulated:
// if needed the bonding energy is accumulated:
if (eflag_either) {
if (eflag_atom) { eatom[i] -= eb_i; }
if (eflag_global) { eng_vdwl -= eb_i; }
@ -233,20 +229,7 @@ void PairSMATB::compute(int eflag, int vflag)
3.0 * x3[itype][jtype] * polyval2)) *
qsiexpq;
}
// if needed the repulsive energy is accumulated:
if (eflag_either) {
if (eflag_atom) {
eatom[i] += aexpp;
if (newton_pair || j < nlocal) { eatom[j] += aexpp; }
}
if (eflag_global) {
if (newton_pair || j < nlocal) {
eng_vdwl += 2.0 * (aexpp);
} else {
eng_vdwl += aexpp;
}
}
}
// calculates the module of the pair energy between i and j
fpair = (Fb * (on_eb[i] + on_eb[j]) + Fr) / dij;
@ -258,10 +241,8 @@ void PairSMATB::compute(int eflag, int vflag)
f[j][1] -= del[1] * fpair;
f[j][2] -= del[2] * fpair;
}
if (vflag_atom) {
ev_tally(i, j, nlocal, newton_pair, 0.0,
0.0, //Energy is tally'd in the other parts of the potential
fpair, del[0], del[1], del[2]);
if (evflag) {
ev_tally(i, j, nlocal, newton_pair, 2.0 * aexpp, 0.0, fpair, del[0], del[1], del[2]);
}
}
}
@ -350,6 +331,15 @@ void PairSMATB::coeff(int narg, char **arg)
if (count == 0) error->all(FLERR, "Incorrect args for pair coefficients");
}
/* ------------------------------------------------------------------------ */
void PairSMATB::init_style()
{
if (force->newton_pair == 0) error->all(FLERR, "Pair style smatb requires newton pair on");
neighbor->add_request(this);
}
/* ----------------------------------------------------------------------
init for one type pair i,j and corresponding j,i
------------------------------------------------------------------------- */

1
src/SMTBQ/pair_smatb.h
View File

@ -35,6 +35,7 @@ class PairSMATB : public Pair {
void compute(int, int) override;
void settings(int, char **) override;
void coeff(int, char **) override;
void init_style() override;
double init_one(int, int) override;
void write_restart(FILE *) override;
void read_restart(FILE *) override;

44
src/SMTBQ/pair_smatb_single.cpp
View File

@ -94,11 +94,7 @@ void PairSMATBSingle::compute(int eflag, int vflag)
int newton_pair = force->newton_pair;
// zero out on_eb
if (newton_pair) {
memset(on_eb, 0, nall * sizeof(on_eb[0]));
} else {
memset(on_eb, 0, nlocal * sizeof(on_eb[0]));
}
memset(on_eb, 0, nall * sizeof(double));
int inum = list->inum;
int *ilist = list->ilist;
@ -136,14 +132,14 @@ void PairSMATBSingle::compute(int eflag, int vflag)
qsiexpq = qsiexpq * qsiexpq;
}
on_eb[i] += qsiexpq;
on_eb[j] += qsiexpq;
if (newton_pair) on_eb[j] += qsiexpq;
}
}
}
// communicate the squared bonding energy between the various bins
comm->reverse_comm(this);
if (newton_pair) comm->reverse_comm(this);
// Support Loop: take the square root of the bonding energy and
// accumulate it in the energy accumulator if needed the store the
@ -158,7 +154,7 @@ void PairSMATBSingle::compute(int eflag, int vflag)
} else {
on_eb[i] = 0.0;
}
//if needed the bonding energy is accumulated:
// if needed the bonding energy is accumulated:
if (eflag_either) {
if (eflag_atom) { eatom[i] -= eb_i; }
if (eflag_global) { eng_vdwl -= eb_i; }
@ -205,20 +201,7 @@ void PairSMATBSingle::compute(int eflag, int vflag)
qsiexpq = x5 * polyval5 + x4 * polyval4 + x3 * polyval3;
Fb = ((5.0 * x5 * polyval4 + 4.0 * x4 * polyval3 + 3.0 * x3 * polyval2)) * qsiexpq;
}
// if needed the repulsive energy is accumulated:
if (eflag_either) {
if (eflag_atom) {
eatom[i] += aexpp;
if (newton_pair || j < nlocal) { eatom[j] += aexpp; }
}
if (eflag_global) {
if (newton_pair || j < nlocal) {
eng_vdwl += 2.0 * (aexpp);
} else {
eng_vdwl += aexpp;
}
}
}
// calculates the module of the pair energy between i and j
fpair = (Fb * (on_eb[i] + on_eb[j]) + Fr) / dij;
@ -230,10 +213,8 @@ void PairSMATBSingle::compute(int eflag, int vflag)
f[j][1] -= del[1] * fpair;
f[j][2] -= del[2] * fpair;
}
if (vflag_atom) {
ev_tally(i, j, nlocal, newton_pair, 0.0,
0.0, //Energy is tally'd in the other parts of the potential
fpair, del[0], del[1], del[2]);
if (evflag) {
ev_tally(i, j, nlocal, newton_pair, 2.0 * aexpp, 0.0, fpair, del[0], del[1], del[2]);
}
}
}
@ -247,7 +228,7 @@ void PairSMATBSingle::compute(int eflag, int vflag)
void PairSMATBSingle::settings(int narg, char **)
{
if (narg > 0) error->all(FLERR, "Illegal pair_style command: smatb accepts no options");
if (narg > 0) error->all(FLERR, "Illegal pair_style command: smatb/single accepts no options");
}
/* ----------------------------------------------------------------------
@ -302,6 +283,15 @@ void PairSMATBSingle::coeff(int narg, char **arg)
if (count == 0) error->all(FLERR, "Incorrect args for pair coefficients");
}
/* ------------------------------------------------------------------------ */
void PairSMATBSingle::init_style()
{
if (force->newton_pair == 0) error->all(FLERR, "Pair style smatb/single requires newton pair on");
neighbor->add_request(this);
}
/* ----------------------------------------------------------------------
init for one type pair i,j and corresponding j,i
------------------------------------------------------------------------- */

1
src/SMTBQ/pair_smatb_single.h
View File

@ -35,6 +35,7 @@ class PairSMATBSingle : public Pair {
void compute(int, int) override;
void settings(int, char **) override;
void coeff(int, char **) override;
void init_style() override;
double init_one(int, int) override;
void write_restart(FILE *) override;
void read_restart(FILE *) override;