Problem with defect count in consecutive cascade simulation

OVITO
1

Hello,
I know this is more of a LAMMPS related question but I’m not sure that if my methodelogy is wrong or my OVITO post processing is wrong. That’s why I’m posting it here.

I am having problem with cascade simulation. I am trying to do consecutive cascade to observe cumulative damage progression in FeNiCr FCC alloy. I’m here to seek help from experts. The methodology I’m following -

  1. The system is relaxed for 100ps via NPT ensemble and a snapshot is taken of the
    relaxed system.
  2. A random atom in the system is selected as PKA and the entire cell is shifted to move the PKA to the center to avoid the damage cascade to reach the system boundary.
  3. 5keV kinetic energy is given to the PKA in a random direction.
  4. During the cascade, the boundary layers of the cell (3 atom layers) are kept in NVT ensemble (Tdamp=0.1ps) to cool down the system back to 300K at the end of each cascade and the interior of the system is kept in NVE ensemble.
  5. During the cascade adaptive timestep is used with parameters Tmin=1e-7, Tmax=0.001ps, Xmax=0.05A and Emax=125eV.
  6. Electron stopping is applied for the atoms with energy greater than 10eV.
  7. I am running each cascade for around 33ps enough to cool down the cell to 300K before the next cascade starts.
  8. After each cascade run the atoms were shifted back over the periodic boundaries to keep the cell cantered at the origin in order to allow analysis, which uses as a reference the original location of the box.
  9. A snapshot of the system is taken as dump as each cascade finishes.
  10. Then steps from 2 to 9 are again repeated and this is done 300 times for 300 consecutive cascades.

The system size is 30ax30ax30a (108000 atoms) enough to contain the cascade within the cell boundaries. ZBL potential is used via hybrid overlay with eam/alloy potential from Bonny et al.

I am using Ovito’s Wigner Seitz analysis modifier the identify the point defects in each cascades using the relaxed system snapshot of the 1st step as reference. According to literature, the defect number should increase linearly at first for 30-35 cascades with a very steep slop. Then it should increase slowly. But I am getting random defect spikes in dump files like 6, 13, 20, 26, 30, 140, 16648, 15158, 7444, 1360, 723, 1593, 222, 82, 70.

To investigate further, I have taken dumps during a cascade. Plotting the defect vs time graph for each cascades gives the expected trend (rapid increase of defects at a very short time, then it reaches to a peak value, then recombination starts and the defects starts to fall and at the end some surviving defects remain).

The problem happens when I shift the cell back to it’s original position. After shifting it back, when compared to the initial relaxed system using Wigner Seitz modifier, the defect shows a spike. I’m using default settings of OVITO for Wigner Seitz analysis.

I have tried to use Berendsen temperature controlling instead of using NVT, have ran each cascades for longer times, used different dt/reset parameters, different temperature damping parameters. I even tried relaxing the system using both NVT and NPT ensemble after each cascades. All showed the same problem.

Has anyone faced this problem before? Any kind of advice from you guys would help a lot.

Thank you,
Tawseef

2

There are multiple things to consider when using the Wigner-Seitz defect analysis.

  • If the simulation cell changes size during the time series, you should use OVITO’s affine mapping option to compensate. In your case, such a size change could arise from the creation of defects during irradiation.
  • It’s important that the atomic positions of non-defective atoms match well between the reference configuration and the configuration of interest. That usually means you need to not only undo all the shifts you applied to your system but also compensate for any drift occurring in your system. This drift could be caused by the momentum created by the external irradiation.
3

I think the drifting was the reason for the unusual defect spikes. I have managed to eliminate the drifting. Thank you very much for your advice.

4

Could you please suggest what is the nature of this drift and how it can be compensated for in such a simulation? The momentum given to our chosen particle is distributed among all atoms and carries the system in a certain direction, right? So we should have followed the LCM from the beginning and created an alpha particle with momentum compensating the momentum of the recoil nucleus?

5

Atomic drift occurs during irradiation when momentum is added to the system, causing all atoms to drift. This can be compensated for either through your simulation software (depending on its capabilities) or by using a single atom as a reference point and maintaining its average position as constant.

Disclaimer: I’m by no means an expert on irradiation simulations and there might be better methods.

6

Thank you for your reply! I tried to compensate for the drift by simply subtracting its vector from the coordinates of all atoms. But it turned out that the number of Frenkel pairs calculated by the W-S method sharply increased to a certain value and remained constant for the whole simulation (2.5 ns)