Key role of asymmetric interactions in low-dimensional heat transport

How a fluctuation relaxes in the equilibrium state is very important, not only in its own right, but also because it governs the transport behavior of a system in a nonequilibrium state. Relaxation of a fluctuation is characterized by the corresponding correlation function, and according to the linear Boltzmann and Enskog equations, if a system is not in a critical region, in general its correlation functions are believed to decay exponentially at long times [1, 2]. Consequently, the integral of a flux autocorrelation function, which appears in the Green–Kubo formula, converges and thus guarantees a size-independent transport coefficient [3]. Until the late 1960s, the viewpoint that correlation functions decay exponentially had been prevailing. However, after Alder and Wainwright numerically evidenced the long-time tail of velocity autocorrelation in gas models in 1970 [4], extensive analytical and numerical studies suggested that in one-dimensional (1D) and two-dimensional (2D) momentum-conserving systems, the heat current autocorrelation function (HCAF) generally decays in power-law manners instead. (See [5–7] for reviews and more literatures.) The HCAF is defined as

$$\begin{eqnarray}&&C(t)\equiv \langle J(0)J(t)\rangle,\end{eqnarray} \tag{ 1 }$$

where $\langle \centerdot \rangle$ denotes the equilibrium thermodynamic average and J(t) is the total heat current. For 1D momentum-conserving systems, a recent analytical study has summarized that $C(t)\sim {{t}^{-\gamma}}$ with $\gamma =1/2$ and $\gamma =2/3$, respectively, for systems with symmetric and asymmetric interparticle interactions [8]. For 2D momentum-conserving systems [6, 9], the decaying exponent is $\gamma =1$. An important consequence of power-law decay is divergence of the heat conductivity for $\gamma \leqslant 1$ following the Green–Kubo formula [3, 6]

$$\begin{eqnarray}&&\kappa =\frac{1}{{{k}_{\text{B}}}{{T}^{2}}\text{d}}\underset{\tau \to \infty}{{\lim}}\,\underset{V\to \infty}{{\lim}}\,\frac{1}{V}{\int}_{0}^{\tau}C(t)\text{d}t,\end{eqnarray} \tag{ 2 }$$

where κ, T, d, and V are, respectively, the heat conductivity, the temperature, the dimension, and the volume of the system, and ${{k}_{\text{B}}}$ is the Boltzmann constant. It implies that low dimensional (i.e. 1D or 2D) materials, such as nanowires and graphene flakes, may possess an anomalous thermal transport property. At present, this viewpoint—that the HCAF of low dimensional momentum-conserving systems generally decays in power-law manners—is a mainstream viewpoint that has also been accepted by experimentalists.

Nevertheless, in a recent numerical study, it has been found that in 1D momentum-conserving lattices with asymmetric interparticle interactions, the heat conductivity may turn out to be independent of the system size [10]. This finding implies that in such a system, the HCAF may decay faster than power-law manners, or in a power-law manner but with $\gamma >1$. Later, faster decay has also been observed in several other 1D momentum-conserving lattices with asymmetric interactions [11–13]. The faster decaying behavior is in clear contrast with existing theories and has significant importance. On one hand, it requires us to revisit the theories developed for low dimensional transport problem during last decades; On the other hand, it implies that realistic low dimensional materials, such as nanowires and graphene flakes, may still follow the well-known Fourier heat conduction law [14–16], because real materials usually show the thermal expansion effect, which is a consequence of asymmetric interparticle interactions. This problem hence deserves careful studies.

