A molecular dynamics simulation of thermalization of crystalline lattice with harmonic interaction

Abstract

Understanding the realization of thermal equilibrium through the thermalization process in a many-body system is a fundamental and complex scientific question, bridging thermodynamics and classical dynamics and connecting to a host of physical phenomena, such as mechanical instabilities in a thermal environment. In this work, based on the harmonic lattice model, we investigate the thermalization process in both velocity and coordinate spaces, by examining microscopic dynamics on the atomic level. We show the distinct relaxation rates of the transverse and longitudinal components of the velocity, reveal the power law governing the nonlinear proliferation of dominant frequencies, and observe the concurrent rapid proliferations of frequencies and topological defects. We also show that the lattice system's persistent out-of-plane deformations exhibit two-stage fluctuation behaviors, characterized by distinct power laws of fractional exponents and associated with the broken up-down symmetry. This work demonstrates the rich dynamics underlying the thermalization process, and advances our understanding on the dynamical adaptations of many-body systems to external disturbances.

Figures (7)

• Figure 1: Illustration of the drum model with the anchored boundary (in orange). (a) The drum model consists of the triangular lattice of linear springs (bonds) and point particles (vertices). The five- and seven-fold disclinations are colored in red and blue, respectively. Disclinations are identified by the Delaunay triangulation of particle configuration projected onto the plane. Note that the geometric bonds established by the Delaunay triangulation algorithm (not shown in the figure) can be distinct from the physical bonds (springs) connecting particles. (b) The contour plot shows the instantaneous out-of-plane deformation in the dynamical evolution of the lattice. $t=5\tau_0$. $v_0=0.5$. $N=316$.
• Figure 2: Characterization of the thermalization process by the relaxation of velocities. (a) and (b) Plots of $\delta f_{\parallel}$ and $\delta f_{\perp}$ (the deviations of $v_{\parallel}$- and $v_{\perp}$-distributions from the Maxwell distributions) versus time at typical values of $v_0$. (c) Distinct relaxation time of $v_{\parallel}$ and $v_{\perp}$ versus $v_0$. (d) The log-log plot of temperature (as derived from the Maxwell distributions) versus $v_0$ show that $T\propto v_0^2$ for both cases of $v_{\parallel}$ and $v_{\perp}$. The sampling is during $t\in [900\tau_0, 1000\tau_0]$ at the time resolution $\Delta t=10\tau_0$. The standard deviations are indicated by the error bars. $N=1406$.
• Figure 3: Thermalization process from the perspective of the nonlinear proliferation of dominant frequencies. (a) Plot of the kinetic energy versus time. (b) The frequency spectra of kinetic energy for the observation time of $t=10\tau_0$ and $t=20\tau_0$. (c) The number of dominant frequencies increases from one (green), two (blue) to three (red) with the extension of the observation time $t$. The empty square and triangle symbols show the sum of the lower frequency in both two- and three-frequency regimes and a multiple of $\omega_{\rm{min}}$ (the lowest curve in gray). (d) The growth of the total number of frequencies [$\Omega(t)$] in the early stage conforms to the power law. (e) and (f) Dependence of the exponent $\alpha$ in the power law of $\Omega(t)$ on $v_0$ for the cases of 2D lattices allowing and restricting out-of-plane deformations, respectively. The arrows indicate the critical value of $v_0$, above which the $\alpha(v_0)$ curve experiences a rapid rise. The statistical analysis is based on 10 independent initial distributions of velocity. The sampling is during $t\in [\tau_0, 100\tau_0]$ at the time resolution $\Delta t$. $\Delta t=0.001\tau_0$ in (c); otherwise, $\Delta t=0.01\tau_0$. $v_0=0.1$ in (a)-(c). $N=847$ in (a)-(d).
• Figure 4: Illustration of the discrete theoretical model for understanding the nonlinear nature of frequency proliferation processes. (a) An example for showing the growth of the number of frequencies in time. The initial state is characterized by the two occupied frequencies (orange dots). The evolution of the state follows the prescribed rules; more information is provided in the text. (b) Log-log plot of the total number of frequencies $\Omega$ versus time at typical values of $n_{\rm{max}}$ and $\gamma$. $n_{\rm{max}}$ is the total number of permissible frequencies, and $\gamma$ is the proportion of the number of frequencies involved in the frequency-proliferation process in each time step. The slopes of the fitting lines are indicated at the typical sites, showing that the discrete model accommodates a range of slope values covering both the nonlinear regimes below and above the linear law. The statistical analysis is based on 20 independent simulations; the two initial frequencies are randomly specified in each simulation. The standard deviations, which are similar for both cases of $n_{\rm{max}}=10,000$ and $n_{\rm{max}}=5,000$, are indicated by the error bars only for the case of $n_{\rm{max}}=10,000$ for the sake of visual clarity.
• Figure 5: Characterization of the in-plane and out-of-plane fluctuations. (a) An instantaneous particle configuration containing topological defects; $v_0=0.512$. The red and blue dots represent five- and seven-fold disclinations. The insets provide a zoomed-in view of the quadrupoles highlighted in green. (b) $\langle N_d\rangle /N$ is the time-averaged total number of disclinations, rescaled by the total number of particles. The data for both cases, with and without out-of-plane deformations, collapse on the same curve at typical system sizes. The arrow indicates the critical value of $v_0$, above which rapid proliferation of topological defects occurs. (c) Plot of the logarithm of $\langle h^2\rangle$ (the mean squared transverse displacement) versus $v_0$. The values of the slopes $k_1$ and $k_2$ of the fitting lines are shown. The cross point of the two fitting lines is located at $v_0 \approx 0.025$. $N=847$. (d) Plot of $|\langle h_{+}\rangle - \langle h_{-}\rangle |$ versus $v_0$ shows the breaking of the up-down symmetry with the increase of $v_0$. $\langle h_{+}\rangle$ and $\langle h_{-}\rangle$ are the magnitude of the transverse deformation along the $z$ and -$z$ directions, respectively. The critical value of $v_0$, above which the up-down symmetry is appreciably broken, coincides with the cross point in (c). The statistical analysis in (b)-(d) is based on the sampling during $t\in [\tau_0, 100\tau_0]$ at the resolution of $\Delta t=2\tau_0$. The standard deviations are indicated by the error bars.