Simulating quartz nanoparticles

GULP
1

Hi, i’m trying to do an md run with this nanoparticle cube of alpha quartz. I use the Sanders potentials. It keeps crashing when i try to run it as it gets too hot, the error message is the mdmaxtemp one.

below is a copy of the options i use. Theres around 2000 atoms altogether, and theyre all inputted manually, using a cif output from vesta and since theres so many I cant list them all here.

Thanks for any help!

#Keywords-----------------------------------------------

md conp

#cell---------------------------------------------------

cell
57.136250 57.136250 57.136250 90.0000 90.0000 90.0000

species 3
Si core 4.0
O core 0.8482
O Shel -2.8482

#potentials-------------------------------------------------

buck
O Shel O Shel 22764.000 0.1490 27.8800 0.0 10.0

buck
Si core O Shel 1283.9073 0.3205 10.6616 0.0 10.0

spring
O 74.9204

#md inputs--------------------------------------------------

ensemble npt 0.05 0.05

integrator leapfrog verlet

temperature 300
pressure 0.0001

equilbration 0.5 ps
production 0.5 ps
timestep 0.001 ps
sample 0.001 ps

shellmass
O 0.1

2

Hi Alex,
Shell model MD is a tricky (and more expensive business) and so you might want to try with a rigid ion force field first to make things easier (and faster). I can see that you are missing the 3-body potential from the Sanders force field and so this will cause you problems and lead to incorrect results as the model was fitted to have this.
There are a few other reasons that will cause the MD to heat up (and blow up!) quite often:

  1. If you put in an unrelaxed structure, particularly where the atoms are a long way from the bulk equilibrium structure, then the first thing that has to happen is to relax the structure. This will dump all the excess potential energy into kinetic energy and the thermostat may not be able to remove it fast enough since it’s not designed to handle this. So make sure you’ve partially optimised your structure before starting MD.
  2. You say you are running a nanoparticle, but you are running NPT. This means your unit cell will collapse and the particles will join together to make a rough 3D structure. Is this what you want? It will also definite put a lot of excess energy into the system (as above).
  3. Your input has 0.5 ps of equilibration before 0.5 ps of production. This is way too short to achieve equilibration from a high energy starting point and so I would expect things to blow up when production starts. This is because velocity rescaling is turned off at that point, which is the main way of trying to stop the velocities from getting too large or small.
    Hope that helps you.
3

Hi Julian,

Thank you for your response! It gave me a lot to think about. I think I’ve gone about the whole thing in an awkward and difficult way, which is causing all the problems in the md run. So I’ve taken a step back from original idea I had.

Again, thank you for the response and sorry for my late one!
Alex