[lammps-users] How can I implement the following Potential for a Cascade Simulation

I am trying to use LAMMPS to simulate a displacement cascade in a uranium dioxide lattice. My difficulty lies in implementing a particular potential for LAMMPS to use. My potential consists of a ZBL term for short-range interactions and a coulombic + Buckingham term for intermediate-range interactions. I have splined these two regions together to create a smooth potential for each of the three interacting atom pairs (U-U, O-O, O-U). For each atom pair, I have an analytic expression for this potential.

For input into LAMMPS, I thought that the best way to implement this potential was to use pair_style hybrid/overlay to superpose (1) pair_style coul/long and (2) a table from my above splined potential from which I have subtracted the coulombic term. My motivation was to recreate my original spline and to preserve long-range coulombic interactions. For each atom pair, I input a separate table. I used 11,000 points for these input tables and I used the spline option for interpolation. One noteworthy issue is that the repulsive ZBL potential leads to a very high gradient for values in my input tables that are close to 0. Because of this, I used 10,000 points for this core ‘high gradient’ region.

I implemented this in LAMMPS and I used the pair_write command to output a table containing the interpolated values of the potential created by the hybrid/overlay described above. The graph of this LAMMPS-generated potential has huge osciallations in the region corresponding roughly to the ZBL core.

To trouble-shoot, I examined each part of the hybrid/overlay potential separately. I discovered 2 Problems.

First, I used pair_write to examine LAMMPS’ interpolations of coul/long. That is, I defined all interactions to be only via coul/long (without inputting my tables) and then I output the interpolated potentials using pair_write. I discovered that the potentials for all three atom pairs (U-U, O-O, O-U) were identical. The tables for each were the same despite my specifying different charges for uranium and oxygen. I see no other possible explanation for this.

Then, I defined all interactions to be strictly according to my input tables (without any coulombic interaction). All 3 tables output by pair_write contain huge oscillations near the core region and do not look anything like a splined interpolation. However, the output tables are reasonably close to my input tables at longer distances. I played with the number and spacing of my table points but this did not change anything. When I changed the interpolating option from spline to linear, it removed the oscillations at short ranges; however, instead, at short ranges, the potential does not decrease exponentially but rather increases in the approximate ZBL region before sharply decreasing to more realistic values at greater distances.

I would like to know

  1. why LAMMPS treats the Coulombic potential of all three atom pairs the same;
  2. how to submit a table that LAMMPS can interpolate correctly. Of course, if there is a better way to approach my problem, then I would be happy to hear this.

Thank you for your attention.

It should be possible to do what you want with pair hybrid/overlay, but
I haven't tried it before. So possibly there are problems in LAMMPS with this.
If you have an analytic expression, you may just want to code it up in
a pair style yourself. It's not that hard to do.

I'll let Paul answer your Q about tables and oscillations at the core. It is
hard for splines to fit regions of hi slope, but there are some tricks he
can tell you. You should be able to get LAMMPS to avoid splining
altogether and just do linear interpolation on your values. If you have
11000 of them, that should be sufficient.

First, I used pair_write to examine LAMMPS' interpolations of coul/long. That is, I defined all interactions to be only >via coul/long (without inputting my tables) and then I output the interpolated potentials using pair_write. I discovered >that the potentials for all three atom pairs (U-U, O-O, O-U) were identical. The tables for each were the same >despite my specifying different charges for uranium and oxygen. I see no other possible explanation for this.

Did you specify different charge values as part of the pair_write
command for each of the 3 cases?
LAMMPS has no way of knowing what charge goes with what atom type, as
charges are
defined on a per-atom basis, not per-type.

Steve

Dear Steve,

Thank you very much for your help. Regarding the coulombic potentials output by pair_write, specifying the atom charges as additional arguments to the pair_write command has corrected the problem. Since the cause of the problem was in my arguments for the pair_write command, I assume that the correct potential was actually computed originally, but I was just not able to output it.

I will wait for Paul’s answer on how to interpolate a table of a “high-gradient” function. If at all possible, I would like to avoid writing a new pair style just for this simulation.

Sincerely,
Milad

Milad,

You can get rid of the problem you’re seeing with interpolating a table of “high-gradient” function by doing the following:

  1. Use a large number of points (10,000 should be good) in your input table. Make sure your points are evenly spaced in r^2 space.
  2. Use the same number of points in your pair table command as are in your input table.
  3. Use linear interpolation.

Following this prescription should force LAMMPS to use only data that you have provided (plus linear interpolation between those points that you provided). The funny splines won’t show up anywhere. Test this by outputting using the pair_write command.

If you have further questions regarding pair_table, I recommend a careful re-reading of the documentation ( http://lammps.sandia.gov/doc/pair_table.html ) and of the previous lammps-users threads on this topic.

Paul