The purpose of this paper is to explore if there exists any link between asymmetric interparticle interactions and the decaying behavior of the HCAF in 1D momentum-conserving lattices. The model we focus on is the Lennard-Jones (L-J) lattices, whose interaction asymmetry degree depends on both the temperature and a pair of system parameters. We show that in wide ranges of temperature and system parameters, this model has a high asymmetry degree and the HCAF decays faster than any power-law manners in the time range we can explore numerically. Nonequilibrium simulations also confirm that it has a size-independent heat conductivity in the corresponding system size range. The Fermi–Pasta–Ulam-α–β (FPU-α–β) model is also studied. This model has asymmetric interactions but faster decay of its HCAF has not been reported yet [13, 17–19]. Our analysis shows that in clear contrast to intuition, the effects of asymmetric interactions, such as thermal expansion, do not vanish in the low temperature regime; rather, it becomes even stronger in some sense. In this regime, fast decay of the HCAF is observed as well. To make a comparison, we also study the HCAF of the Fermi–Pasta–Ulam-β (FPU-β) model of symmetric interactions in the low temperature regime. It is confirmed again that under proper conditions, asymmetric interactions are closely related to fast decay of the HCAF. These evidences suggest that certain link exist between fast decay of the HCAF and asymmetric interactions, and this finding could be a helpful clue for revealing the mechanism(s) of the former. As the two lattice models of asymmetric interactions (the L-J model and the FPU-α–β model) are typical, we suspect such link may also be found in other asymmetric lattices with proper parameter values. Nevertheless, we have to emphasize that at present exact conditions for observing such link are not clear yet, and it would be risky to interpret it as a cause-effect relationship. These problems deserve more efforts and should be clarified in future.

We adopt molecular dynamics studies for our aim here. At present numerical analysis plays an important role in studying heat conduction properties of low dimensional systems. Although the problem can be dealt with analytically, certain approximations and assumptions have to be resorted to. For example, in analytical studies, it has been generally assumed that all slow variables of relevance for the long-time behavior of hydrodynamics and the related time correlation functions are the long-wavelength Fourier components of the conserved quantities' densities [8]. For this reason, particular attention must be paid when analytical results are compared with simulation and experiment results [17]. As to experimental studies, despite the fact that in recent years it has become technically feasible to measure the heat conductivity of low dimensional materials, the accessible precision is still far away for drawing convincing conclusions. In addition, realistic materials studied in laboratories, such as nanowires and graphene flakes, may not be genuine 1D and 2D objects, in view of their possible transverse motions. Therefore, how to control or evaluate the effects of the transverse motions turns out to be a new experimental challenge. In contrast, the molecular dynamics method does not suffer from any of these problems.

To perform a molecular dynamics simulation, one has to be satisfied with a finite system and evolve it for a finite time. In this work, we consider the system size (i.e. the number of particles) from $N\sim {{10}^{3}}$ to 10⁵ and the time range of the HCAF up to a time scale ${{t}_{\text{num}}}$. For $t<{{t}_{\text{num}}}$, our simulation results of the HCAF are reliable, in the sense that their errors are negligible and they are free of the finite size effect. Our computing algorithms and resources allow us to make ${{t}_{\text{num}}}\sim {{10}^{3}}$ for the L-J lattices and ${{t}_{\text{num}}}\sim {{10}^{4}}$ for the FPU-α–β lattices. To make the results accurate enough, large ensemble average is taken, and to make sure they do not suffer from the finite size effect, a careful analysis following [20] is done. To be precise, our numerical results of the HCAF are accurate for $t<{{t}_{\text{num}}}$ and $N>4\times {{10}^{3}}$. Our aim in this work is trying to find in this time range (i.e. $t<{{t}_{\text{num}}}$), what characteristics the HCAF may have, their possible relation to the interaction symmetry, and their implications to the heat conduction properties. One may wonder whether these characteristics may sustain to longer times, and whether the power-law tail of  ∼  t^−2/3 predicted by the analytical studies [8, 21, 22] in the thermodynamic limit will take over. This is interesting and important, and deserves further studies in future, probably with more powerful computing resources. In this work we do not attempt to address this question.

In the following we shall first specify the simulation setup of our numerical study. Then, with a measure of asymmetry degree introduced, the lattice model of Lennard-Jones potential will be analyzed. After that we shall turn to the asymmetric FPU-α–β model, and discuss to what a large system size the Fourier heat transport can be expected. A brief summary and discussions on some related problems will be presented in the last section.

