about the fix reax/c/species command in lammps

Dear Lammps Community,

Please forgive me to again report a strange phenomenon resulted from fix reax/c/species command.

I built an bicrystal alloy model consisting of total 3010 atoms.

By the atom type setting command: set group alloy type/fraction 1 0.82 12963 type/fraction 2 0.13 23486 type/fraction 3 0.05 23486,

I intend to set 82% Fe atom, 13% Cr atom and 5% Al atom.

And the log.lammps shows the following reasonable results:

2454 settings made for type/fraction

397 settings made for type/fraction

143 settings made for type/fraction

However, the fix reax/c/species command output a different result: Fe2118Cr251Al624.

Although I enabled the reax/c with kokkos package by the following changes, the same problem occurred again.

e.g. pair_style reax/c into pair_style reax/c/kk & fix reax/c/species into fix reax/c/species/kk

My lammps information is as follows:

Lammps version: 31Mar17;

Packages: Kokkos, User-reaxc, User-misc, molecule, qeq, reax, replica, rigid, kspace, manybody

Makefile.ubuntu: KOKKOS_DEVICES=OpenMP;

Running command: mpirun -np 4 -k on -sf kk -in in.mini

How do you like such problem?

Any advices will be highly appreciated! I look forward to your positive response,

Thanks in advance,

Best regards,

Xiaolong Liu from Korea Aerospace University.

please produce and post a simple and fast to run input example that quickly reproduces the issue, so that somebody else can reproduce and look into it.
by preference report this in our github issue tracker at: https://github.com/lammps/lammps/issues
examples for how to submit this kind of issue are #523, #500, #472

thanks,
axel.

Thanks Professor Axel,

I posted the following information (input script, log.lammps, output by fix reaxc/species)
on the GitHub as you advised,

Hope for the solutions,

Best regards,

Xiaolong

Part 1 The Input Script

# 1-Initialization

units real

atom_style charge

atom_modify id yes map hash first alloy

boundary p p p

newton on

dimension 3

# 2-Creation of Fe82Cr13Al5 bicrystal

lattice bcc 2.866
region box block 0 15 -10 10 0 5 side in
create_box 3 box

region left block 0 15 -10 0 0 5 side in
lattice bcc 2.866 orient x 5 -9 0 orient y 9 5 0 orient z 0 0 1
create_atoms 3 region left

region right block 0 15 0 10 0 5 side in
lattice bcc 2.866 orient x 5 9 0 orient y -9 5 0 orient z 0 0 1
create_atoms 3 region right

group alloy region box
group grain1 region left
group grain2 region right

set group alloy type/fraction 1 0.82 12963 type/fraction 2 0.13 23486 type/fraction 3 0.05 23486 # Fe75Cr20Al5

mass 1 55.85000
mass 2 52.00000
mass 3 26.98000

change_box all x scale 1.2 y scale 1 z scale 1

# 3-Force field

pair_style reax/c lmp_control safezone 4.0

pair_coeff * * ffield.reax Fe Cr Al

fix QEQ all qeq/reax 1 0.0 10.0 1.0e-6 reax/c

fix_modify QEQ energy yes

delete_atoms overlap 0.35 grain1 grain2 compress yes

fix BR all box/relax iso 0.0

neighbor 2.5 bin

neigh_modify every 1 delay 0 check yes once no cluster no exclude none page 100000 one 2000 binsize 0.0

# 4-Compute & thermo setting

compute total_potential_energy all pe

compute total_kinetic_energy all ke

compute reax all pair reax/c

thermo 250

thermo_style custom step time atoms temp press pe ke etotal epair

thermo_modify lost error lost/bond error norm yes flush yes line multi format float %10.4g temp thermo_temp press thermo_press

dump myimageSCC all image 500 dump.*.jpg type type atom yes size 1024 1024 view 30 60 center d 0.5 0.5 0.5 up 0 0 1 zoom 1.0 box yes 0.02 axes yes 0.05 0.1 subbox no 0.0 shiny 1.0

dump_modify myimageSCC append no element Fe Cr Al first yes flush yes thresh none acolor * red/blue/green adiam 1 2.48 adiam 2 2.50 adiam 3 2.86 backcolor black boxcolor white

# 5-Minimization

min_style cg

min_modify dmax 0.1 line quadratic

minimize 1.0e-4 1.0e-6 100000 1000000

unfix BR

# 6-Equilibration

velocity all create 600 13425 loop local

fix NVT1 all nvt temp 600 600 100 drag 0.2

fix SP all reax/c/species 1 1 100 mole.out element Fe Cr Al

timestep 0.20

restart 100 restart.*.stgb1

run 1000

unfix NVT1

Part-2 log.lammps

LAMMPS (31 Mar 2017)
OMP_NUM_THREADS environment is not set. Defaulting to 1 thread. (…/comm.cpp:90)
using 1 OpenMP thread(s) per MPI task

1-Initialization

units real

atom_style charge

atom_modify id yes map hash first alloy

boundary p p p

newton on

dimension 3

2-Creation of Fe82Cr13Al5 bicrystal

lattice bcc 2.866
Lattice spacing in x,y,z = 2.866 2.866 2.866
region box block 0 15 -10 10 0 5 side in
create_box 3 box
Created orthogonal box = (0 -28.66 0) to (42.99 28.66 14.33)
2 by 2 by 1 MPI processor grid

region left block 0 15 -10 0 0 5 side in
lattice bcc 2.866 orient x 5 -9 0 orient y 9 5 0 orient z 0 0 1
Lattice spacing in x,y,z = 3.89719 3.89719 2.866
create_atoms 3 region left
Created 1505 atoms

region right block 0 15 0 10 0 5 side in
lattice bcc 2.866 orient x 5 9 0 orient y -9 5 0 orient z 0 0 1
Lattice spacing in x,y,z = 3.89719 3.89719 2.866
create_atoms 3 region right
Created 1505 atoms

group alloy region box
3010 atoms in group alloy
group grain1 region left
1520 atoms in group grain1
group grain2 region right
1520 atoms in group grain2

set group alloy type/fraction 1 0.82 12963 type/fraction 2 0.13 23486 type/fraction 3 0.05 23486 # Fe75Cr20Al5
2454 settings made for type/fraction
397 settings made for type/fraction
143 settings made for type/fraction

mass 1 55.85000
mass 2 52.00000
mass 3 26.98000

change_box all x scale 1.2 y scale 1 z scale 1
orthogonal box = (-4.299 -28.66 0) to (47.289 28.66 14.33)
orthogonal box = (-4.299 -28.66 0) to (47.289 28.66 14.33)
orthogonal box = (-4.299 -28.66 0) to (47.289 28.66 14.33)

3-Force field

pair_style reax/c lmp_control safezone 4.0

pair_coeff * * ffield.reax Fe Cr Al

fix QEQ all qeq/reax 1 0.0 10.0 1.0e-6 reax/c

fix_modify QEQ energy yes

delete_atoms overlap 0.35 grain1 grain2 compress yes
Neighbor list info …
update every 1 steps, delay 10 steps, check yes
max neighbors/atom: 2000, page size: 100000
master list distance cutoff = 12
ghost atom cutoff = 12
binsize = 6, bins = 9 10 3
3 neighbor lists, perpetual/occasional/extra = 2 1 0
(1) command delete_atoms, occasional
attributes: full, newton on
pair build: full/bin/atomonly
stencil: full/bin/3d
bin: standard
(2) pair reax/c, perpetual
attributes: half, newton off, ghost
pair build: half/bin/newtoff/ghost
stencil: half/ghost/bin/3d/newtoff
bin: standard
(3) fix qeq/reax, perpetual, copy from (2)
attributes: half, newton off, ghost
pair build: copy
stencil: none
bin: none
Deleted 15 atoms, new total = 2995

fix BR all box/relax iso 0.0

neighbor 2.5 bin

neigh_modify every 1 delay 0 check yes once no cluster no exclude none page 100000 one 2000 binsize 0.0

4-Compute & thermo setting

compute total_potential_energy all pe

compute total_kinetic_energy all ke

compute reax all pair reax/c

thermo 250

thermo_style custom step time atoms temp press pe ke etotal epair

thermo_modify lost error lost/bond error norm yes flush yes line multi format float %10.4g temp thermo_temp press thermo_press

dump myimageSCC all image 500 dump.*.jpg type type atom yes size 1024 1024 view 30 60 center d 0.5 0.5 0.5 up 0 0 1 zoom 1.0 box yes 0.02 axes yes 0.05 0.1 subbox no 0.0 shiny 1.0

dump_modify myimageSCC append no element Fe Cr Al first yes flush yes thresh none acolor * red/blue/green adiam 1 2.48 adiam 2 2.50 adiam 3 2.86 backcolor black boxcolor white

5-Minimization

min_style cg

min_modify dmax 0.1 line quadratic

minimize 1.0e-4 1.0e-6 100000 1000000
Neighbor list info …
update every 1 steps, delay 0 steps, check yes
max neighbors/atom: 2000, page size: 100000
master list distance cutoff = 12.5
ghost atom cutoff = 12.5
binsize = 6.25, bins = 9 10 3
2 neighbor lists, perpetual/occasional/extra = 2 0 0
(1) pair reax/c, perpetual
attributes: half, newton off, ghost
pair build: half/bin/newtoff/ghost
stencil: half/ghost/bin/3d/newtoff
bin: standard
(2) fix qeq/reax, perpetual, copy from (1)
attributes: half, newton off, ghost
pair build: copy
stencil: none
bin: none
Per MPI rank memory allocation (min/avg/max) = 516.4 | 517.7 | 518.8 Mbytes
---------------- Step 0 ----- CPU = 0.0000 (sec) ----------------
Step = 0 Time = 0 Atoms = 2995
Temp = 0 Press = 5.947e+04 PotEng = -73.63
KinEng = 0 TotEng = -73.63 E_pair = -73.63
---------------- Step 1 ----- CPU = 4.8355 (sec) ----------------
Step = 1 Time = 1 Atoms = 2995
Temp = 0 Press = 5.902e+04 PotEng = -73.63
KinEng = 0 TotEng = -73.63 E_pair = -73.63
Loop time of 4.83944 on 4 procs for 1 steps with 2995 atoms

46.4% CPU use with 4 MPI tasks x 1 OpenMP threads

Minimization stats:
Stopping criterion = energy tolerance
Energy initial, next-to-last, final =
-73.6304892688 -73.6304892688 -73.6348942056
Force two-norm initial, final = 120792 120044
Force max component initial, final = 110248 109434
Final line search alpha, max atom move = 9.07043e-10 9.92614e-05
Iterations, force evaluations = 1 1

MPI task timing breakdown:
Section | min time | avg time | max time |%varavg| %total

*Dear Lammps Community,*

Please forgive me to again report a strange phenomenon resulted from fix
reax/c/species command.

I built an bicrystal alloy model consisting of total 3010 atoms.

By the atom type setting command: set group alloy type/fraction 1 0.82
12963 type/fraction 2 0.13 23486 type/fraction 3 0.05 23486,

I intend to set 82% Fe atom, 13% Cr atom and 5% Al atom.

And the log.lammps shows the following reasonable results:

  2454 settings made for type/fraction

  397 settings made for type/fraction

  143 settings made for type/fraction

However, the fix reax/c/species command output a different result:
Fe2118Cr251Al624.

Although I enabled the reax/c with kokkos package by the following
changes, the same problem occurred again.

​yes, because the problem is *not* with fix reax/c/species but falls into
the PEBCAC category, i.e. it resides between the chair and the computer.

the "fraction" is always applied to the group of the set command (which is
all atoms). so with the second type/fraction set will *also* consider
atoms, that have already been changed to type 1. you can confirm, that fix
reax/c/species is working correctly, by issuing the following 3 commands
and checking the log file:

group Fe type 1
2131 atoms in group Fe
group Cr type 2
254 atoms in group Cr
group Al type 3
625 atoms in group Al

to explain. the first type/fraction command does change 2454 atoms from
type 3 to type 1, the second command changes 397 atoms to type 2, however
quite a few of them were of type 1 and not type 3. the last command makes
things even worse, as you now change 5% of all 3010 atoms to type 3, and -
again - it is most likely that you change atoms of type 1 or 2.

the correct way to go about this, is the following:

set group all type/fraction 1 0.82 12963

​# now there are about 18% atoms left as type 3
# you want a ratio of 13:5 Cr to Al, i.e. a fraction of 13/18 = 0.72222
needs to be converted from Cr to Al
​# make a group that only contains type 3 atoms.

group CrAl type 3
set group CrAl type/fraction 2 0.722222222 23486

group Fe type 1
2454 atoms in group Fe
group Cr type 2
39
​4​
atoms in group Cr
group Al type 3
16
​2​
atoms in group Al

​# that is: 81.5% Fe, 13.1% Cr and 5.4% Al​

or even better, do it the other way around, i.e. create all atoms as type
1, then convert 18% of those to type 2 and 27.7777% of the type 2 to type 3.
with that i get:

set group all type/fraction 2 0.18 12963
Setting atom values ...
  508 settings made for type/fraction
group CrAl type 2
508 atoms in group CrAl
set group CrAl type/fraction 3 0.2777777777 23486
Setting atom values ...
  149 settings made for type/fraction

group Fe type 1
2502 atoms in group Fe
group Cr type 2
359 atoms in group Cr
group Al type 3
149 atoms in group Al

​-> 83.1% Fe, 11.9%​ Cr, 5% Al

mystery solved. LAMMPS is working fine, but your input script needs an
upgrade. :wink:

axel.

Dear Prof. Kohlmeyer,

Thanks so much and it really works!

Thank you very very much for you continuous help me, a LAPPMS beginner!

Best regards,

Xiaolong