[lammps-users] Observation of C2 molecules resulting from LAMMPS/ReaxFF simulations

Dear LAMMPS users:

I am writing to inquire about some issues we faced in our MD simulations using LAMMPS/Reax. To verify proper implementation of ReaxFF in LAMMPS, we have simulated methane/O2 and propene/O2 molecular systems. These systems are chosen the same as those reported in the paper published by Chenoweth and et al. [1]. In their paper, they had used serial version of ReaxFF developed by Dr. van Duin. Reaction pathways obtained from our simulations showed some outstanding discrepancies compared to corresponding ones found in the mentioned paper. We generated 10 unique propene/O2 configurations in order to conduct 10 independent simulations. In contrast to Ref. [1], C2 species was observed from all unique propene/O2 configurations. It should be noted here that the production of C2 had already been observed from our simulation of phenol monomer decomposition at two temperatures of 3000 K and 3400 K. However, this species is not reported in the experimental works as a decomposition product of phenol monomer.
Nielson and et al. [2] have discussed about inclusion of a new partial energy contribution, E(C2), in the ReaxFF definition of total energy of the system. It is done to destabilize C2 molecule. They argue that during the re-parameterization of the hydrocarbon ReaxFF description (ReaxFF-CH) with the QC data for all-carbon clusters, they found that the energy terms (excluding (EC2)) were sufficient to capture the relative energy and geometry of most of these clusters, the only exception being the C2 molecule. ReaxFF-CH erroneously predicts that the two carbons form a very strong (triple) bond, while in fact the triple bond would get de-stabilized by terminal radical electrons, and for that reason the carbon-carbon bond is not any stronger than a double bond.
I have asked about this issue for Dr. van Duin, developer of ReaxFF, and here is his response:
“I believe that the LAMMPS code does not recognize the C2-correction
that my standalone code employs. Because of that, the C2 molecule
becomes excessively stable.

Here are my questions:

Has anyone observed the same thing in MD simulations combustion/decomposition of hydrocarbons as we did? If not, what could be the reason for our observation of C2? Could it be because of the lack of EC2 in the Reax forcefield?

Thanks a lot for all your help in advance.

Best regards.

References
1. Chenoweth, K, van Duin, A. C. T. and Goddard, W. A., “A ReaxFF Reactive Force Field for Molecular Dynamics Simulations of Hydrocarbon Oxidation,” Journal of Physical Chemistry A, 112, 2008, pp.1040-1053.
2. Nielson, K.D., van Duin, A. C. T., Oxgaard, J., Deng, W.Q., and Goddard, W. A., “Development of the ReaxFF Reactive Force Field for Describing Transition Metal Catalyzed Reactions, with Application to the Initial Stages of the Catalytic Formation of Carbon Nanotubes,” Journal of Physical Chemistry A, 109, 2005, pp. 493-499.

Farshad Saberi

This is a question for Aidan.

Steve

After consulting with the ReaxFF authors, they found the C2 problem was due
to it being zeroed out in ffield.reax.cho. See diff below. I just checked
this in, and it should be posted on the website shortly.

--- ffield.reax.cho (revision 5213)
+++ ffield.reax.cho (working copy)
@@ -5,7 +5,7 @@
    26.5405 !p(coa2)
     1.5105 !p(trip4)
     6.6630 !p(trip3)
- 0.0000 !kc2
+ 70.0000 !kc2
     1.0588 !p(ovun6)
     4.6000 !p(trip2)
    12.1176 !p(ovun7)

The ReaxFF authors have found an additional change that is required for the ffield.reax.cho version of ReaxFF, related to triple bond stabilization. Normally this is applied to all triple bonds, but for ffield.reax.cho version, it should only be applied to CO triple bonds i.e. in carbon monoxide. I have added an argument to pair_style reax to toggle this behavior.:

reax args: hbcut hbnewflag tripflag precision

pair_style reax 10.0 0 1 1.0e-6

It is checked in now and will appear in the next patch. You will know you have the new version when your old script breaks. See doc page for more info. Don’t forget to recompile lib/reax. Please let me know how it works for you.

Aidan

Dear Aidan,

I really appreciate all of your efforts and help. I am trying it and will let you know.

Thanks

Farshad Saberi

Dear Aidan,

I'd like to appreciate again for your help! After simulating several hydrocarbon systems with the modified ReaxFF, I no longer observed the
production of C2 species.

Thanks.

Farshad