We consider 1D lattices with the nearest neighboring coupling whose Hamiltonian can be written as

$$\begin{eqnarray}&&H=\underset{i}{\sum}\,\frac{p_{i}^{2}}{2{{m}_{i}}}+U\left({{x}_{i}}-{{x}_{i-1}}-a\right),\end{eqnarray} \tag{ 3 }$$

where p_i, x_i, m_i, and U represent, respectively, the momentum, the position, the mass of the ith particle and the potential between two neighboring particles. For all the models we study in the following, we assume that all the component particles are identical and have unit mass; i.e. m_i  =  1. The lattice constant, a, is set to be unity (a  =  1) so that the system length L equals the particle number N. Without confusion, in the following we shall use N to indicate the system length as well.

In our simulations for the HCAF, the periodic boundary conditions are imposed, and the total momentum of the system is set to be zero. Note that in 1D lattices, if there is no steady motion (i.e. the total momentum is zero), then the heat current equals the energy current [6]. We consider the total heat current defined as $J\equiv {\sum}_{i}{{\overset{\centerdot}{{x}}\,}_{i}}\frac{\partial U}{\partial {{x}_{i}}}$ [23]. To numerically measure the heat current in the equilibrium state, the system is first evolved from an appropriately assigned random initial condition for a long enough time (>10⁸ in our simulations) to ensure that it has relaxed to the equilibrium state, then the total heat current at ensuing times is measured. The energy density, or the energy per particle, denoted by $\langle E\rangle$, is determined by the initial condition and is conserved during the simulation. At the low temperature regime, by using the Viral theorem, we have ${{k}_{\text{B}}}T\approx \langle E\rangle \approx 2\langle U\rangle$, where $\langle U\rangle$ is the averaged potential energy per particle, or the average potential energy between two neighboring particles, at the equilibrium state. Note that this relation holds only when the harmonic term dominates the potential; hence only applies at the low temperature regime. The Boltzmann constant, ${{k}_{\text{B}}}$ is set to be unity (${{k}_{\text{B}}}=1$) throughout.

In order to check the heat conductivity computed based the Green–Kubo formula and the HCAF, denoted by ${{\kappa}_{\text{GK}}}$, we shall compute the heat conductivity directly with a nonequilibrium setup as well. In doing so we impose the fixed boundary conditions; Namely, we introduce two auxiliary particles, the zeroth and the (N  +  1)st, connect them to the first and the Nth particle, and fix their positions at x  =  0 and x  =  N, respectively. Two Langevin heat bathes [7] with two different temperatures ${{T}_{+}}=T+\Delta T$ and ${{T}_{-}}=T-\Delta T$ are applied to the first and the Nth particle. The system is evolved then, and after the stationary state is established, the heat conductivity is evaluated by assuming the Fourier law $j=-{{\kappa}_{\text{neq}}}\text{d}T/\text{d}x$, where $j\equiv \langle J\rangle /N$ is the heat current and $\text{d}T/\text{d}x\approx \left({{T}_{+}}-{{T}_{-}}\right)/N$ is the temperature gradient along the lattice. The former can be measured directly, and the latter can be set by specifying $\Delta T$. ${{\kappa}_{\text{neq}}}$ denotes the heat conductivity computed with this nonequilibrium setup.

The L-J potential has been widely adopted in modeling realistic materials. It is asymmetric with respect to the equilibrium point, and our study has shown that in 1D lattices with L-J potentials, the HCAF can decay faster than power-law manners [11]. This finding is in clear contrast to the well known theoretical prediction of the power-law decay in the thermodynamic limit.

The L-J potential has the form of $U(x)=\left[{{(x+1)}^{-m}}-2{{(x+1)}^{-n}}+1\right]$. It involves a pair of parameters, m and n, that control the asymmetry degree of the potential. Without loss of generality, we fix m  =  2n so that the minimum of U(x) is fixed at x  =  0. The potential is asymmetric with respect to the equilibrium point (see figure 1(a)). In order to compare the asymmetry degree for different (m, n) and at different temperatures, we introduce the following measure of asymmetry degree:

