Simple Harmonic Oscillator Model of O 2 Molecule in Vacuum: A Classical Molecular Dynamics Study

The simple harmonic oscillator model allows a basic understanding of all processes and can be used to analyse optical vibrational modes and electronic transitions in atoms, molecules and crystals, in order to derive general properties of harmonic generation to all orders. In particular, we are done to investigate single, ten, and twenty harmonic oscillator of O 2 molecule using GROMACS 4.55. Time is the important factor in simulation because the best result of simulation can be obtained by increasing the number of steps and by decreasing the timesteps. Since the properties only depend on time and not on the specific microscopic model, they can also be adopted for the quantum-mechanical description by using in the classical molecular dynamics.


Introduction
Computation and simulation are top role in modern science.They are bridge between theoretical and experimental methods to obtain data in science.The harmonic oscillator is extremely useful in chemistry as a model for the vibrational motion in a diatomic molecule (Guillaume L, 2003), (J.J. Rael et al, 2000) and (Charles B et al, 1990).The atoms are viewed as point masses which are connected by bonds which act approximately like springs obeying Hooke's law (Ugural, et al, 2003), (Keith., 1971) and (Simo, J. C, 1998).To get started with it, a simple model of bond vibration of diatomic molecule O 2 can be done.They are many approaches to investigate bond vibration (F. lachello, 1998), (Philip M, 2002) and (O. S. van Roosmalen., 1998).
One of them is harmonic oscillator.Here, oxygen atom is a type to describe the atoms in the diatomic molecule.Two oxygen atoms are treated as a classical object which are connected by a spring.The ilustration and force formulation of Hooke's law of harmonic oscillator which working in the spring is in the equation 1 and Figure 1.
IJCR-Indonesian Journal of Chemical Research p. ISSN: 2354-9610 e. ISSN: 2614-5081 Vol.3 No.1, Hal. 8-16 9 Where F is the resulting force, x is the displacement of the mass from equilibrium (x = rr q ) and k is the force constant.

Molecular dynamics simulation is based in
Newton's equations of motion (N.Goga et al, 2015).The simulations use classical mechanical descriptions of energy as a motion.

Single harmonic oscillator of O 2 molecule in vacuum
The energy can be regarded as vibrational energy only because rotational and translational motions do not occur.The data of potential energy, kinetic energy, total energy, average distance, and bond length of simulation are shown in Table 1.The motion motion O 2 molecule is periodic.Because the given force constant is quite high, and the period is quite small.The other analysis with different time (2-6 ps) simulation is showed in Table 2 and    4. It shows that the high force constant is also high in the average and tot-Drift energy.
The effect of mass can be seen in Table 5.
The mass affect the frequency but it does not have effect on the total energy nor on the drift.It shows that they differ in the magnitute of the probability.The distribution of large system also shows that the peak of single oscillator is not the same with ten and twenty oscillator.

Figure 1 .
Figure 1.Ilustration of a simple harmonic oscillator single harmonic oscillator simulation, one molecule O 2 was put in a vacuum 3D box.The mass of an atom (m O ) was 16 amu.The distance at which the potential was minimum (b 0 ) was set to be 0.11 nm, the force constant that contribute to the vibration (k) was 500 kJ.mol -1 .nm -2 In GROMACS, nuclei and electron are usually treated as an entities which is described by lennard-Jones potential as follows in equation 2, U (r) = 4 ԑ [ (σ/r) 12 -(σ/r) 6 ](2)ԑ is depth of potential well, σ is the finite distance at which the inter-particle potential is zero, and r is distance between two particles.(C6 1.0e-6, C12=1.0e-10).Molecular dynamics simulation will be performed with 10000 steps which 0.001 ps timesteps.The center of mass is not fixed in a certain coordinate to study the effect of vibration to translational motion.Periodic boundary condition is not applied.If it will be applied, the size of 3D box has to be increased to make sure that there is no interaction among O 2 molecule.

Force
constant is proportional to the total energy.The force constant is not the only parameter that affect the period.The other one is the mass of the atoms involved in the motion.The larger force constant applied, the larger total energy is showed in O 2 molecule.The investigations are continued by the modification of timesteps to know the maximum timestep (algorithm is stable).The timesteps of 0.01ps is good to represent the energy and bond length average but the timesteps of 0.2 ps is bad to represent them because the average bond lenght differs and the standard deviation is high.

Figure 2 .
Figure 2. Potential, kinetic and total energy as a function of time (a) and at 2-6 ps (b)

Figure 3 .
Figure 3. (a) Average distance of single harmonic oscillator of O 2 molecule, and (b) Distribution of harmonic bond distance as a function of bond length

Figure
Figure 2 (a).In order to improve energy conservation by reducing time, the timestep try to be reduced but it makes a problem which is higher computational cost.In this experiment 0.001 ps is enough to represent the system and conservation energy that can be maintained.The result of reducing timestep can be seen in Table3.Natural frequency is a function of force constant and mass of atom.Many variations of natural frequency can be made by manipulating force constant and mass atom.

Figure 4 .
Figure 4. (a) Total energy as a function of time and (b) at velocity=0 the graph does not look like Boltzman distribution, It is only vibration can occur and the atom spends more time in maximum and minimum bond length.The initial atomic coordinate play an important role because it determines the amplitude and energy of vibration but the frequency of vibration still the same.Weakly coupled harmonic oscillator of O 2 molecule in vacuum The distribution of harmonic bond length of ten and twenty harmonic oscillators at constant energy are given in Figure 6 (a) and (b).The distribution is actually the same between ten and twenty oscillators.Comparing Figure 6 (a) and (b) to distribution at single oscillator (Figure 2 b).

Figure 6 .
Figure 6.Distribution of ten (a) and twenty (b) harmonic oscilators as a function of bond lenght (nm)

Figure 7 .
Figure 7. Distribution of twenty oscillators at constant temperature as a function of bond lenght (nm)

Figure 8 .
Figure 8. Distribution at 10 oscillator constant energy (black), 20 oscillator constant energy (red), 20 oscillator at constant temperature (green) as a function of bond distance.

Table 1 .
Calculation result of single harmonic oscillator O 2 molecule (time =10 ps)

Table 3 .
Total energy, average, error estimate, RMSD and Tot-

Table 4 .
Total energy and total energy drift in different force constant Variational of initial bond length ro change the total energy.The |ro-bo| considered as the amplitude of vibration.

Table 5 .
Total energy and total energy drift in different mass