Hi All,
I’ve got a simulation running which contains H, C, O, Fe in the gas phase. My gas is unusual in that H is very abundant. Because I don’t know yet which interatomic potential to use for my system, I figured a good place to start would be the ReaxFF potentials that came with LAMMPS: ffield.reax.Fe_O_C_H based on Aryanpour, Van Duin, Kubicki, Phys Chem A 2010, 114, 6298.
After running for a few picoseconds at 300K, I find that the oxygen forms primarily H3O with the hydrogen, not H2O like I would expect. I did some literature searching and found modified ReaxFF potentials which should improve H3O, HO fidelity (Zhang, Van Duin, J Phys Chem B 2017, 121, 24). I transferred the terms for H, O, as well as the angles, bonds and torsions involving just H and O. However, the result appears to be the same.
Question: is my composition far enough out of the expected range for which this ReaxFF potential was calibrated, so that I am not justified in using it? If so, is there some interatomic potential I should be using for this type of system instead? I did some poking around on OpenKIM, and the NIST repository, but I didn’t see anything obvious that had all four elements and was for gas phase (except for LJ potentials which I doubt will suffice here). Do I need to look at making my own potentials for this system starting with DFT? Simply run longer? Shorter timestep (using 0.25 fs)?
I’m attaching the script in case there is something obvious to someone else that may be throwing my results off.
Cheers,
Zack
box.lmp (98.7 KB)
ffield.reax.Fe_O_C_H (12.9 KB)
in.ism (1.1 KB)
When I look at the abstracts of the two publications you mention, both are mentioning the liquid phase and Iron and Iron oxide and thus they are quite far away from the - rather unusual - system that you are describing. So I doubt that either is a suitable potential for your use case. ReaxFF potentials can handle a variety of compounds but they also are not very transferable, i.e. you not just need the elements, but also the same conditions and types of compound.
I doubt that there are any ready to use force fields unless you go to very “generic” (and thus much less accurate) ones. Learning how to parameterize force fields is a non-trivial operation. Most certainly not something that is suitable for a beginner in the field. My suggestion would be to look for a method that is more flexible, which usually comes with a higher computational codes like some semi-empirical quantum software or even some (pseudo potential) DFT tools. At the very least, it should give some some initial insight into the behavior of small systems under the given conditions and you may use that later to review and improve any more empirical calculations.
Also, gas phase with H, C, O and Fe at 300K sounds very unusual. One would expect rather higher temperatures… or else everything would condense.
Hi Axel,
Thank you for the hints, notably about how 300K is low temp. I therefore tried spending more time equilibrating at higher temp (1000K) before bringing the temperature back down and the ratio of H2O/H3O improved, which led me to suspect it may be that my starting conditions were sketchy: namely I start with H atoms, not H molecules, so an H atom colliding with an H2O will produce an H3O.
I made a new simulation where I fill the box with H, ran it until the H atoms reacted into H2 gas and then added the other atoms. After a second equilibration it behaved much better. Most of the O is now bound into water. It seems just finding the right initial conditions already made a huge difference.
It is good to know about transferability for ReaxFF. After sending out e-mails, I have been made aware of some ReaxFF potentials used for catalysts and so simulate the gasses next to the metals. I’m hopeful these are enough to bring my quality over the threshold of usefulness. My fingers are crossed.
Cheers,
Zack