$$\begin{eqnarray}&&\eta \equiv \frac{\text{d}}{\text{d}\langle U\rangle}\left({{x}_{+}}-|{{x}_{-}}|\right),\end{eqnarray} \tag{ 4 }$$

where x₊ and x₋ are, respectively, the right- and the left-side zero point of $U(x)-\langle U\rangle =0$. Note that ${{x}_{+}}-|{{x}_{-}}|$ represents expansion and η is equivalent to the thermal expansion coefficient upon a factor of the heat capacity [24]. As η captures and reflects thermal effects of asymmetric interactions, it is a natural and physically meaningful choice to measure the asymmetry degree. As figure 1(b) shows, the asymmetry degree of the L-J potential increases as the average potential energy $\langle U\rangle$, and for a fixed $\langle U\rangle$ value it increases as the parameter pair (m, n) varies from (m, n)  =  (12, 6) to (2, 1). We shall focus on the case of (m, n)  =  (12, 6), the most frequently adopted parameters in literatures, but also discuss other values of (m, n) when it is in order. In figure 2 we show the temperature and the sound speed, v_s, of the L-J lattices at various values of $\langle E\rangle$, obtained by using the fluctuating hydrodynamics theory. The simulation results of the HCAF are shown in figure 3.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

First of all, we have to perform a finite size effect analysis to make sure that the simulation results of the HCAF are accurate and free of the finite size effect in the time range $t<{{t}_{\text{num}}}$. In doing so we have several options of the analysis method developed recently [20, 25, 26], but because those provided in [25, 26] are in the frequency representation, we adopt the one in the time representation [20] that severs our aim here more conveniently. We find that in 1D L-J lattices, the finite size effect is weak, which is a very favorable property. To show this, in figures 3(a) and (b) the HCAF at different system sizes are compared for (m, n)  =  (12, 6) and (m, n)  =  (2, 1), respectively. The energy density $\langle E\rangle$ is fixed at $\langle E\rangle =0.5$, which corresponds to a temperature of $T\approx 0.55$ and an average potential energy per particle $\langle U\rangle \approx 0.23$. It can be seen that in both cases, all the curves of C(t) collapse onto one perfectly upon scaling by the system size N. This observation explicitly suggests that the finite size effect is negligible for t  <  10³ in the case of (m, n)  =  (12, 6) and for $t<5\times {{10}^{3}}$ in the case of (m, n)  =  (2, 1), given that $N>4\times {{10}^{3}}$. Or, equivalently, these results suggest that ${{t}_{\text{num}}}\approx {{10}^{3}}$ in the former case and ${{t}_{\text{num}}}\approx 5\times {{10}^{3}}$ in the latter. Note that the oscillating tails around zero for t  >  10³ appear in figure 3(a) (also in figures 3(c) and (d)) are due to uncertainty of statistical average of J(0)J(t) (see equation (1)), and the negative parts of the oscillations are not shown as the vertical axis is in the logarithmic scale.

It shows clearly in figure 3(a) that the HCAF decays faster than any power-law manners for $t<{{t}_{\text{num}}}$. The inset presents the data in the semi-log scale, from which one can see that the decay manner is already very close to an exponential one. Figure 3(b) shows the results for (m, n)  =  (2, 1), and C(t) curves can be regarded as decaying exponentially with certainty. It can be found from figure 1(b) that the L-J potential with (m, n)  =  (2, 1) has a relatively higher asymmetry degree; We thus conjecture the higher the asymmetry degree is, the closer the decaying behavior of the HCAF tends to be exponential. To check this conjecture, in figure 3(c) we compare the HCAF for various (m, n) values. It shows that, as (m, n) varies from (12, 6) to (2, 1), i.e. as the asymmetry degree increases (see figure 1(b)), the C(t) curve in the semi-log scale becomes straighter and straighter (see the inset of figure 3(c)), suggesting that the decaying behavior indeed tends to be exponential progressively.

