Relationship between nonlinearities and thermalization in classical open systems: The role of the interaction range

We discuss results on the dynamics of thermalization for a model with Gaussian interactions between two classical many-body systems trapped in external harmonic potentials. Previous work showed an approximate power-law scaling of the interaction energy with the number of particles, with particular focus on the dependence of the anomalous exponent on the interaction strength. Here we explore the role of the interaction range in determining anomalous exponents, showing that it is a more relevant parameter to differentiate distinct regimes of responses of the system. More specifically, on varying the interaction range from its largest values while keeping the interaction strength constant, we observe a crossover from an integrable system, approximating the Caldeira-Leggett interaction term in the long-range limit, to an intermediate interaction range in which the system manifests anomalous scaling, and finally to a regime of local interactions in which anomalous scaling disappears. A Fourier analysis of the interaction energy shows that nonlinearities give rise to an effective bath with a broad band of frequencies, even when starting with only two distinct trapping frequencies, yielding efficient thermalization in the intermediate regime of interaction range. We provide qualitative arguments, based on an analogous Fourier analysis of the standard map, supporting the view that anomalous scaling and features of the Fourier spectrum may be used as proxies to identify the role of chaotic dynamics. Our work, that encompasses models developed in different contexts and unifies them in a common framework, may be relevant to the general understanding of the role of nonlinearities in a variety of many-body classical systems, ranging from plasmas to trapped atoms and ions.

Figure 1

Scaling exponent $\alpha$ versus $\lambda$, evaluated, as described in Ref. [17], by fitting the interaction energy at long times to power law in the number of harmonic oscillators, $N_A=N_B=50,100,200,400$. Both short- and long-range interactions do not result in anomalous scaling for the exponent $\alpha$, while at intermediate values the scaling exponents are anomalous. The curves presented in each plot differ by the value of $\gamma$, showing a propensity to a broader and deeper anomalous region in the $\lambda$ space with increasing $\gamma$. Also, the two plots show a dependence on the relative values of the trapping frequency, (a) being obtained for ${\omega }_B/{\omega }_A=144/89$, (b) for the case of identical trapping frequencies, ${\omega }_B/{\omega }_A=1.$

Figure 2

Evolution of a particle in phase space (position in abscissa, momentum in ordinate) for various values of the interaction range $\lambda$, and the same initial conditions for all particles in the system, with interaction strength $\gamma =1$. The initial condition in phase space for the selected particle is $(0.13455,0.14265)$. The presence of rare interactions between particles is evidenced in (a) corresponding to a quasilocal interaction. Progressive increases by two orders of magnitude of $\lambda$ correspond to plots (b) and (c) in which more frequent interactions appears especially in the center of the trapping potential, which the latter situation corresponding to large energy exchange, while (d) shows a case of very large $\lambda$ for which nonlinearities in the interaction potentials are not enough to create energy exchange. For comparison, in the latter case we also report, in red (innermost ellipse, thin line), the trajectory corresponds to a pure harmonic motion of the particle, starting from the same initial conditions but without interparticle interactions (closed dynamics).

Figure 3

Dependence of the dynamics of thermalization on the interaction range $\lambda$ for two systems made of 500 particles each. (a) Shown are the curves of the dependence of the inverse temperatures versus time for the same initial inverse temperatures of 0.2 (dashed lines) and 2.0 (continuous lines) in arbitrary units, and same coupling strengths $\gamma =1$. The final value of the equilibrium inverse temperature, as well as the thermalization timescales, strongly depend on $\lambda$ in the long-range limit. (b) Plot of the final inverse temperatures, after ${10}^5$ time steps, versus the interaction range $\lambda$ for three different values of the interaction strength $\gamma$. This shows clearly that thermalization occurs, for the chosen timescale, only for small $\lambda$ and large $\gamma$.

Figure 4

