Dear all,

I need some clarification in the c++ file of pair_style coul/wolf as I would like to modify this a bit and use a modified potential.

A part from the program is given below,

for (ii = 0; ii < inum; ii++) {

i = ilist[ii];

qtmp = q[i];

xtmp = x[i][0];

ytmp = x[i][1];

ztmp = x[i][2];

jlist = firstneigh[i];

jnum = numneigh[i];

qisq = qtmp*qtmp;
e_self = -(e_shift/2.0 + alf/MY_PIS) * qisq*qqrd2e;

if (evflag) ev_tally(i,i,nlocal,0,0.0,e_self,0.0,0.0,0.0,0.0);

for (jj = 0; jj < jnum; jj++) {

j = jlist[jj];

factor_coul = special_coul[sbmask(j)];

j &= NEIGHMASK;

delx = xtmp - x[j][0];

dely = ytmp - x[j][1];

delz = ztmp - x[j][2];

rsq = delx*delx + dely*dely + delz*delz;

if (rsq < cut_coulsq) {

r = sqrt(rsq);

prefactor = qqrd2e*qtmp*q[j]/r;

erfcc = erfc(alf*r);
erfcd = exp(-alf*alf

*r*r);

**v_sh = (erfcc - e_shift**

*r) * prefactor;*

**dvdrr = (erfcc/rsq + 2.0****alf/MY_PIS * erfcd/r) + f_shift;**

**forcecoul = dvdrr**

*rsq*prefactor;**if (factor_coul < 1.0) forcecoul -= (1.0-factor_coul)*prefactor;**

**fpair = forcecoul / rsq;**

**f[i][0] += delx*fpair;**

**f[i][1] += dely*fpair;**

**f[i][2] += delz*fpair;**

**if (newton_pair || j < nlocal) {**

**f[j][0] -= delx*fpair;**

**f[j][1] -= dely*fpair;**

**f[j][2] -= delz*fpair;**

**}**

**if (eflag) {**

**ecoul = v_sh;**

**if (factor_coul < 1.0) ecoul -= (1.0-factor_coul)*prefactor;**

**} else ecoul = 0.0;**

**if (evflag) ev_tally(i,j,nlocal,newton_pair,**

**0.0,ecoul,fpair,delx,dely,delz);**

}

}

}

My questions are

1)What is **v_sh** calculating? Is that the force term?

- And in the next step
**dvdrr,****erf****cc**is divided by**rsq**but the next step is divided by**r**.

3)What is the importance of **factor_coul** in this program?

Kindly clarify this. Thanks a lot in advance.