[lammps-users] Recoil atom simulation question

Hi,

I am trying to model the effect of a PKA atom on a FCC metal. For the moment I decided to fix the box T to 300 K w/ a npt ensemble. Then thermostat to 300 K (Nosé-Hoover or Langevin) just the borders of the box while the cascade does its evolution.

Dumping the results I noticed that if I unfix the first npt fix before choosing the recoil atom and setting its velocity, the atom does not move and the cascade does not occur. It looks like the atom needs a temperature fixed before its velocity (recoil energy) is set. Is it normal?

If so, I need to run an npt fix on the whole box and a nvt/langevin thermostat on the borders simultaneously. The whole box does have a fixed T then. In this case, does the fixed temperature of the atoms interfere with the evolution of the cascade? I mean, if I set the npt fix (whole box) to run every, say, 100 timesteps, does the cascade get its evolution just for 100 timesteps maximum and then the npt fix just “freezes” all the atoms back to 300 K?

Thank you.

Cheers,

Stefano

Dear Stefano,

You should not do an NPT of the region where the cascade is progressing - the temperature control will affect the cascade evolution. However, as you mention, the borders can be temperature controlled. PFA a sample input for simulating a collision cascade simulation in FeCr where the borders are temperature controlled. You will have to give the appropriate Vx, Vy and Vz corresponding to the PKA energy. I hope this helps.

Regards
Manoj

in.pkaSample (2.19 KB)

Hi. Please see my point-wise reply between your text below:

  1. You used the npt ensemble before the cascade simulation and you a unfixed it (one of my mistakes). Still in such fix, I saw you set a temperature damping parameter of 0.1 timesteps even tough it is suggested to stay on the order of 100 timesteps. Is your choice an active part of the modeling of the cascade? Why did you chose a short “Tdamp”?

The timestep is femo-seconds and the Tdamp is 0.1 pico-seconds (= 100 timesteps)

  1. There is an nve ensemble in fixed for the central box. Keeping constant the energy in the central region doesn’t “amplify” the cascade effect in time?

Why should it? an NVE ensemble does not actively change either the velocities or the positions.

  1. I noticed you used a compute nfrenkel and a compute Frenkel/local. However I don’t see such compute commands in lammps manual (is that something you developed by yourself?)

Oops sorry. That is something we developed ourselves where we calculate the Frenkel pairs created during the cascade and output only their co-ordinates - helpful in reducing MD output. The extra time taken to calculate the defects during the run is lesser than the time spent by the program to write out all the position co-ordinates for post-processing.

  1. The choice of the “neigh_modify delay 5” is still an active choice of the model? Is that useful in a cascade simulation?

I believe it is useful since the neighbors in a collision cascade can change and the neighbor list must be updated.

Please note that if you are carrying out a high energy PKA simulation (few tens of keVs) you might want to use the “fix ttm” or the “fix electron/stopping” commands in LAMMPS. That is also not included in the example file I sent you.

Regards
Manoj

Hi Manoj,

Thank you very much for your help and your input example. It is exactly what I needed. I was indeed doing something wrong with the thermostats. I’d like to ask you another one or two questions if you don’t mind:

  1. You used the npt ensemble before the cascade simulation and you a unfixed it (one of my mistakes). Still in such fix, I saw you set a temperature damping parameter of 0.1 timesteps even tough it is suggested to stay on the order of 100 timesteps. Is your choice an active part of the modeling of the cascade? Why did you chose a short “Tdamp”?
  2. There is an nve ensemble in fixed for the central box. Keeping constant the energy in the central region doesn’t “amplify” the cascade effect in time?
  3. I noticed you used a compute nfrenkel and a compute Frenkel/local. However I don’t see such compute commands in lammps manual (is that something you developed by yourself?)
  4. The choice of the “neigh_modify delay 5” is still an active choice of the model? Is that useful in a cascade simulation?

Thank you again.

Stefano

Thank you Manoj,

I found your suggestions and explanations very useful. I have just one last question that is a problem I found running my “complete” model with your implementations as well.

I am trying to run some “multiple cascade” simulations on a NiFe alloy with PKA energies of 5 keV (1500 cascades – about 0.5 dpa, box of 30x30x30 – 108’000 atoms). I am using the following algorithm:

Loop pka 1500

Read_restart restart_file.*

Set velocity (300 K)

npt ensemble (300 K, 1 bar)

Run 100 ps (100’000 timesteps)

Minimize

Unfix npt

Hi.

50000 steps are indeed too much for a 5 keV collision cascade.

regarding abnormally high defects in some cases at 5 keV, I’m not sure why that is happening. Are you using stiffened potentials? I do notice that you are using variable time-step. A plot of the repulsive part of the inter-atomic potential may give you the answer.

Cheers!
Manoj

Hi,

I am not using stiffened potential (I am using this potential in particular: 2013–Bonny-G-Castin-N-Terentyev-D–Fe-Ni-Cr. However, I am running an alloy/overlay command as such:

pair_style hybrid/overlay eam/alloy zbl 0.35 0.5

pair_coeff * * eam/alloy /FeNiCr_Bonny_2013_ptDef.eam.alloy Fe Ni Cr

pair_coeff 1 2 zbl 26 28

neighbor 0.3 bin

neigh_modify delay 5

I’m not sure on how to plot the repulsive part of the potential in this case.

Thank you very much.

Cheers,

Stefano