[lammps-users] Help on the output file of 'fix reax/bonds' command

LAMMPS LAMMPS Mailing List Mirror
1

Dear Lammps-er,

I am trying to write some code to analyze the output from command ‘fix reax/bonds’ and want to make sure of the two things:

(1) When Lammps is writing the output from ‘fix reax/bonds’ command, it does not write each line in the order of increasing atom ID, but instead scans the regions and writes the output from +X to -X, +Y to -Y, +Z to -Z. Is that correct?

(2) Each output line in the body part contains the information in the order of:

id type nb id_1…id_nb mol bo_1…bo_nb abo nlp q

Could you please show me how to determine the number of “bo_nb” there? Also, I do not know what ‘mol’ stands for.

Part of my output from ‘fix reax/bonds’ is in the section below.

Thanks a lot in advance.

Best,
Paul

NCSU, Raleigh

2

This is a question for Aidan. But I presume there is no
particular order to the atoms written to the file. Each
line has the atom ID, so that should be all you need.

Steve

3

The list of atoms is not sorted and the order can be different on different
timesteps. Instead, they are identified by id. I attach a c code that
analyzes the LAMMPS reax/bonds file. If you find it useful (and tell me), I
will add it to the LAMMPS tools directory.

The statement " The format of the output file should be self-explantory. "
is our way of saying the format is not documented. If you want to know what
this stuff means, you have to understand the arrays ia and iag:

lib/reax/cbkia.blk: $/cbkia/ ia(nat,mbond+3),iag(nat,mbond+3)

which are populated in lib/reax/reax_connect.F.

Aidan

mol_fra.c (24.3 KB)

4

Dear Aidan and Steve,

Thanks for the explanations and the code. I have not tried it out yet, but I will let you know whether it is working.

I am also attaching my fortran code that can be used to collect the bond connectivity information from the output of ‘fix reax/bonds’. It is not a general code for any system, but people can easily modify and take advantage of it, if they want to.

Best wishes
Paul

bondConnectCheck_reaxLammps.f90 (7.63 KB)

5

Thanks! I put mol_fra.c and bondConnectCheck.f90 in a new directory
tools/reax. However, when I tried bondConnectCheck.f90 on attached file
(from examples/reax), it gave an error. It would be nice to have something
that works in a simple case. Could you try it yourself and see what is
wrong?

Thanks,

Aidan

bonds.reax.tatb (71 KB)

6

I found the problem in bondConnectCheck.f90. I see now that it reads the
bond file and echoes it to output, but does not do any analysis. I think you
will find mol_fra.c more useful.

7

Hi Aidan,

Sorry that I did not put it clear in my previous email, but as I mentioned, my code is not a general one, and I use it to subtract the bond connectivity of NH3 molecules only.

That is why in the beginning I/O section, I have the followings:

open (unit=20, file=‘N129.txt’, status=‘unknown’)
open (unit=21, file=‘N133.txt’, status=‘unknown’)
open (unit=22, file=‘N137.txt’, status=‘unknown’)
open (unit=23, file=‘N141.txt’, status=‘unknown’)
open (unit=24, file=‘N145.txt’, status=‘unknown’)
open (unit=25, file=‘N149.txt’, status=‘unknown’)
open (unit=26, file=‘N153.txt’, status=‘unknown’)
open (unit=27, file=‘N157.txt’, status=‘unknown’)

The numbers after the letter “N” in the file names here are the index of N atoms in my initial data file (which is used in the command of “read-data” command of input script).

The basic idea of my code is to scan the output from ‘fix reax/bonds’ line by line and subtract all the connectivity information for all nitrogen atom. If the number of bonds of nitrogen atom is not equal to 3, in my case, it implies the happening of reaction at that step. I then write the image index (which is corresponding to the simulation step or time, if you read that number together with the ‘thermo’ frequency and ‘run’ iterations).

I checked your file (“bonds.reax.tatb”) and it looks good to me. So, for example, if you want to extract the connection for the atom with id 100, you need to change the constant definition in the code,

image = 1
headline = 7
tailline = 1
natom = 384

You totally have two images in the file, the first is at t=0, and the second is from one iteration, that is why “image = 1” in the code. Then you may want to change the I/O section, and corresponding “IF” commands there. I was using something like "if (id .eq. 129) then ", but you need change all those IFs (there are several IFs and you have to change all) into “if (id .eq. 100) then” for your case.

I got the code running correctly for your “bonds.reax.tatb” file after those changes. You may also need to delete the last blank line there in your file. The code reads the file as repeated combinations of “Head( 7 lines)+Body(384 lines ,in your case)+tail(1 line)”. There is no blank line in between.

I may update that code later so that we can use it for general purpose.

Best wishes,
Paul

8

Yeah, I realize that ^-^.

That fortran code only meets with my personal interest and needs a lot of improvement for general purpose. Your c code is much more useful in that sense.

Paul