Trying to make the simulation periodic and avoid mismatched atoms

[Dalia]

Hello,

I wrote this code to calculate the stacking fault energy for NiCrFe alloy, when I change the line from p p s to p p p and visualize in ovito it shrinks like that attached, how can I make it all periodic? Also how to avoid grey atoms resulting. Thank you so much!

------------------------ INITIALIZATION ----------------------

units metal
dimension 3
boundary p p s
atom_style atomic
variable latparam1 equal 3.52 # Adjusted lattice parameter for Ni

variable xdim equal {latparam1}*sqrt(6)/2*10 variable ydim equal {latparam1}sqrt(2)/210

----------------------- ATOM DEFINITION ----------------------

lattice fcc {latparam1} region 1 block 0 {xdim} 0 {ydim} 0 20 region 2 block 0 {xdim} 0 {ydim} 20 40 region whole block 0 {xdim} 0 {ydim} 0 200 units box create_box 3 whole # Specify 3 atom types lattice fcc {latparam1} orient x 1 1 2 orient y -1 1 0 orient z -1 -1 1

Create atoms for NiCrFe 60:20:20 alloy in region 1

create_atoms 1 region 1
set type 1 type/fraction 2 0.20 12345
set type 1 type/fraction 3 0.20 12346 # Adjusted to maintain the 60:20:20 ratio

Create atoms for NiCrFe 60:20:20 alloy in region 2

create_atoms 1 region 2
set type 1 type/fraction 2 0.20 12347
set type 1 type/fraction 3 0.20 12348 # Adjusted to maintain the 60:20:20 ratio

----------------------- FORCE FIELDS -----------------------

pair_style eam/alloy
pair_coeff * * /home/dalia/FeCrNi_d.eam.alloy Ni Cr Fe

------------------------- SETTINGS --------------------------

group top region 1
group bot region 2

Compute the initial energy

compute peratom all pe/atom
compute eatoms all reduce sum c_peratom

Dump to comp for Ovito post processing

dump 1 all custom 1 dump.comp.* id type xs ys zs c_peratom fx fy fz

thermo 1
thermo_style custom step pe c_eatoms

Minimize energy before displacement

fix 1 all setforce 0 0 NULL
min_style cg
minimize 1e-10 1e-10 1 1

Store the initial energy

variable initial_energy equal c_eatoms
run 0
variable stored_initial_energy equal ${initial_energy}

------------------------- Displacement -----------------------

displace_atoms bot move -1.0 0.0 0.0 units box

Minimize energy after displacement

fix 1 all setforce 0 0 NULL
min_style cg
minimize 1e-10 1e-10 1 1

Store the final energy

variable final_energy equal c_eatoms
run 0
variable stored_final_energy equal ${final_energy}

Calculate the stacking fault energy in eV/Å^2 and convert to mJ/m^2

variable area equal {xdim}*{ydim}
variable diff_energy equal “v_stored_final_energy - v_stored_initial_energy”
variable stacking_fault_energy equal “v_diff_energy / v_area * 16021.7733” # Conversion factor from eV/Å^2 to mJ/m^2

Print the energies and stacking fault energy

print “Initial energy: {stored_initial_energy} eV" print "Final energy: {stored_final_energy} eV”
print “Stacking fault energy: ${stacking_fault_energy} mJ/m^2”

SIMULATION DONE

print “All done”

[hothello]

This is your first post, so Please Read This First: Guidelines and Suggestions for posting LAMMPS questions. Also, please edit your post and enclose the input script in the preformatted text environment.

Said so, there is nothing wrong in your simulation. Ovito just wraps atoms at the periodic boundaries, which means that the grey atoms on top of the box have z < zlo.
If you want to visualise the slab without breaking, you need to either decrease the zlo value, or shift the whole system upwards, or use the unwrap modifier in Ovito.

[Dalia]

Thank you so much! I tried decreasing zlo to 90 but I still get grey atoms as shown attached, also, I noticed that when I lower zlo to 40 or 50, the stacking fault comes out a negative very small number, why is that?

[Dalia]

I’m sorry it doesn’t let me edit the post here’s the formatted code:

%%writefile NiCrFe.in
# This LAMMPS input script calculates the stacking fault energy of a NiCrFe 60:20:20 alloy.
# The script initializes the system with metal units, 3 dimensions, and periodic boundaries
# in the x and y directions, while being shrink-wrapped in the z direction. It sets the
# lattice parameter and calculates the dimensions of the simulation box.

