corrections to documentation and implementation of pair styles lj/mdf, buck/mdf, and lennard/mdf

[Examples:]

pair_style lj/mdf 2.5 3.0

pair_coeff * * 1 1

pair_coeff 1 1 1 1.1 2.8 3.0 3.2 :pre

pair_coeff * * 1.0 1.0

pair_coeff 1 1 1.1 2.8 3.0 3.2 :pre

pair_style buck 2.5 3.0

pair_coeff * * 100.0 1.5 200.0

pair_coeff * * 100.0 1.5 200.0 3.0 3.5 :pre

pair_style lennard/mdf 2.5 3.0

pair_coeff * * 1 1

pair_coeff 1 1 1 1.1 2.8 3.0 3.2 :pre

pair_coeff * * 1.0 1.0

pair_coeff 1 1 1021760.3664 2120.317338 3.0 3.2 :pre

[Description:]

:c,image(Eqs/pair_mdf-4.jpg)

The following coefficients must be defined for each pair of atoms

types via the pair_coeff command as in the examples above, or in the

data file or restart files read by the "read_data"_read_data.html or

"read_restart commands"_read_restart.html, or by mixing as described

below:

Either the first two or all of the following coefficients must be

defined for each pair of atoms types via the pair_coeff command as

in the examples above, or in the data file read by the

"read_data"_read_data.html. The two cutoffs default to the global

values and epsilon and sigma can also be determined by mixing as

described below:

epsilon (energy units)

sigma (distance units)

:line

For the {buck/mdf} pair_style, the potential energy, {E(r)}, is the

standard Buckingham potential:

standard Buckingham potential with three required coefficients.

The two cutoffs can be omitted and default to the corresponding

global values:

:c,image(Eqs/pair_mdf-5.jpg)

\rho (distance units)

C (energy-distance^6 units)

r_m (distance units)

r_{cut}$ (distance units) :ul

r_{cut} (distance units) :ul

:line

For the {lennard/mdf} pair_style, the potential energy, {E(r)}, is the

standard 12-6 Lennard-Jones written in the $A/B$ form:

standard 12-6 Lennard-Jones written in the A/B form:

:c,image(Eqs/pair_mdf-6.jpg)

The following coefficients must be defined for each pair of atoms

types via the pair_coeff command as in the examples above, or in the

data file or restart files read by the read_data or read_restart

commands, or by mixing as described below:

data file read by the read_data commands, or by mixing as described below.

The two cutoffs default to their global values and must be either both

given or both left out:

A (energy-distance^12 units)

B (energy-distance^6 units)

[Mixing, shift, table, tail correction, restart, rRESPA info]:

For atom type pairs I,J and I != J, the epsilon and sigma coefficients

and cutoff distance for all of the lj/cut pair styles can be mixed.

and cutoff distances for the lj/mdf pair style can be mixed.

The default mix value is {geometric}.  See the "pair_modify" command

for details.

All of the {lj/cut} pair styles support the

"pair_modify"_pair_modify.html shift option for the energy of the

Lennard-Jones portion of the pair interaction.

The {lj/cut/coul/long} and {lj/cut/tip4p/long} pair styles support the

"pair_modify"_pair_modify.html table option since they can tabulate

the short-range portion of the long-range Coulombic interaction.

All of the {lj/cut} pair styles support the

"pair_modify"_pair_modify.html tail option for adding a long-range

tail correction to the energy and pressure for the Lennard-Jones

portion of the pair interaction.

All of the {lj/cut} pair styles write their information to "binary

restart files"_restart.html, so pair_style and pair_coeff commands do

not need to be specified in an input script that reads a restart file.

The {lj/cut} and {lj/cut/coul/long} pair styles support the use of the

{inner}, {middle}, and {outer} keywords of the "run_style

respa"_run_style.html command, meaning the pairwise forces can be

partitioned by distance at different levels of the rRESPA hierarchy.

The other styles only support the {pair} keyword of run_style respa.

See the "run_style"_run_style.html command for details.

for details. The other two pair styles buck/mdf and lennard/mdf do not

support mixing, so all I,J pairs of coefficients must be specified

explicitly.

None of the lj/mdf, buck/mdf, or lennard/mdf pair styles supports

the "pair_modify"_pair_modify.html shift option or long-range

tail corrections to pressure and energy.

These styles write their information to "binary restart

files"_restart.html, so pair_style and pair_coeff commands do not need

to be specified in an input script that reads a restart file.

These styles can only be used via the {pair} keyword of the "run_style

respa"_run_style.html command.  They do not support the {inner},

{middle}, {outer} keywords.

:line

    memory->destroy(rhoinv);

    memory->destroy(buck1);

    memory->destroy(buck2);

    memory->destroy(offset);

  }

}

  memory->create(rhoinv,n+1,n+1,"pair:rhoinv");

  memory->create(buck1,n+1,n+1,"pair:buck1");

  memory->create(buck2,n+1,n+1,"pair:buck2");

  memory->create(offset,n+1,n+1,"pair:offset");

}

/* ----------------------------------------------------------------------

void PairBuckMDF::coeff(int narg, char **arg)

{

  if (narg != 5 && narg != 7) error->all(FLERR,"Incorrect args for pair coefficients");

  if (narg != 5 && narg != 7)

    error->all(FLERR,"Incorrect args for pair coefficients");

  if (!allocated) allocate();

  int ilo,ihi,jlo,jhi;

  buck1[i][j] = a[i][j]/rho[i][j];

  buck2[i][j] = 6.0*c[i][j];

  if (offset_flag && (cut[i][j] > 0.0)) {

    double rexp = exp(-cut[i][j]/rho[i][j]);

    offset[i][j] = a[i][j]*rexp - c[i][j]/pow(cut[i][j],6.0);

  } else offset[i][j] = 0.0;

  cut_inner[j][i] = cut_inner[i][j];

  cut_inner_sq[i][j] = cut_inner[i][j]*cut_inner[i][j];

  cut_inner_sq[j][i] = cut_inner_sq[i][j];

  rhoinv[j][i] = rhoinv[i][j];

  buck1[j][i] = buck1[i][j];

  buck2[j][i] = buck2[i][j];

  offset[j][i] = offset[i][j];

  // compute I,J contribution to long-range tail correction

  // count total # of atoms of type I and J via Allreduce

  if (tail_flag) {

    int *type = atom->type;

    int nlocal = atom->nlocal;

    double count[2],all[2];

    count[0] = count[1] = 0.0;

    for (int k = 0; k < nlocal; k++) {

      if (type[k] == i) count[0] += 1.0;

      if (type[k] == j) count[1] += 1.0;

    }

    MPI_Allreduce(count,all,2,MPI_DOUBLE,MPI_SUM,world);

    double rho1 = rho[i][j];

    double rho2 = rho1*rho1;

    double rho3 = rho2*rho1;

    double rc = cut[i][j];

    double rc2 = rc*rc;

    double rc3 = rc2*rc;

    etail_ij = 2.0*MY_PI*all[0]*all[1]*

      (a[i][j]*exp(-rc/rho1)*rho1*(rc2 + 2.0*rho1*rc + 2.0*rho2) -

       c[i][j]/(3.0*rc3));

    ptail_ij = (-1/3.0)*2.0*MY_PI*all[0]*all[1]*

      (-a[i][j]*exp(-rc/rho1)*

       (rc3 + 3.0*rho1*rc2 + 6.0*rho2*rc + 6.0*rho3) + 2.0*c[i][j]/rc3);

  }

  return cut[i][j];

}

double PairLJ_AB_MDF::init_one(int i, int j)

{

  if (setflag[i][j] == 0) {

    aparm[i][j] = mix_energy(aparm[i][i],aparm[j][j],

                               bparm[i][i],bparm[j][j]);

    bparm[i][j] = mix_distance(bparm[i][i],bparm[j][j]);

    cut_inner[i][j] = mix_distance(cut_inner[i][i],cut_inner[j][j]);

    cut[i][j] = mix_distance(cut[i][i],cut[j][j]);

  }

  if (setflag[i][j] == 0) error->all(FLERR,"All pair coeffs are not set");

  cut_inner_sq[i][j] = cut_inner[i][j]*cut_inner[i][j];

/* ---------------------------------------------------------------------- */

PairLJMDF::PairLJMDF(LAMMPS *lmp) : Pair(lmp) {}

PairLJMDF::PairLJMDF(LAMMPS *lmp) : Pair(lmp) {

  writedata = 1;

}

/* ---------------------------------------------------------------------- */