As figure 1(b) shows, for a given parameter pair (m, n), the asymmetry degree is controlled by the average potential energy per particle $\langle U\rangle$: it increases dramatically as the latter. As the temperature monotonically increases with $\langle U\rangle$, it implies that the temperature could also be correlated to the decaying behavior of the HCAF; i.e. we can expect that for the L-J model at a higher temperature, the average potential energy per particle is higher and consequently, potential's asymmetry degree would be higher (see figure 1(b)). In figure 3(d), temperature dependence of the decaying behavior of the HCAF is studied for (m, n)  =  (12, 6) with $\langle E\rangle =0.5$, 1, and 2. For these $\langle E\rangle$ values, the corresponding temperature is $T\approx 0.55$, 1.2, and 2.4 (see figure 2(a)), and the corresponding $\langle U\rangle$ is about 0.2, 0.4, and 0.8. Again, the C(t) curve in the semi-log scale (see the inset of figure 3(d)) becomes straighter and straighter and at $\langle E\rangle =1$, it has already become a perfect exponential function. This again confirms that the decaying behavior of the HCAF is correlated to the asymmetry degree.

Following the Green–Kubo formula, the fact that C(t) decays faster than the power law of $\sim {{t}^{-1}}$ suggests that the integral of C(t) is convergent. In addition, for $N>4\times {{10}^{3}}$, the fact that C(t)/N curves for different system sizes agree with each other as shown in figures 3(a) and (b) implies that the thermal conductivity does not depend on the system size. In figure 4, we present the heat conductivity given by the Green–Kubo formula. Note that for obtaining the heat conductivity at a given system size N by using the Green–Kubo formula (denoted by ${{\kappa}_{\text{GK}}}(N)$), a practical procedure is to truncate the time integration [6, 27] up to the time ${{\tau}_{tr}}(N)=N/\left(2{{v}_{s}}\right)$ (v_s is the sound speed of the system) [20]; i.e.

$$\begin{eqnarray}&&{{\kappa}_{\text{GK}}}(N)=\frac{1}{{{k}_{\text{B}}}{{T}^{2}}}{\int}_{0}^{{{\tau}_{tr}}(N)}\frac{C(t)}{N}\text{d}t.\end{eqnarray} \tag{ 5 }$$

Zoom In Zoom Out Reset image size

It can be seen from figure 4 that ${{\kappa}_{\text{GK}}}$ increases as N for $N<4\times {{10}^{3}}$ but gets saturated for larger N, which is in good consistence with the results of C(t)/N (see figure 3(a)).

To further verify the convergence of the heat conductivity, we compute it directly by using the nonequilibrium simulations. The temperatures of the two Langevin heat bathes applied to the two ends of L-J lattices are T₊   =  T  +  0.15 and T₋  =  T  −  0.15 with T  =  0.55. The heat conductivity measured, i.e. ${{\kappa}_{\text{neq}}}$, are also presented in figure 4 for a close comparison. Again, the heat conductivity measured in this way tends to saturate and the value it tends to is the same as that obtained by integrating the HCAF with the Green–Kubo formula.

To summarize, in the time range accessible to numerical studies, the HCAF of 1D L-J lattices has been found may decay faster than power-law manners, and the decaying behavior is correlated to the asymmetry degree of interactions. In addition, the system may follow the Fourier law at the corresponding system sizes.

The potential of the 1D FPU-α–β model is $U(x)={{x}^{2}}/2-\alpha {{x}^{3}}/3+\beta {{x}^{4}}/4,$ where the two parameters α and β determine, respectively, the asymmetric and the symmetric nonlinear term. As our aim is to investigate the effects induced by the asymmetric term, we fix $\beta =1$ throughout. Figure 5(a) shows the potential for three different values of α, and figure 5(b) shows the corresponding result of the asymmetry degree.

Zoom In Zoom Out Reset image size