Contour plots for the scaling exponent $\alpha$ [left plots, (a) and (c)] and the thermalization parameter $\overline{\epsilon }\left(t_f\right)$ [right plots, (b) and (d)] in the $\lambda \text{−}\gamma$ plane. The plots are taken after evolving the system for ${10}^5$ time steps (top plots) and for ${10}^6$ time steps (bottom plots). Notice the consistent behavior of the scaling exponent near 2 in the right-bottom part of the related plot, where the Caldeira-Leggett approximation is expected to hold. The initial inverse temperatures are 2.0 and 0.2, resulting in an initial temperature contrast of $\epsilon \left(0\right)=18/11\simeq 1.636$.

Figure 5

Fourier analysis of the thermalization dynamics for $N=400$ and $\gamma =1$. The three-dimensional cascade plots on the left of each panel show the evolution of the FFT spectrum of the interaction energy versus time. The FFT is constructed by considering time series samples of 2048 time steps. The plots on the right of each panel show the FFT at two specific time windows, the initial segment (first 2048 time steps, continuous black lines), and at the end of the full time period considered (last 2048 time steps, dashed red lines). The cases shown are for (a) $\lambda ={10}^{-3}$, (b) ${10}^{-2}$, (c) ${10}^{-1}$, (d) 1, (e) ${10}^2$ while, for contrast to the case (e), the Fourier spectrum of the interaction energy when the particles of each bath move harmonically with equal amplitudes at the frequencies of $2{\omega }_A/2\pi$ and $2{\omega }_B/2\pi$ is shown in (f). (e) corresponds to a case where the interaction term closely resembles the Caldeira-Leggett form, for which we do not expect the composite system to thermalize.

Figure 6

Comparison between the Fourier transforms at initial (continuous lines) and final (dashed lines) times in the two cases of (a) the interaction Hamiltonian based on Eq. (1), and (b) on its quadratic approximation as in Eq. (7) for $\lambda =10$ and values of $\gamma ={10}^{-4}$ (blue), $\gamma ={10}^{-2}$ (red), and $\gamma =1$ (black). In these plots we have considered time series of 8192 time steps in computing the FFT, resulting in a higher-frequency resolution with respect to the plots in Fig. 5.

Figure 7

Dependence of the scaling exponent $\delta$ upon the interaction range $\lambda$. For each value of $\gamma$ and $\lambda$, the FFTs are constructed using the same procedure as indicated in Fig. 5, using 20 sequences of length 2048 time steps at the end of the simulation. A range of frequencies is uniformly chosen from 1/8 to 3/8 of the Nyquist frequency and the FFT values in this range are fitted with a power-law function to extract the exponent $\delta$ for each sequence. The 20 values of $\delta$ are then analyzed to get both the mean values and error bars quoted in the figure. We have verified that varying the range of frequencies or number of segments considered does not affect the qualitative features of the plot of our analysis. Three regimes are clearly visible corresponding to a nonmonotonic behavior, with constant $\delta$ at low $\lambda$, then a sudden decrease until $\delta$ reaches values of about $-2$, followed by an increase at large $\lambda$ to a new plateau corresponding to $\delta \simeq -0.8$. The data also show the dependence for three distinct values of $\gamma$, most notably an increase of the value of $\lambda$, while increasing $\gamma$, for which $\delta$ reaches its minimum value.

Figure 8

Dependence of the Fourier scaling exponent $\delta$ for the simple standard map dynamics upon the stochasticity parameter $K$. For each $K$, a long time series (${10}^5$ time steps) of the momentum is generated for a set of initial 40 conditions. The ensemble averaged momentum is used for constructing the FFTs and an analogous procedure as adopted for generating Fig. 7 is followed to extract the exponent $\delta$ and its error bars for each sequence. These values are plotted against $1/K$ tto allow for a better inferential match with the case of the Gaussian interaction potential, since the role of $K$ in the latter system is played by $\gamma /{\lambda }^2$. Note that large $K$ corresponds to global chaos in the standard map while for values $1\lesssim K\lesssim 4$ stable regions coexist with chaos. The two curves labeled broad (circles) and narrow (diamonds) correspond to initial conditions drawn from the entire domain and from a limited (10% of the allowed domain) area centered at a randomly selected point in phase space. Beyond the threshold for global chaos, the behavior is independent of the choice.