Dear Lammps users,

I am trying to select a certain number of atoms as PKAs using the following commands:

group select empty

group PKAstore empty

variable numPKA equal count(PKAstore)

variable n loop 500

variable check atom grmask(Si,FREE)

label jump

variable select equal (random(1,atoms,4224):%.0f)
variable a equal (v_check[v_select])

if “v_check[v_select] == 1” then &

“group PKAstore id {select}"
if "{numPKA}<10” then &

“next n” &

“jump SELF jump”

variable n delete

However, I ran into the following error:

ERROR: Indexed per-atom vector in variable formula without atom map. (…/variable.cpp:4604)

Last command: variable a equal $(v_check[v_select])

I have looked up LAMMPS manual and added “atom_modify map yes id yes”in my input, yet the error still exists. I don’t know how to fix this atom mapping problem, and if there is a better way to select a certain number of atoms randomly within a specified group and region?

Many thanks for any suggestion!

Best wishes,

Dingyu

The whole in script is as follows:

variable T1 equal 300.0 # K

variable RDN equal 4928459

variable PROD equal 10000 #100ps test

variable dt1 equal 0.1

variable P1 equal 1.0

variable Tpka equal 55830.0

units real

boundary p p p

atom_style charge

atom_modify map yes id yes

read_restart PHPS50_slab_Eq300K.restart

pair_style reaxff lmp_control safezone 32 mincap 2000 minhbonds 1000

pair_coeff * * ffield N Si H

fix 1 all qeq/reaxff 1 0.00 10.00 1e-6 reaxff

neighbor 2.0 bin

neigh_modify every 10 delay 0 check yes

balance 1.1 shift z 10 1.05

timestep ${dt1}

reset_timestep 0

region BOTTOM block INF INF INF INF INF -16.0 units box

region FREE block INF INF INF INF -16.0 23.0 units box

region SURF block INF INF INF INF 10.0 23.0 units box

group BOTTOM region BOTTOM

group FREE dynamic all region FREE

group SURF region SURF

group FREEstatic region FREE

group N type 1

group Si type 2

group H type 3

group FREESi intersect Si FREEstatic

group SURFSi intersect Si SURF

#randomly select a Si PKA

group select empty

group PKAstore empty

variable numPKA equal count(PKAstore)

variable n loop 500

variable check atom grmask(Si,FREE)

label jump

variable select equal (random(1,atoms,4224):%.0f)
variable a equal (v_check[v_select])

if “v_check[v_select] == 1” then &

“group PKAstore id {select}"
if "{numPKA}<10” then &

“next n” &

“jump SELF jump”

variable n delete

velocity BOTTOM set 0.0 0.0 0.0

fix FREEZE BOTTOM setforce 0.0 0.0 0.0

fix spec all reaxff/species 1 1 100 species_0_1000_ps.out element N Si H &

delete species.del masslimit 0 47

fix bond all reaxff/bonds 1000 bonds_0_1000_ps.reax

compute mytemp FREE temp

compute_modify mytemp dynamic/dof yes

thermo_style custom step time atoms temp c_mytemp pe ke evdwl ecoul etotal press xlo xhi ylo yhi zlo zhi cpuremain

thermo 1000

dump 1 all custom 10000 PHPS50_slab_rad.lammpstrj id type q x y z vx vy vz

dump_modify 1 sort id

#setup PKA

velocity select set 0.3388362 0.3388362 0.3388362 units box

fix 2 all nve

fix 3 all dt/reset 1 1.0E-3 1.0E-1 0.2 units box

restart 1000000 PHPS50_slab_rad.restart

run ${PROD} # 100ps