Lammps Cascade Simulations

LAMMPS
1

Dear experts,

I am running consecutive PKA cascades on the same simulation box with periodic boundaries on all sides to accumulate damage. The choice of PKA atom is randomized and uniform irradiation is also ensured. There are thermostats on thin boundary layers on all sides maintained at 300 K. PKA energy is like 1 keV and box size is like 20x20x20. Time step is 1 fs with adaptive modification. A Restart file is written after every cascade is run, which is read for the next cascade to begin. Each cascade is sufficiently allowed to evolve so that any additional energy in the system is minimum. Everything is running fine, but when the point defects are analyzed (using the OVITO - Wigner-Seitz method) there are oscillations in the numbers of point defects, i.e., increasing slowly in a linear graph then suddenly becomes very high number (10 to 100 times a value) reaches a maximum value and then decreases to a very low value - this pattern is continuing with the increasing number of cascades. The overall number of vacancies / interstitials as a result do not increase. I do not find any explanation for this.

Please help me.

Thank you for your time.
Regards.

2

In general existing defects can anneal out if additional energy (subsequent cascades)
is added to the system. Since it sounds like you are running a series of independent
simluations, I suggest you isolate the one that is causing a dramatic increase
in defect count. Verify that the initial and final snapshot are what you observe in the OVITO plot.
Check the logfile for that run, see if there is anything anomolous. Visualize just that simulation. Etc.

Steve

3

Thank you for the suggestion. I shall do as you suggested and give the feedback.
Regards,
uas

4

For example in the 16 th. cascade the number of vacancies is 30 (increasing slowly till this) then in the 17 th. cascade the number goes up to 240 and then in the next cascade it is something around 3000, then 3540 and then again around 34. When this happens for a long series of about 1200 successive cascades it becomes really difficult to again identify the problematic cases and run them individually from the restart files.

I tried running these cases separately, but sometimes just got this behavior repeated without any remedy and sometimes got the numbers as expected.

There are no differences observed (apart from the variables that are estimated using a random number seed) in the log files between two runs that produced an expected number and that gave an anomalous one.

Your suggestions will be very helpful.

Thank you.

5

What degree of radiation is it that you want to represent? How large is your system?
If you are modeling a large number of incidents, at some point your overall structure is no longer easily detectable and that may be throwing off the analysis.
1200 incidents on the same system sounds quite excessive to me for typical size MD systems. Even 100.
Are you certain that you are following the correct simulation protocol?
Could it be that you should be making the restart before selecting the PKA?

6

My simulation cell has 50 unit cells on all three sides. PKA energy is 1 keV. To get at least 0.1 dpa therefore I need to simulate in the order of 4000 successive cascades. As I mentioned before, I am doing an uniform irradiation simulation over the system. A restart binary is prepared after equilibration. Then after every cascade a restart binary is made to start the next cascade from it. A PKA is selected from the configuration in the restart file to begin the cascade.

7

What do you mean by you repeat an individual case. Sometimes it repeats the behavior
in the successive simulatoin, sometimes it does not. Are you saying an individual simulation
is not repeatable, running on the same # of processors? Why is that?

Steve

8

No. Let me clearly state it again. I am saving the restart files. So when I run just the nth. cascade only from the saved (n-1)th. saved data (the # of processors and all other setup is unchanged) and visualize the outcome, sometimes the number of interstitials/vacancies come in the correct order as expected from the numbers obtained by visualizing the result of the previous, i.e., (n-5)th., (n-4)th, …, (n-1)th. cascades. But sometimes the anomalously blowing up of the number of vacancies / interstitials by two to three orders of magnitude again occurs even in scenario of running only the problematic case alone.

9

By “repeatable” I meant if you run a single simulation multiple times (one that produces a zillion vacancies) do you always get the same answer. Running it each time by itself, indepedent of the original sequence of simulations. If not, I’d like to understand why. If yes, I suggest doing thermo output every step in that run and looking at it to see if there any anomalies. And are you saying
that if you use OVITO on only the output from that single simulation the vacancy count goes from 240 → 3000 ?

A related Q. Is OVITO or your definition of a vacancy robust to the entire box of atoms shifting in x,y, or z. If your box is periodic on all sides, how are you preventing the entire system from drifting in the direction of momentum you’ve added? And does the 3000 vacancy snapshot look normal or odd?

10

Yes, it is independent. However, if I get one or two orders of difference in the number when running the same nth. cascade from the (n-1)th. cascade system history stored in restart file, then I do not know what to tell about that.

Yes, after increasing slowly the numbers abruptly increase in a few consecutive cascades and then again abruptly decreases to the initial values. This thing repeats like a pattern throughout the number of cascades I run.
In OVITO, the Wigner-Seitz analysis is used.

“If your box is periodic on all sides, how are you preventing the entire system from drifting in the direction of momentum you’ve added?”
I do not know how to do/ensure that. If you mean to use any particular command or fix, please suggest me then I shall definitely look into that.
I have tested with normal x, y, z, scaled and unwrapped coordinates.
Number of particles in the box are conserved.
“And does the 3000 vacancy snapshot look normal or odd?”
When the numbers are in this kind of sequence: 3, 5, 5, 8, 13, 15, 90, 600, 2000, 3400, 1600, 50, 8, 10, 7, 16, 26, 40, 42, 38, 280, 1900, 5040, 50000, 6040, 2342, 54, 20, 32, … then it is odd. When the numbers are abruptly high suddenly I see the box is full with vacancies and interstitials.
Everything else but the point defect numbers are coming correctly.

11

ok - thanks for the info

You’ve (mostly) convinced me that your anomalous defect counts are not
coming from:

a) a bias you mistakenly created by running consecutive sims one after the other
b) issues with the sequence of consecutive frames you are feeding to OVITO,
e.g. making it think atoms displaced further than they did

So a natural Q is how do you know anything is wrong? Maybe these highly
varying stats are the nature of the physics?

I think your sequence of numbers: 3,5,5,8,etc are the defect
counts after successive cascades.

I suggest picking an anomalous pair, like when 15 → 90 due to a
single knock-on simulation. Take the 15-defect state as a starting
point. Run 10 or 100 one-cascade sims from that starting point. With
different random atoms chosen as the knock-on, different random
directions for the kick direction, etc. Maybe that’s as simple as
setting a different RNG seed for each of the 100 sims.

Then look at the defect stats for the 100 outcomes. Do 10% have huge
increases and 90% not? Or vice versa?

Do the same for a starting state like 1600 → 50 where the defect
count drops dramatically.

If big changes never happen, then it seems likely there is in fact
something wrong about the way you are running N simulations one after
the other.

If big changes happen occasionally, then maybe that’s the way the
physics works, and your “odd” sequence of defect counts is just
sampling from those stats. E.g. it seems like you might observe a big
difference if you happen to pick a knock-on direction that shoots an
atom directly at one it’s nearest neighbors.

And also (I’ve already suggested this), pick a one-cascade simulation
where there is a big increase in defect count and examine its logfile
thermo output (with every-step output turned on) and verify there is
nothing screwy about it. E.g. maybe your adaptive timestep settings
are off, leading to huge spikes in temperature or energy.

Also, the typical way to ensure no drift of a simulation box you are
adding a big momentum spike to (repeatedly), is to freeze a plane of
atoms (or even a single atom), say the bottom z-layer or a bottom
corner atom. Presumably your cascade is entirely contained in the
middle of the box. Or more mildly, attach a spring to the bottom
plane or bottom corner atom to keep it roughly in the same position.

Steve

12

Many thanks for these suggestions and explanations. I shall try to implement these and give you the feedback.

Regards.