From figure 5(b) it can be seen that the FPU-α–β model is quite different from the L-J model: Though in general the asymmetry degree increases as α, it generally decreases and tends to zero as the average potential energy $\langle U\rangle$ increases. In addition, for larger α the asymmetry degree may depend on $\langle U\rangle$ in a nonmonotonic manner as in the case of $\alpha =1.5$. At the first glance, these results seem to contradict to our intuition; e.g. when the temperature tends to zero, $\langle U\rangle$ decreases accordingly, and the potential energy represented by the quadratic term in the potential function would become dominant. For this reason, one may expect that the physical properties of the FPU-α–β model would converge to those of a harmonic lattice. However, we would like to point out that, for some nonlinear effects, such as thermal expansion, the cubic term in the potential function is always the leading term; hence whether a nonlinear effect is non-negligible in the low temperature limit depends on if it could stand out from the linear effect. As far as thermal expansion is concerned, in the low temperature limit the thermal expansion coefficient tends to be a temperature independent constant proportional to g/c² with g and c being, respectively, the coefficient of the cubic and the quadratic term of the Taylor expansion of the given potential function (see [24]). For the FPU-α–β model $g/{{c}^{2}}=4\alpha /3$; i.e. the larger is α, the larger is the thermal expansion coefficient, which is in good consistence with the results of the asymmetry degree given in figure 5(b). (Note that for the L-J model, the thermal expansion coefficient does not tend to zero in the low temperature limit either, and the value of g/c² for various parameter (m, n) is also in consistence with the results given in figure 1(b).)

If the decaying behavior of the HCAF is correlated to the asymmetry degree, then based on figure 5(b), we could expect that faster decay may be observed in the low temperature regime where the asymmetry degree is high enough. In the high temperature regime, the asymmetry degree tends to be low, which is consistent with the fact that as the temperature increases, the symmetric quartic term becomes overwhelmingly dominating so that the expansion quantity ${{x}_{+}}\,-\,|{{x}_{-}}|$ tends to be a constant, inducing in turn a decreasing thermal expansion rate. As a result, one may expect instead slow decay of the HCAF.

Figure 5(b) also tells that the asymmetry degree is relatively strong only at low temperatures, in particular for the potential with a smaller value of α. To reveal the effects of the asymmetric potential, we first investigate the HCAF at a very low energy density of $\langle E\rangle =0.05$. For this energy density, we have $T\approx 0.05$ and $\langle U\rangle \approx 0.025$. In figure 6(a) we show the HCAF for $\alpha =0.5$, 1, and 1.5 with N  =  8192; it can be seen that for $t<{{t}_{\text{num}}}=2\times {{10}^{4}}$, before the fluctuations (due to uncertainty of statistical average) set in, C(t) decays more rapidly—in an exponential-like way—than any power-law manners in all three cases. In addition, we have also verified that for all these cases C(t)/N does not change as N as long as N  >  10³.

Zoom In Zoom Out Reset image size

Now let us turn to the case of a high temperature: $\langle E\rangle =0.8$, which corresponds to $T\approx 0.9$ and $\langle U\rangle \approx 0.3$. Compared with the low temperature case, it can be seen from figure 5(b) that the asymmetry degree has greatly diminished. In this case, the HCAF does show a power-law tail with the decaying exponent close to  −2/3 (see figure 6(b)) as predicted by the self-consistent mode coupling theory [21, 22].

Figure 6(b) also suggests that the HCAF undergoes an initial faster decaying stage which lasts longer and longer with the increase of α. This is a signal that the asymmetric potential term even plays a role in the high temperature regime; i.e. a higher asymmetry degree can maintain a longer initial faster decaying stage. It can be expected that if the temperature is increased further, the effects induced by the asymmetric potential will become even more unnoticeable due to the rapid decrease of the asymmetric degree. This has been verified by our simulations.

Concerning the results given in figure 6(a), one may wonder if there exists a crossover from the exponential-like decay to the power-law decay even at low temperatures, but the time corresponding to the crossover point, denoted by t ^*, is so long that has gone beyond the scope of our simulations. To study this question, we computed the HCAF at several low temperatures with $\alpha =1$ and N  =  8192 (see figure 7(a)). It can be seen that such a crossover may exist and t ^* may increase very fast as the temperature decreases.