# Regions are defined for the atoms, and the lattice orientation is specified. Atoms are
# created and distributed in two regions with a 60:20:20 ratio of Ni, Cr, and Fe. The
# initial energy of the system is computed and minimized, and this energy value is stored.

# Atoms in the bottom region are then displaced, and the energy of the system is minimized
# again. The final energy is stored, and the difference between the final and initial
# energies is used to calculate the stacking fault energy, which is then converted to 
# mJ/m^2. Finally, the script prints the initial energy, final energy, and stacking fault
# energy and then indicates that the simulation is done.

# ------------------------ INITIALIZATION ----------------------
units metal
dimension 3
boundary p p s
atom_style atomic
variable latparam1 equal 3.52  # Adjusted lattice parameter for Ni

variable xdim equal ${latparam1}*sqrt(6)/2*10
variable ydim equal ${latparam1}*sqrt(2)/2*10

# ----------------------- ATOM DEFINITION ----------------------
lattice fcc ${latparam1}
region 1 block 0 ${xdim} 0 ${ydim} 0 20 
region 2 block 0 ${xdim} 0 ${ydim} 20 40
region whole block 0 ${xdim} 0 ${ydim} 0 200 units box
create_box 3 whole  # Specify 3 atom types
lattice fcc ${latparam1} orient x 1 1 2 orient y -1 1 0 orient z -1 -1 1

# Create atoms for NiCrFe 60:20:20 alloy in region 1
create_atoms 1 region 1
set type 1 type/fraction 2 0.20 12345  # Type 2 is Cr
set type 1 type/fraction 3 0.20 12346  # Type 3 is Fe

# Create atoms for NiCrFe 60:20:20 alloy in region 2
create_atoms 1 region 2
set type 1 type/fraction 2 0.20 12347  # Type 2 is Cr
set type 1 type/fraction 3 0.20 12348  # Type 3 is Fe

# Set remaining fraction for Ni (Type 1)
set type 1 type/fraction 1 0.60 12349  # Type 1 is Ni

# ----------------------- FORCE FIELDS -----------------------
pair_style eam/alloy
pair_coeff * * /home/dalia/FeCrNi_d.eam.alloy Ni Cr Fe

# ------------------------- SETTINGS --------------------------
group top region 1
group bot region 2

# Compute the initial energy
compute peratom all pe/atom
compute eatoms all reduce sum c_peratom

thermo 1
thermo_style custom step pe c_eatoms

# Minimize energy before displacement
fix 1 all setforce 0 0 NULL
min_style cg
minimize 1e-10 1e-10 1000 1000

# Dump the configuration after minimization before displacement
dump dump_before all custom 1 dump.before id type xs ys zs c_peratom fx fy fz
run 1
undump dump_before

# Store the initial energy
variable initial_energy equal c_eatoms
run 0
variable stored_initial_energy equal ${initial_energy}

# ------------------------- Displacement -----------------------
displace_atoms bot move -1.0 0.0 0.0 units box

# Minimize energy after displacement
fix 1 all setforce 0 0 NULL
min_style cg
minimize 1e-10 1e-10 1000 1000

# Dump the configuration after minimization after displacement
dump dump_after all custom 1 dump.after id type xs ys zs c_peratom fx fy fz
run 1
undump dump_after

# Store the final energy
variable final_energy equal c_eatoms
run 0
variable stored_final_energy equal ${final_energy}

# Calculate the stacking fault energy in eV/Å^2 and convert to mJ/m^2
variable area equal ${xdim}*${ydim}
variable diff_energy equal "v_stored_final_energy - v_stored_initial_energy"
variable stacking_fault_energy equal "v_diff_energy / v_area * 16021.7733"  # Conversion factor from eV/Å^2 to mJ/m^2

# Print the energies and stacking fault energy
print "Initial energy: ${stored_initial_energy} eV"
print "Final energy: ${stored_final_energy} eV"
print "Stacking fault energy: ${stacking_fault_energy} mJ/m^2"

# SIMULATION DONE
print "All done"

[Dalia]

@hothello Hello, thank you so much for your help! I tried shifting the whole system upwards but as I originally have it in two Ni regions, the type 1 region which is in red appear also on top with an interface with the type 2 region in blue which is wrong not FCC nor a stacking fault, here are visualizations before and after the shift respectively:

After the shift: