I agree with Stéphane's opinion as well. Make sure that the system is well equilibrated at 300K then gradually increase the temp, probably with long periods of holding the temp steady at intermediate stages. Once you have a stable NVT ensemble at 600K then switch to to NVE.

Thanks: truly I think the answers are good. As Stéphane says as an NVE there is no thermostat. During the NVE, I use a first stage of equilibration in which I controlled the temperature with just a 'scale' directive: answering Raiker I'm not sure if this is a particular thermostat: I have not found data about this. Finally Stéphane's idea (also thanks to Pavel and David) seems quite good: my equilibration at 600 K started from the minimised structure, but I didn't make a previous equilibration at 300 K. What I -instead- tried is to increase more and more the number of steps of the equilibration stage. Perhaps this option of making a two-stage equilibration (300 K, 600 K) through a NVT ensemble will give a better result.

I have no experience in solid state Molecular Dynamics simulations, but it appears obvious to me that the temperature is not constant if you don't apply a temperature coupling to the system.

You are simulating a microcanonical ensemble instead of a canonical one.

Does the system include only the zeolite or also gas or solvent molecules are there?

If there are some special reason for not using temperature coupling for the production run, I would suggest to slowly equilibrate to the desire temperature (using a T_coupling algorithm) and once the systems appears to be stable switch off the temperature coupling. In the equilibration-phase I would also use pressure coupling, because an increase/decrease of temperature would cause an increase/decreas of volume of your system.

Thanks Andrea. Prior to the MD simulation I have made a molecular mechanics study where I have obtained accurate values for the unit cell size at the 3 different temperatures I'm studying: 300, 450 and 600 K. This shows, by the way, the usual negative thermal expansion of many zeolites, and hence the cell volume slightly shrinks as the temperature increase. With the obtained unit cell size, I fix it in the NVE MD and I add a certain loading of hydrocarbon molecules inside. I prefer to have the energy constant and check that the temperature remains fix (i.e. use the microcanonical ensemble) rather than using the canonical ensemble (NVT) where the temp. is fixed and I should check that the energy of the system remains constant.

I'm following your suggestion of slow equilibration, in relation with other answers. If not, I may try the NVT option.

Hi,

I agree with Abel and Reichert answers, But I'd add that if the framework is fixed, it is difficult to maintain a stable temperature, because the partile momentum is not conserved. The reason why the proble appears at 600 K and not at 300 K should probably be looked for in the dynamics of the adsorbed molecules. It is likely that at 300 K the motion of the molecules is alost confined at the adsorbing sites, while at 600 K their mobility is sufficiently high to make the collisions with the framework more frequent. It is this kind of collisions which do not conderve the particle momentum. Therefore, in my opinion, if slow gradual heating does not solve the problem, a thermostat must be applied. It can be provided by the flexible framework or by a Nosé-hoover or Berendsen thermostat with suitable parameters.

Caro Giuseppe, grazie per la risposta.

I probably forgot to specify that all my MD are with fully flexible zeolite and hydrocarbon. Currently I'm making the simulations that you and the others have so kindly suggested and before the weekend I hopefully will report here an update.

So, finally it has worked! It may have been a combination of fixes, but mainly I would say two: (1) Stepwise increase of temperature as separate NVT runs until a smooth thermalisation at my desired temperature (600 K) has been obtained; and (2) Using a better optimised unit cell: initially I had considered an estimation of the unit cell experimentally obtained taking into account the negative thermal expansion: yes, very accurate, BUT not predicted by the force field in use. Hence I was using a unit cell size not sufficiently close to the minimum, and I'm sure everybody understand it is not wise to make a long MD with a too constrained unit cell: instabilities start to appear and some bond distances/angles, etc may start to depart from reasonable conformations, all of these resulting in energy drifts, etc.

I have to say that the temperature is not perfectly constant with the current setting: it goes from 600 to 625 K in 10 nanoseconds. But, considering the difficulty of the setup I am satisfied. Further, the test was made with the smallest number of hydrocarbon molecules inside the zeolite: I hope with more molecules temperature will remain even more close to 600 K.

Again, thanks to everybody! :-)

Hi German - one additional thing that I didn't see appear yet in the answers above, is that at 600 K, you will have adsorption and desorption of your hydrocarbons much more frequently. These processes are exothermic/endothermic, so each time this occurs, heat is released/absorbed, and therefore your temperature cannot remain constant - which is perfectly correct. Besides, if you run your simulation for a very long time, you will average out over all the adsorption/desorption events (if in dynamic equilibrium), and you will have some average temperature (i.e., constant over long time scales, not over short time scales).

Hi everybody. I agree that single adsorption/desorption event cause a transfer between kinetic and potential energy, changing instantaneous temperature, but adsorbed molecules are termostatted by the heath bath ensured by the flexible framework. We demostrated more than twenty years ago that the adsorbed molecules subsystem (J. Phys.Chem.96, 1482 (1992)) acts ad a small canonical enesemble, whose instantanoeous temperature fluctuaitions are LARGER than those of a NVE enesemble, as expected on the basis of statistivcal mechanics theory.(Lebowitz et al, Phys Rev. 153, 250 (1967)). Avrerage temperature, as remarked by Eric, usually remains constant, in the condidered case, if adsorption/desorption event numebrs is in turn constant, on avreage, over sufficiently long times.

I would also add another paper by Giuseppe also on the topic :

"A comment on the flexibility of framework in molecular dynamics

simulations of zeolites"

Microporous and Mesoporous Materials 2009, 125, 160–168

doi: 10.1016/j.micromeso.2009.03.032

In the NVE is expected that the temperature will fluctuate. Temperature depends on the velocities (cell size may not matter that much at the initial phase) , even when you initially adjust the velocities, so that the initial T is 600 K, this is just the instantaneous T, and T is a macroscopic properties and as such only its average value is meaningful. So, if you don't have a way how to keep close to the required T, it will fluctuate or even drift from the initial value. This is more problematic for high temperatures and especially heterogeneous system.

You should try some rescaling of the velocities until the system is equilibrated around 600 K, then turn it off, or the best way is to use thermostat, NH could be a good choice. You may also try to run long run NVT at 600 K until system is equilibrated, then switch to NVE (if for some reason you like the production run in NVE).

Hello everyone! I am also interested in your question because I meet the same problem now. I run my system in 300K without any problem but when I increase it up to 473K, time step 2fs and run within 20ns. It dies around 16ns. I have stable temperature by NVT with 5ns first and run NVE later. The point is if I only run 10ns, it did not show any error but my molecules diffuse not well. So do I need to add more time simulation in NVT to make sure that my system is heat slowly? Thank you very much!