Zoom In Zoom Out Reset image size

Nevertheless, we would like to show that, even though the HCAF turns out to have a power-law tail, in practice the faster decay of the HCAF for $t<{{t}_{\text{num}}}$ as we have verified via simulations can still guarantee an effective constant heat conductivity even to a macroscopic system size [28], as long as the faster decaying stage lasts sufficiently long. To show this, we express the heat conductivity for a large system size N with the Green–Kubo formula and divide the integral of C(t) into two parts³ i.e.

$$\begin{eqnarray}&&\kappa =\frac{1}{{{k}_{\text{B}}}{{T}^{2}}}\left[{\int}_{0}^{{{t}^{\ast}}}\frac{1}{N}C(t)\text{d}t+{\int}_{{{t}^{\ast}}}^{{{\tau}_{tr}}(N)}\frac{1}{N}C(t)\text{d}t\right],\end{eqnarray} \tag{ 6 }$$

and in the time range of the second integral, i.e. $\left({{t}^{\ast}},{{\tau}_{tr}}(N)\right)$, we assume that the HCAF decays as $C(t)\sim {{t}^{-2/3}}$. Taking the case of $\alpha =1.5$ as shown in figure 6(a) as an example, where one can see that the initial rapid decay lasts at least up to ${{t}_{\text{num}}}=2.5\times {{10}^{4}}$ and $C\left({{t}_{\text{num}}}\right)/N$ drops to about 10⁻⁶. To estimate the upper bound of the second integral, we assume t ^* to be ${{t}_{\text{num}}}$. Then the first integral ${\int}_{0}^{{{t}^{\ast}}}\left(C(t)/N\right)\text{d}t\approx 4$ and the second integral ${\int}_{{{t}^{\ast}}}^{{{\tau}_{tr}}}\left(C(t)/N\right)\text{d}t\approx 2.6\times {{10}^{-3}}{{t}^{1/3}}|_{{{t}^{\ast}}}^{{{\tau}_{tr}}}$. It follows that the contribution of the power-law tail (given by the second integral) to the heat conductivity is not comparable to that of the faster decaying part until the system size reaches up to N  =  10⁹, given that [6, 20, 27] ${{\tau}_{tr}}(N)\approx N$. In other words, if we suppose that the average distance between two neighboring particles is one angstrom, then the heat conductivity would keep in effect unchanged over a wide system size range from about one micrometer to ten centimeters. So as long as the faster decay stage $\left(0,{{\tau}_{e}}\right)$ lasts long enough, as evidenced in figure 6(a), even though we assume that the HCAF has a power law tail, its influence to the heat conduction can still be safely neglected.

Above analysis raises a relevant question: If these properties remain in the FPU-β model. In figure 7(b), we show the HCAF at several temperatures for the FPU-β model with N  =  8192. Roughly, its behavior looks similar as in the FPU-α–β model in that there is a fast decay stage followed by a slow power-law tail. However, qualitative difference exists: for FPU-β model, the decaying behavior of the HCAF is of power law rather than exponential-like, and it is slower than $C(t)\sim {{t}^{-1}}$ throughout (see figure 7(b)). As a result, the thermal conductivity for FPU-β model would increase with the system size significantly.

In summary, we have numerically studied the HCAF in the 1D L-J lattices. It has been shown that the HCAF generally decays faster than power-law manners. In addition, we have observed a correlation between the decaying behavior of the HCAF and the interaction asymmetry degree: As the system parameters (m, n) change from (m, n)  =  (12, 6) to (2, 1), the asymmetry degree increases, and meanwhile the HCAF tends to decay more and more exponentially. In particular, for (m, n)  =  (2, 1), the HCAF shows clear signal of exponential decay. On the other hand, the asymmetry degree increases as the average potential energy per particle, or equivalently the temperature. Again, the HCAF tends to decay exponentially as the temperature increases.

