[lammps-users] Doubts on .trj file output by reax/c

Dear lammps users,
I have some doubts about .trj file output after running a simulation with reax/c pair style. First, the timestep in .trj file is shown as 250fs, whereas I have specified it as 500fs. Second, it shows the ensemble as NVE, but I have used NVT. Third, it shows all box dimensions and temperature as zero. I am using this .trj file and fix reax/c/bonds to get the bond information for post-processing. Is the information from .trj file correct? Also is it the best method to extract bond information from reax/c simulation?

Regards
Abhiram

pe.trj (738 KB)

you should use a regular dump file. the trajectory functionality in the USER-REAXC code is unmaintained and untested and will be removed soon since there are known bugs.

Thanks for your reply. For reax/c pair style, I am using atom style “charge”. So we cannot use “compute bond/local”. For using a regular dump file to get bond information, did you mean to dump per-atom information from “fix reax/c/bond”?

Regards
Abhiram

Thanks for your reply. For reax/c pair style, I am using atom style “charge”. So we cannot use “compute bond/local”.

the problem is not so much the atom style (you can run reax/c calculations also with atom style full) but the fact that in reaxff simulations there are no explicit bonds defined in the input (which is what compute bond/local requires).

fix reax/c/bonds can output its information to a file directly. please see its documentation.

axel.

Yes reax/c/bond gives me a bond.reaxc file as output. But my doubt here is that in order to get bond information, I should write a separate script or Matlab code combining bond.reaxc and a dump file containing atomic coordinates? Because the output file from fix reax/c/bond is not in a readable format by post-processing software like VMD or OVITO. Am I right here?

Regards
Abhiram

i don’t know about OVITO, but in VMD you cannot have dynamic bond information from a dump file anyway. bond information is stored only globally from reading a topology file or from using a distance/covalent radius based heuristic.

VMD is assuming that the bond topology is fixed like it is required in CHARMM force field simulations.
So you will have to do custom tcl scripting to update the bond visualization from and external source anyway, or use the simplified atom distance based dynamic bond representation.

For how OVITO deals with dynamic bonds, please see in the OVITO documentation or their support forum/mailing list.

Just for clarity, “bonds” in ReaxFF are not binary on/off, but you have a bond order parameter that can change dynamically from 0 to N. so the fix reax/c/bonds output uses bond order ranges/cutoffs to detect bonds.

axel.