[lammps-users] Reax-C and Reax-Fortran implementations

Hello,
I have been testing the various energy contributions computed by the C and Fortran versions of reax that are implemented within LAMMPS (REAX and USER-REAXC packages). I have two observations to share:

  1. I see disagreement between the two versions in the VdW energy. It seems that the C version includes the VdW inner wall term, while the Fortran version does not. If we are following Chenoweth paper that is cited in the LAMMPS documentation, then the VdW inner wall term should not be included in the equations (i.e., the Fortran version is consistent, the C version is not). If we are following another reference for this term, could you point me to the correct reference?

Chenoweth, van Duin and Goddard, Journal of Physical Chemistry A, 112, 1040-1053 (2008)

  1. With regards to the hydrogen bond energy, it seems that we have again deviated from the Chenoweth paper where the hydrogen bond term uses a sin4 term, rather than a sin8 term used in Chenoweth. The C version uses the sin4 term. The Fortran has both the sin4 and sin**8 terms coded, where the functional form is controlled by the lhbnew flag in reax_poten.F as written below:

2227 if (lhbnew .eq. 0) then
2228 ehbh=(1.0-exphu1)dehb(ityhb)exphu2sin2sin2sin2sin2
2229 else
2230 ehbh=(1.0-exphu1)dehb(ityhb)exphu2sin2sin2
2231 endif

The lhbnew flag is set within the lammps/src/REAX/pair_reax.cpp file (ihbnew) and set to 1, indicating that the sin**4 functional form should be used. Here is the code from pair_reax.cpp:

48 cutmax = 0.0;
49 hbcut = 6.0;
50 ihbnew = 1;
51 iprune = 4;
52 ihb = 1;
53 chpot = 0;

However, when passing the variable back to the reax library, somehow it is reset back to lhbnew=0, thereby defaulting to the sin**8 functional form. Perhaps I have missed where this is set back to 0 in the LAMMPS source, but it seems that this is not the desired outcome. Could you please provide comments as to whether you observe similar behavior?

Thanks,
–Jim

These are questions for Aidan ...

Steve

Jim,

Thanks for these observations. I am reminded of a joke that Irish people
like to say is told by English people about Irish people.

Victorian English gentleman: Stationmaster, why is the clock on the
eastbound platform ten minutes behind the clock on the rightbound platform?

Irish Stationmaster: Sure, if they both gave the same time, we wouldn't be
needing two clocks.

The stationmaster would understand the value of having two ReaxFF
implementations.

1. The modified inner-wall is used with some versions of ReaxFF, notably the
one Sergey Zybin used in his study of PETN. It is not coded in pair_style
reax. It is in pair_style reax/c, but I think it reverts to the standard
inner wall if not specifically requested in the ffield.reax file.

2. pair_style reax defaults to ihbnew=1, and this is correctly documented in
pair_reax.html.

Perhaps I have missed where this is set back to 0 in the LAMMPS source, but

The file example/reax/in.reax.tatb explicitly sets ihbnew=0 in the
pair_Style command:

pair_style reax 10.0 0 1.0e-6

pair_style reax/c only supports sin**4, which corresponds to ihbnew = 1.

Perhaps this default should be reverted to ihbnew=0, to match Chenoweth et
al.? Also, the differences between Chenoweth et al. and pair_style reax/c
should be documented?

I am hoping to soon have a more extensive set of examples that will help
guide people through some of these variations.

Aidan