We have also studied the decaying behavior of the HCAF in the FPU-α–β model. While the asymmetry degree increases with the asymmetry parameter α, its dependence on temperature is quite different from that of the L-J model: in the high temperature regime, the asymmetry degree decreases monotonically as the average potential energy per particle increases, in consistence with the fact that the symmetric quartic term in the potential becomes dominating. We show that in the low temperature regime the exponential-like fast decay of the HCAF can last for a significantly long time though a power-law tail may follow.

One important question that cannot be definitely answered via numerical simulations only is whether the power-law tail would show up if the system size is out of the scope accessible to simulations. Our observation is that there are not any signals of the power-law tail in the L-J model with all the system parameters investigated. However, as to the FPU-α–β model, we cannot exclude such a possibility based on our present simulation results.

Nevertheless, we have shown that for the FPU-α–β model at a low temperature, even if the HCAF has a power-law tail eventually, effective constant heat conductivity can still be expected if only the initial faster decay lasts long enough. In particular, if the HCAF keeps decaying faster over more than three orders, the contribution of the assumed power-law tail to the heat conductivity can be neglected even when the system reaches a macroscopic scale. Namely, the slow tail of the HCAF may not necessarily imply abnormal transport in practice, hence should be carefully analyzed when theoretical predictions based on the power-law tail are applied to experiments. (It is worth noting that though for the FPU-β model there is also a fast decaying period in low-temperature regime, the decaying rate is slower than  ∼1/t). Based on these analysis, we conclude that with proper system parameters, both the L-J model and the FPU-α–β model (at a sufficient low temperature) may have the normal heat conduction property given by the Fourier law. As realistic materials usually contain asymmetric interactions, we think this result may be significant to applications as well.

In a recent work, heat conduction of 1D L-J lattices was studied [13] at different particle densities. Though the authors did not show how the HCAF decays, they provided the convergent heat conductivity at low temperature regime by using the Green–Kubo formula, which implies that the HCAF decays faster than  ∼1/t. This is consistent with our earlier conclusion on this model [11]. In that work [13], the authors also conjectured that the lattices with the parabolic and/or quartic confining potential, to which the Fermi–Pasta–Ulam model belongs, should all exhibit anomalous heat transport. Our study suggests that their conjecture should be studied further because as we have shown here, the heat conduction behavior of the FPU-α–β model can also be normal. In addition, normal heat conduction has also been evidenced in other lattices with confining potential [10–12], such as the piecewise parabolic potential [12]. Therefore, more efforts are still needed to unveil the exact conditions under which the HCAF decays faster.

Finally, we emphasize that asymmetric interactions are not a sufficient condition for the faster decay of the HCAF and normal heat conduction, just as suggested by the results of the FPU-α–β model at high temperatures. Asymmetric interactions may not be a necessary condition for the faster decay of the HCAF and normal heat conduction, either. To this end, the 1D rotator lattice [29, 30] is the only known example. This model has normal heat conductivity under certain conditions, but its interaction is symmetric. However, recent studies suggest that the angle variables of this model do not constitute a conserved field [31, 32], implying that this model is in effect irrelevant with the subject we discuss here. It is worth noting that if this example is excluded, then all the 1D momentum-conserving systems where the HCAF has been reported to decay faster up to now have asymmetric interactions [10–13]. Recently, two more model systems have also been added to this category: 1D hard-point gas and 1D Toda lattice with alternative masses [33]. The faster decay of the HCAF and normal heat conduction are found in both models in certain range of the system size when the mass ratio tends to unity, i.e. when the systems approach their integrable limit. Taking into account all these studies, the faster decay of the HCAF seems to be a collective consequence of some different factors to which asymmetric interactions belong.

Very useful discussions with S Lepri, R Livi, and A Politi are gratefully acknowledged. This work is supported by the National Natural Science Foundation of China (Grants No. 11335006, No. 11275159, No. 10805036, and No. 11535011).