The number of atoms created in the box is larger than input + PotEng: Nan - NaYF4

LAMMPS LAMMPS Beginners
1

Hi all,
I am Raj Kumar Ramamoorthy, an experimental nanoscience researcher working on the nucleation and growth of nanoparticles. To complement my experimental work, I am beginning to use the LAMMPS simulation tool. My initial project involves simulating the stable polymorph crystal structure of NaYF4, for which I have constructed a cubic phase with Na at the corners, Y at face-centered positions, and F at tetrahedral sites.

3d NaYF4 minization

#–Initialization
units metal
#The units command also sets the timestep size and neighbor skin distance to default values for #each style: For style real these are dt = 1.0 femtoseconds and skin = 2.0 Angstroms.
dimension 3
boundary p p p
atom_style atomic

#–Variable
variable latpara equal 5.47
variable material string NaYF4

#-- Atom definition and creation of simulation box

lattice custom ${latpara} &
a1 1 0 0 &
a2 0 1 0 &
a3 0 0 1 &
basis 0 0 0 &
basis 0.5 0.5 0 &
basis 0.5 0 0.5 &
basis 0 0.5 0.5 &
basis 0.25 0.25 0.25 &
basis 0.75 0.25 0.25 &
basis 0.25 0.75 0.25 &
basis 0.75 0.75 0.25 &
basis 0.25 0.25 0.75 &
basis 0.75 0.25 0.75 &
basis 0.25 0.75 0.75 &
basis 0.75 0.75 0.75 &

region box block 0 2 0 2 0 2

unit cell

create_box 3 box

3 kinds of atoms

create_atoms 1 box &
basis 1 1

fill with 1 atoms # no & is needed

create_atoms 2 box &
basis 2 2 &
basis 3 2 &
basis 4 2

fill with 2 atoms

create_atoms 3 box &
basis 5 3 &
basis 6 3 &
basis 7 3 &
basis 8 3 &
basis 9 3 &
basis 10 3 &
basis 11 3 &
basis 12 3

fill with 3 atoms

mass 1 22.9898
mass 2 88.9059
mass 3 18.9984

pair_style buck 2.5
pair_coeff 1 1 0 1 0
pair_coeff 1 2 0 1 0
pair_coeff 1 3 1094.2 0.2709 0
pair_coeff 2 2 0 1 0
pair_coeff 2 3 1568.9 0.2950 0
pair_coeff 3 3 1127.7 0.2753 15.83

#minimize 1.0e-4 0.001 10 100

energy force iteration

velocity all create 300 1029329 dist gaussian

fix 1 all nvt temp 300 300 100

thermo 5
thermo_style custom step temp pe
run 10
unfix 1

While the above code runs without errors, I am having difficulty interpreting the log file and output.

  1. My expectations were to generate 8 unit cells, resulting in 64 Na, 48 Y, and 64 F atoms. However, the simulation yields 96 atoms for each element.
  2. Furthermore, the potential energy at the 0th step is NaN, and the temperature tends towards NaN during the simulation.

I would be grateful for any guidance on these outputs to refine my simulation scripts.

2

Hello,

Please format your input using the preformatted text option – it makes the input much more readable. Additionally, you’re not following the forum rules. For example, you didn’t provide the LAMMPS version you’re using.

That can happen when your potential is incorrect or when there is overlap between atoms.
My advice is to print the positions of your atoms using:

write_dump all atom dump.lammpstrj
run 0

Then, double-check your structure using VMD or Ovito to make sure everything looks as expected.

Simon

3

Thank you for the prompt response. My apologies for the rule violation. I am running LAMMPS version (29 Aug 2024 - Update 2). Please see below for the formatted script.

# 3d NaYF4 minization

#--Initialization
units		metal 
#The units command also sets the timestep size and neighbor skin distance to default values for #each style: For style real these are dt = 1.0 femtoseconds and skin = 2.0 Angstroms.
dimension	3
boundary	p	p	p
atom_style	atomic

#--Variable
variable latpara equal 5.47
variable material string NaYF4

#-- Atom definition and creation of simulation box

lattice custom ${latpara} &
a1 1 0 0 & 
a2 0 1 0 &
a3 0 0 1 &
basis 0 0 0 &
basis 0.5 0.5 0 &
basis 0.5 0 0.5 &
basis 0 0.5 0.5 &
basis 0.25 0.25 0.25 &
basis 0.75 0.25 0.25 &
basis 0.25 0.75 0.25 &
basis 0.75 0.75 0.25 &
basis 0.25 0.25 0.75 &
basis 0.75 0.25 0.75 &
basis 0.25 0.75 0.75 &
basis 0.75 0.75 0.75 &

region	box block 0 2 0 2 0 2
# unit cell

create_box	3 box 
# 3 kinds of atoms

create_atoms	1 box &
basis 1 1 

# fill with 1 atoms # no & is needed

create_atoms	2 box &
basis 2 2 &
basis 3 2 & 
basis 4 2
# fill with 2 atoms

create_atoms	3 box &
basis 5 3 &
basis 6 3 &
basis 7 3 &
basis 8 3 &
basis 9 3 &
basis 10 3 &
basis 11 3 &
basis 12 3 
# fill with 3 atoms


mass 1 22.9898
mass 2 88.9059
mass 3 18.9984


pair_style	buck 2.5
pair_coeff	1 1 0 1 0 
pair_coeff	1 2 0 1 0
pair_coeff	1 3 1094.2 0.2709 0
pair_coeff	2 2 0 1 0
pair_coeff	2 3 1568.9 0.2950 0
pair_coeff	3 3 1127.7 0.2753 15.83

#minimize	1.0e-4 0.001 10 100	
# energy force iteration

velocity all create 300 1029329 dist gaussian

fix		1 all nvt temp 300 300 100

thermo		5
thermo_style custom step temp pe
write_dump all atom dump.lammpstrj
run 0
#run		1 start 0 stop 10
unfix 		1

As suggested by you, I’ve examined the positional coordinates of the Na, Y, and F atoms. Interestingly, all atoms, regardless of their type, are occupying identical positional coordinates. I believe the issue might stem from atom overlap. In this situation, what is the best way to overcome this? Specifically, how can I remove selected atoms from the simulation box?

4

A cutoff of 2.5 Angstrom is quite unreasonable for metal units. Where did you take it from?

When I add the command:

delete_atoms overlap 0.2 all all

It deletes 192 atoms, so those are all on top of the remaining 96.

I don’t quite understand why you have to have three create_atoms commands. Shouldn’t there be just one? You should re-read the documentation for that command very carefully.

5

As has been stated many times before. You should not do this on your own, but need to find an experienced collaborator that can train and advise you. There is more to correct and meaningful simulations, than sticking some commands into a software and have it (or us) tell you what is wrong. You would not give me access to your experimental equipment without proper training by a practitioner, do you?

6

Thank you for the suggestions. I agree that an experienced collaborator would be very helpful. I was initially overambitious to think I could do this alone. Until I find someone, I will continue to learn and advance in this direction.

Since I’m working with an ionic system, I’ve restricted the cutoff to 2.5. Your comments on “delete_atoms” and “create_atoms” as well with “set type/fraction” command, I am now able to develop an improved script and move forward with the minimization.

7

That makes physically no sense at all. The cutoff will be determined by the force field. 2.5 angstrom is shorter than your effective particle diameter!
And if you have an ionic system, you likely will need to use a long-range coulomb solver.

8

Thanks for your feedback. I will include a longer/long-range coulomb cutoff distance.