15
$\begingroup$

I completed an ab initio molecular dynamic (AIMD) simulation in VASP for a hypothetical material. The simulation was done for 5ps with 1fs steps and the temperature was set to 300K using Nosé–Hoover thermostat. Now I've extracted the total energy variation w.r.t time and could not decide if the material is potentially stable or not from the plot. AIMD plot

Does this support that the material could be stable since the total energy converges around $\ce{-800 eV/cell}$ after 1ps? My concern is there is a huge variation in total energy during the 1st pico second, what does that mean? And what other tests should I do to investigate stability/instability of materials using AIMD simulations (e.g. variation of bond length over time)?

Improve this question
$\endgroup$
6
  • 2
    $\begingroup$ +1. Interesting question, but what type of stability are you looking for? Are you looking for dynamical stability (the structure of the material stays the same), chemical stability (segregation into different compounds), or something else? $\endgroup$
    – ProfM
    Jul 3, 2020 at 7:16
  • 1
    $\begingroup$ I am looking for both. Does the material maintain its structural integrity and does it not decompose to different compounds. What other tests should I do in addition to this energy curve to inspect this? $\endgroup$
    – Achintha Ihalage
    Jul 3, 2020 at 7:59
  • 2
    $\begingroup$ Did you relax the structure first, before starting MD? If you didn't, then what you see in the MD is just the structure relaxation effect. If you did, then the MD appers to have found a close-lying lower minimum. $\endgroup$
    – Susi Lehtola
    Jul 3, 2020 at 10:53
  • $\begingroup$ I did relax the structure first and used the optimised structure to start the MD simulation. By "close-lying lower minimum" do you mean the initial structure breaks at this temperature? $\endgroup$
    – Achintha Ihalage
    Jul 3, 2020 at 11:41
  • $\begingroup$ @AchinthaIhalage your energy lowers by 600 eV in the simulation, which is huge. How big is the system, i.e. what is the energy change per atom? $\endgroup$
    – Susi Lehtola
    Jul 4, 2020 at 9:29

1 Answer 1

Reset to default
13
$\begingroup$

This is a description of some methods available to study dynamical and chemical stability, not limited to MD:

Dynamical stability. The question here is: given a compound in a particular structure, is this structure dynamically stable or is there a lower-energy structure that can be accessed without jumping over an energy barrier? From a potential energy surface (PES) point of view, the question is: is this structure at a local minimum of the PES, or is it at a saddle point? The most straight-forward way of answering this question is to perform a phonon calculation. If all phonon frequencies are real, then the system is dynamically stable (at a local minimum), but if there are some imaginary frequencies, then the system may be dynamically unstable (at a saddle point). So my immediate suggestion would be to do a phonon calculation (which is also much cheaper than an MD simulation). If all frequencies are real, then you are done, the structure is dynamically stable. If you encouter imaginary frequencies, then the structure may or may not be dynamically stable. The structure definitely sits at a saddle point of the PES, but entropic contributions may help stabilize it at finite temperature (you are looking at the free energy surface rather than at the PES). A well-known example of this latter case is that of perovskites, which undergo a series of temperature-induced structural phase transitions. In this case, phonons are of no use (the relevant terms are anharmonic), and MD may come in handy. My suggestion would then be to start with the lower energy structure that the imaginary phonons take you to, and then run and MD simulation to see if the structure "symmetrizes" from the local minimum to the saddle point.

Chemical stability. The question here is: given a compound with a particular stoichiometry, will is decompose into other compounds? MD is not really a viable way to answer it in general; instead, you probably want to go down the route of structure prediction. What you would do is to generate multiple candidate structures at different stoichiometries, and then construct a Hull diagram to investigate if your structure will stay put or decompose. MD can simply not access the long time scales that are necessary to see something like a phase decomposition.

Improve this answer
$\endgroup$
5
  • 1
    $\begingroup$ Thanks for this answer! The calculation is indeed for a potential perovskite. I am trying to grasp what you meant by "jumping over an energy barrier". Does this mean, for a dynamically stable material, the energy should not deviate too much during a MD simulation, in contrary to my figure? I used a small supercell (2x2x2) for this simulation because of the computational power. $\endgroup$
    – Achintha Ihalage
    Jul 3, 2020 at 11:38
  • $\begingroup$ Secondly, after MD I noticed that the symmetry of the original structure lowers from Fm-3m to P1. I guess this basically means that the material is unstable? By the way I followed VASP MD documentation to get these results. $\endgroup$
    – Achintha Ihalage
    Jul 3, 2020 at 11:44
  • 2
    $\begingroup$ A structure is dynamically stable if it sits at a local minimum. When this is the case, it will have real phonons. However, it may be that there is a lower energy minimum elsewhere that you could reach by jumping an energy barrier. If this is the case, then phonons are of no use. MD may work, but it is very possible that the time scales are too long for MD. If this was the situation, I would be tempted again to go for structure prediction... $\endgroup$
    – ProfM
    Jul 3, 2020 at 13:46
  • 1
    $\begingroup$ If the symmetry lowers, then this indeed may indicate that you were not at a local minimum of the PES, and the MD is taking you there. However, if you go from Fm-3m to P1 (which essentially means no symmetry whatsoever), then my suspicion is that this may be a "sampling" issue with MD. You may have to be a bit more lenient with your symmetry tolerances to make sure that you have P1 symmetry rather than a higher symmetry which your MD ensemble is not large enough to capture. $\endgroup$
    – ProfM
    Jul 3, 2020 at 13:50
  • $\begingroup$ Hello @ProfM, great answer as always. About chemical stability: how would you create a hull diagram at a temperature above $0$ K? What if a particular phase or stoichiometry becomes stable at a higher T? $\endgroup$
    – Hitanshu Sachania
    Sep 27, 2021 at 19:55

You must log in to answer this question.

Not the answer you're looking for? Browse other questions tagged .