Potential kinetic energy and temperature

[Hypatio]

I have a working non-equilibrium molecular dynamics model of MgO. I first find the equilibrium positions of the atoms in the lattice and then introduce a temperature by assigning initial velocities in random directions for each atom. An equation relating atomic velocity and temperature is

$mv^2=3k_BT$
where k_b is the Boltzmann constant.

So for example, oxygen with mass=2.66*10^-26 kg, and with T=2000K each oxygen is assigned a velocity v=1766 m/s.

Running the simulation, however, I observe that as the atoms move out of their T=0K equilibrium positions their velocity decreases until they begin bouncing off zero as the atoms bounce off one another. When the system is in a steady state of vibration the average temperature and atomic velocity in the lattice becomes about 1/2 of the initially assigned temperature, ~880 m/s and ~1000 K in this case.

This seems to be a transfer of kinetic energy into potential energy from moving out of the bottom of the potential energy well, but I'm not sure? Is this associated with a particular concept or phenomenon? Latent heat? The decrease in velocity is roughly 1/2 for all temperatures so I don't think it is related to specific heat. Where is the energy going? Is the actual temperature of the lattice 1000K or 2000K? How do I determine temperature in such a case?

 

[Useful nucleus]

Assigning one value for the velocity for all atoms of one type is not reasonable. You'd better create a velocity distribution (Maxwell-Boltzmann).
Regarding what happens later on in your simulation, you need to tell us more about the time step of your simulation. Is it short enough so as not to blow up the simulation. What ensemble are you simulating? If you're simulating a microcanonical ensemble, is the total energy conserved? etc...

 

[Hypatio]

Hi UN,

I would have created a velocity distribution, but it appears that from the initial condition of constant temperature an appropriate distribution quickly occurs in the model.

The timestep is lower than a point of instability, although I think it is blowing up for "high" temperature simulations (>2500 K). I am currently using 1.0 femtoseconds. Significantly decreasing the timestep doesn't result in any clear changes for lower temperature simulations but I will try it for higher temperature.

I believe my simulation is essentially microcanonical. The lattice is isolated and free to expand or contract from vibration, although I will later simulate temperature gradients to evaluate thermal conductivity. It is a 20x20 (400 atom) two-dimensional lattice. I am not certain that energy is conserved although there is no reason to doubt it unless my interpretation of the temperature change is incorrect. The only thermal energy lost should be that related to the momentum of the entire lattice, which is small.

I think that the average velocity of the lattice becomes 1/2 the initial velocity because the potential energy increases when atoms are displaced from their equilibrium positions. It is always 1/2 because the given velocity is the maximum and the minimum is zero at the point of elastic rebound against adjacent atoms. But what does this mean for temperature?

The total energy associated with each atom should be the sum of it's kinetic and potential energy, but is the potential energy thermal? Does it contribute to temperature?

 

[Timo]

Yes, the potential energy is "thermal". What you do is starting in a positional configuration corresponding to T=0 and a kinematic configuration corresponding to T=2000 K (or whatever value). When your system equilibrates, some of the kinetic energy goes into potential energy (-> equipartition theorem), resulting in a configuration representative for some temperature 0<T<2000 K. The most obvious solution to your problem seems to be coupling a thermostat to your system. That's a dirty trick from a purist's perspective (but few people a purists nowadays, as dirty tricks get you more publications than formally proper work), but as I understood it you ultimately don't want to work in a microcanonical ensemble, anyways.

 

[Hypatio]

So perhaps, due to equipartition, I am seeing a kinetic 'cooling' on the order of 1/2 because I have 2 degrees of freedom in my two-dimensional model? And the 2k_bT/3 law cannot describe the temperature of my system. Furthermore, this should be related to specific heat. But then why am I getting similar values (1/2 initial temperature) for all temperatures? Is there something in my model that prevents a reasonable approximation of the strong temperature-dependence of specific heat? Maybe I need a larger lattice in which phonons develop?

 

[Timo]

My guess (without having seen your simulation): You do not really see cooling. You see equilibration. Assuming the deviations from the ground state positions are not too small, then you can Taylor expand around those. This will give you potential energy terms that are harmonics springs. The number of degrees of freedom (in the sense of the equipartition theorem) for the potential energy therefore equals the number of degrees of freedom in the kinetic energy for this approximated Hamiltonian => half of the energy (that was initially only in the kinetic terms) goes into the potential energy. That's irrespective of temperature (at least for as long as the ansatz with the Taylor expansion holds), the spatial dimension of your system, its size, and the specific heat.

 

[Hypatio]

Well, how is this related to the thermodynamic temperature then? It seems that the equation

[itex]\frac{1}{2}mv^2=\frac{3}{2}k_b T[\itex]

cannot then be true for a particle because half the energy is converted into potential energy.