Hi,

I am trying to modify based on pair coul/debye, its potential is

But maybe the following red r2inv should be rinv? It might be a mistake or I understand any place wrongly?

double PairCoulDebye::single(int i, int j, int itype, int jtype,

double rsq, double factor_coul, double factor_lj,

double &fforce)

{

double r2inv,r,rinv,forcecoul,phicoul,screening;

r2inv = 1.0/rsq;

r = sqrt(rsq);

rinv = 1.0/r;

screening = exp(-kappa*r);

forcecoul = force->qqrd2e * atom->q[i]*atom->q[j] *
screening * (kappa + rinv);
fforce = factor_coul*forcecoul * r2inv;

phicoul = force->qqrd2e * atom->q[i]*atom->q[j] * rinv * screening;
return factor_coul*phicoul;

}

Best Regards,

NIE PIN(Ms)

Nanyang Technology University, School of Physical and Mathematical Science

Tel: +65 84522947

email: [email protected]…2775…