# Fermi Pasta Ulam systems (FPU): mathematical aspects

Post-publication activity

Curator: Bob Rink

## Introduction

The Fermi Pasta Ulam (FPU) chain is a Hamiltonian model in Statistical Mechanics, that describes a nonlinear beaded string as well as a one dimensional crystal. In 1955, one of the first computer simulations in Theoretical Physics revealed its surprising recurrent behavior. This marked the starting point of a revolutionary development in Nonlinear Science.

The FPU chain was introduced by physicist and Nobel prize winner Enrico Fermi, computer expert and physicist John Pasta and mathematician Stan Ulam, in a classified scientific report with the title "Studies of nonlinear problems" (Fermi et al. 1955). Employed in the Manhattan project in Los Alamos and having access to the supercomputer MANIAC-1, these scientists decided, together with computer programmer Mary Tsingou, to simulate "a one-dimensional continuum ... with forces acting on the elements of this string." It thus seems that the model was originally meant as a discretization of a nonlinear wave equation, but nowadays crystals and DNA strands are also often modeled by the FPU chain.

The aim of the 1955 numerical experiment was to investigate the statistical properties of the chain, and in particular the question how fast a many particle system reaches thermal equilibrium, as predicted by the Boltzmann-Gibbs theory. The latter theory is strongly based on the ergodic hypothesis, which claims that many particle systems behave more or less randomly. Highly surprisingly, the simulations of Fermi, Pasta and Ulam revealed that, in the regime explored by the authors, the FPU chain does not obey the Boltzmann-Gibbs laws: it turns out that an initial long wave excitation does not properly spread to short waves within observable times. Even worse: instead of being ergodic, at low energy the FPU chain displays strong recurrent behavior, which seems to prevent it from ever reaching thermal equilibrium.

These phenomena are now called the FPU problem or FPU paradox and they continue to inspire important discoveries in Nonlinear Science. Highlights include the KAM theory of quasi-periodic motion, the discovery of integrable systems and solitons and many revolutionary developments in Chaos Theory.

## The model

The FPU chain is a mechanical model that consists of identical atoms exerting forces on their nearest neighbors, according to the second order ordinary differential equations

$$\frac{d^2 q_j}{dt^2} = F(q_{j+1}-q_j) - F(q_j-q_{j-1}) \ , \ j = \ldots, -1, 0, 1, \ldots\ .$$

The local "stress-strain relation" $$F$$ is assumed to have a Taylor expansion

$$F(x) = x + \alpha x^2 + \beta x^3 + \ldots \ .$$

In the literature one encounters finite and infinite chains and various types of boundary conditions.

## Hamiltonian systems

The differential equations for the $$q_j(t)$$ are an example of a Hamiltonian system

$$\frac{dq_j}{dt} = \frac{\partial H}{\partial p_j}\ , \ \frac{dp_j}{dt}=-\frac{\partial H}{\partial q_j} \ ,$$

where the Hamiltonian function is the sum of kinetic and potential energy$H(q, p) = \sum_j \frac{1}{2} p_j^2 + W(q_{j+1}-q_j) \ .$

Here, $$p_j = \frac{dq_j}{dt}$$ denotes the particle momentum and the function $$W$$ is an anti-derivative of $$F\ .$$

Important properties of a Hamiltonian system are that its solutions conserve the total energy $$H(q,p)$$ and that its flow is both symplectic and volume-preserving.

One should note that local potentials $$V(q_j)$$ are absent in the original FPU model. In this respect, the FPU chain differs considerably from other one-dimensional conservative interacting particle systems, such as the linear Klein-Gordon chain and its nonlinear variants, the sine-Gordon chain and the $$\phi^4$$ chain. The nonlinear Klein-Gordon chains nevertheless share much of their qualitative behavior with the FPU chain.

## Equipartition of energy

The widely accepted ergodic hypothesis says that after excitation, a many particle system visits every part of the energy surface $$H(q, p)=e$$ equally often. More precisely: for every integrable function $$f = f(q, p)$$ one has that

$$\lim_{T\to \infty} \frac{1}{T} \int_0^T \! f(q(t), p(t))\ dt = \frac{1}{|H^{-1}(e)|}\int_{H^{-1}(e)}\! f(q, p)\ d_{n-1}(q, p)\ ,$$

for almost every initial condition $$(q(0), p(0))$$ at the energy level $$e\ .$$ In the FPU experiment, this hypothesis was tested for a special set of observables, the so-called "normal mode energies". These are introduced by a Fourier transform, viewing the solution $$q(t)$$ as a sum of waves. For instance, in the case of a finite chain with fixed boundary conditions and $$n-1$$ moving particles, one writes

$$q_j = \sqrt{\frac{2}{n}}\ \sum_{k=1}^{n-1} \sin\left(j k \pi / n\right) Q_k \ .$$

In the literature, the $$Q_k$$ are called normal mode coordinates. In terms of these new variables, the Hamiltonian function is

$$H(Q, P) =\sum_{k} \frac{1}{2}(P_k^2 + \omega_k^2 Q_k^2) + \alpha H_3(Q) + \beta H_4(Q) + \ldots \ ,$$

where the frequencies $$\omega_k$$ are given by the dispersion relation

$$\omega_k =2 \sin (k\pi/2 n)\ .$$

We note that in the harmonic chain ($$F(x)=x$$) every normal mode variable $$Q_k(t)$$ simply performs a sinusoidal oscillation with frequency $$\omega_k\ .$$ Thus the solution $$q(t)$$ is simply a superposition of noninteracting waves. Moreover, each "normal mode energy"

$$E_k = \frac{1}{2}(P_k^2 + \omega_k^2 Q_k^2)$$

is a constant of motion, so that the solutions of the harmonic chain move on low-dimensional invariant tori. In particular, the harmonic chain can never be ergodic: only a nonlinear FPU chain can be.

Indeed, a nonlinearity will couple the normal modes and this may allow energy to be transferred between them. Moreover, because $$\omega_k \approx \frac{k\pi}{n}$$ is nearly linear for $$k \ll n\ ,$$ one expects that resonance causes this energy transfer to be particularly strong between acoustic (=long, that is small wave number $$k$$) waves. Hence, the energy in an acoustic excitation should quickly be equipartitioned among the all modes, causing the chain to attain what physicists refer to as a thermal equilibrium.

The numerical integrations described in the Los Alamos report were performed on FPU chains with 32 or 64 particles and with fixed boundary conditions, at moderate energy levels and with acoustic initial excitations (typically exciting the longest wavelength/lowest frequency mode only).

The astonishing outcome of the computer simulations of Fermi, Pasta and Ulam was that, within the computation time available to them, no energy equipartition was observed: energy that was initially put in one long wave, was shared by only a few normal modes. Moreover, within a rather short time nearly all the energy in the system returned to the initial wave. This recurrent behavior has later been confirmed in simulations with much larger numbers of particles and on much longer time-scales, and one is led to believe that at energies that are small enough, rather than being ergodic, the FPU chain behaves more or less quasi-periodically, as if it were an integrable Hamiltonian system.

On the other hand, as will be explained below, some very recent numerical simulations reveal that in the "thermodynamic limit" (that is for a large numbers of particles $$n$$ and total energy proportional to $$n$$), the FPU chain may actually be ergodic. To observe the convergence of time-averages though, it appears that one needs to integrate the equations of motion for unexpectedly long times. This phenomenon now goes under the name "metastability".

## Integrability

Quasi-periodic behavior is typical for "completely integrable" Hamiltonians. These are Hamiltonians $$H=H_1$$ for which there exist $$n-1$$ other Hamiltonians $$H_2, \ldots, H_n$$ such that every $$H_i$$ is a constant of motion for the Hamiltonian flow of every $$H_j \ .$$ It should be noted that such integrable Hamiltonians are extremely exceptional.

The theorem of Liouville-Arnol'd states that if a level set $$H^{-1}(h):=\{ (q,p)\ | H_1(q,p) = h_1, \ldots, H_n(q,p)=h_n \ \}$$ is compact, then each of its connected components is diffeomorphic to a torus $$\mathbb{T}^n = (\mathbb{R}/\mathbb{Z})^n \ .$$ Moreover, there exist canonical "action-angle" coordinates $$(\phi, a)\in \mathbb{T}^n\times U\ ,$$ with $$U$$ an open subset of $$\mathbb{R}^n \ ,$$ around this torus. In these action-angle coordinates, the Hamiltonian differential equations of $$H$$ take the form

$$\frac{d\phi_j}{dt} = \omega_j(a) \ , \ \frac{da_j}{dt} = 0 \ .$$

Solutions of this differential equation are simply rectilinear$\phi_j(t)=\phi_j(0)+\omega_j(a(0)) \cdot t, a_j(t)=a_j(0) \ .$ Depending on $$\omega(a) \ ,$$ they are periodic or quasi-periodic. In particular, a completely integrable Hamiltonian system is far from ergodic on its phase space $$\mathbb{T}^n\times U \ .$$

## Integrable approximations

The observation of recurrent and possibly quasi-periodic motion in the FPU experiment suggests that the Fermi Pasta Ulam system, if not integrable, can be at least well-approximated by an integrable one.

The most straightforward example of such an integrable approximation is the harmonic FPU chain, with $$F(x)=x \ .$$ The integrals of motion are the normal mode energies $$E_k\ .$$ It is known though that small perturbations of harmonic oscillators can be extremely chaotic: the harmonic approximation does not explain the FPU paradox.

However, some more interesting integrable approximations of the FPU chain were discovered already in the 1960s. We will discuss the most relevant ones below.

### The Toda lattice

The stress-strain relation $$F(x) = e^x$$ also yields an integrable Hamiltonian. The resulting chain is called the Toda lattice. In this special case the Hamiltonian equations can be rewritten as a differential equation for a matrix $$A=A(t)$$ of the form

$$\frac{dA}{dt} = AB-BA\ ,$$

for some matrix $$B=B(A,t)\ .$$ The Toda equations are said to admit a Lax-pair. Since for every $$B\ ,$$ the matrix $$AB-BA$$ is tangent to the conjugacy class of $$A\ ,$$ one concludes that all $$A(t)$$ are conjugate: there exists a family $$C=C(t)=C(t, A(0))$$ such that $$A(t) = C(t)A(0)C(t)^{-1}\ .$$ In particular, the eigenvalues of $$A$$ are constants of motion.

### Birkhoff normal forms

In certain cases, for instance for the finite FPU chain with fixed boundary conditions, an integrable approximation is found by the method of Birkhoff normal forms. This is a concept borrowed from the bifurcation theory of low-dimensional Hamiltonian systems and reminiscent of the "method of averaging". These methods are described in Rink 2006.

### The Boussinesq and Korteweg-De Vries equations

Since Fermi, Pasta and Ulam studied the slow evolution of a long wave with small amplitude, it is reasonable to make an Ansatz for a continuum approximation. Depending on the physical limit one wants to model, one may postulate for example that

$$q_j(t) = \varepsilon u(\varepsilon j, \varepsilon t)\ ,$$

where $$\varepsilon$$ is a small parameter and $$u=u(x, \tau)$$ is a smooth function. Taylor expansion of the equations for $$u$$ with respect to $$\varepsilon$$ then gives the perturbed wave equation

$$u_{\tau\tau} = u_{xx} + \varepsilon^2\left( 2 \alpha u_x u_{xx} + \frac{1}{12} u_{xxxx} \right) +\mathcal{O}(\varepsilon^4)\ .$$

Without the higher order terms, this is the well-known Boussinesq equation. The further, somewhat unmotivated, Ansatz

$$u_x(x,\tau) = U(x + \tau, \varepsilon^2 \tau) = U(X, T)$$

for a slowly modulated unidirectional pressure wave, produces a perturbed Korteweg-De Vries (KdV) equation

$$U_T = \alpha UU_X + \frac{1}{24} U_{XXX} + \mathcal{O}(\varepsilon^2)\ .$$

The KdV equation is famous for supporting a family of localized waves, the so-called solitons,

$$U(X,T)=c S(\sqrt{c}(X-cT))\ ,$$

where $$S$$ is the profile function of the soliton and the free parameter $$c>0$$ measures the wave speed, amplitude and inverse square width simultaneously. The profile $$S$$ is easily found by solving a certain conservative ODE and has the property that $$S(Z)\to 0$$ exponentially as $$Z\to \pm \infty\ .$$

In a numerical study in 1965, Zabusky and Kruskal (Zabusky et al. 1964) discovered that solitons with different speeds hardly interact when they meet: they only shift their phases . Apparently solitons can be "nonlinearly superposed". Indeed, it was later proved by Lax that the KdV equation is an infinite dimensional integrable system: it has a Lax pair and infinitely many symmetries, integrals and action-angle variables.

In fact, it was later found by Zakharov that the Boussinesq equation is integrable as well, but until now its solutions are analytically less understood. On the other hand, while the above derivation of the Boussinesq equation clearly provides an approximation to the FPU equations of motion, the KdV approximation was considered questionable among mathematicians for a long time, in spite of the apparent existence of solitary waves in the FPU chain found by Friesecke and Wattis (Friesecke et al. 1994). For more general acoustic waves, the KdV approximation was justified only recently by Bambusi and Ponno (Bambusi et al. 2006), who derived from FPU a system of two KdV equations by an infinite dimensional variant of the method of averaging. This implies that the KdV approximation of FPU is justified on intermediate time scales.

## The KAM theorem

One wonders if the existence of the above integrable approximations has implications for the dynamics of the FPU chain. A partial answer may be given by the Kolmogorov-Arnol'd-Moser (KAM) theorem, which states that under Kolmogorov's nondegeneracy condition, many of the quasi-periodic motions of an integrable Hamiltonian $$H$$ persist under small Hamiltonian perturbations. Kolmogorov's condition requires, roughly speaking, that every invariant torus $$\mathbb{T}^n\times \{a\}$$ of the integrable system has a different frequency $$\omega(a)\ .$$ So the theorem indeed does not hold for the harmonic chain for which $$\omega$$ is a constant function.

The tori that survive the perturbation are the most irrational ones. A typical KAM statement is for instance that tori for which the frequency satisfies the Diophantine inequalities

$$|\langle \omega, k \rangle | > \frac{\gamma}{(1+|k|)^{\tau}} \ \mbox{for all} \ k \in \mathbb{Z}^n\ ,$$

will survive any smooth perturbation smaller than order $$\gamma^2\ .$$ Under the further condition that $$\tau>n-1\ ,$$ the set of Diophantine $$\omega$$'s has large measure, so in nearly-integrable systems one finds many quasi-periodic motions.

The few solutions away from the KAM tori nevertheless need not stay close to an unperturbed torus, but can drift from one torus to another by a process called "Arnol'd diffusion". On the other hand, Nekhoroshev's theorem states that this diffusion is extremely slow if $$\omega$$ is a convex function of $$a\ .$$ An example of such a convex integrable system is the Toda lattice.

Using the Toda lattice or Birkhoff normal forms as integrable approximations, the existence of quasi-periodic KAM tori in the FPU chain has recently been proved, although at low energy only, see Rink 2006 and Henrici et al. 2008.

Unfortunately, KAM estimates hopelessly deteriorate as the number of particles $$n$$ grows. And in spite of the fact that KAM theorems exist for perturbations of the KdV equation, these do not apply to the FPU situation. Altogether, the physically relevant cases are still far from being understood by KAM theory.

The infinite-dimensional version of KAM-theory (KAM-tori are infinite-dimensional and stable) is found in the equation with FPU-nonlinearity by L. D. Pustyl'nikov in 1994 (See Reference below).

## Lyapunov modes

Although integrable approximations are important for our understanding of the FPU paradox, the FPU chain also exhibits aspects of nonintegrable behavior and chaos. For example, Zabusky, Sun and Peng (Zabusky et al. 2006), in a numerical investigation of the occurrence and interaction of solitons in FPU chains, have found both chaos and unstable soliton interactions in the FPU chain.

To understand the absence, emergence and presence of chaos in the FPU chain, scientists have been very interested in particular periodic solutions of the FPU equations, as well as their stability. Of particular interest are low-frequency Lyapunov modes or "q-breathers", see Flach et al. 2008. On the one hand these solutions may act as organizing centers for destabilization and chaos. On the other hand, they lie extremely close to the long wave initial conditions of the original FPU experiment, and thus considerably influence it.

## Recent developments: the metastability scenario

From a physical point of view, an important development is the discovery of the metastability phenomenon. The idea of metastability was proposed in Fucito et al. 1982 as a possible explanation for the absence of equipartition in nonlinear particle chains, but until recently the phenomenon had never been confirmed either mathematically or numerically. The metastability scenario predicts that the FPU chain will reach a thermal equilibrium after very long time scales, but before it does so, it spends a long time in a so-called metastable state in which energy is not equipartitioned. In other words: the FPU chain does thermalize, but only very slowly.

Our first theoretical understanding of metastability came with a paper of Bambusi and Ponno in 2006 (Bambusi et al. 2006) that provided the justification of the KdV approximation. This leads to a rigorous lower bound for the equipartition time-scale and explains the formation of a metastable state, i.e. a non-equilibrium state at intermediate time-scales. See also Benettin et al. 2009.

In view of a thermodynamical interpretation, the most important question is whether the metastability scenario persists in the thermodynamic limit, i.e. in the limit $$n\to\infty$$ at fixed nonzero energy per particle $$\varepsilon=E/n\ .$$ If it does, then the FPU paradox is relevant, and possibly paradoxical, for thermodynamics. Various recent numerical experiments seem to confirm that $$\varepsilon$$ is the parameter that governs the formation of the metastable state, see also Berchialla et al. 2003, Benettin et al. 2009.

Berchialla et al. (Berchialla et al. 2004) provided numerical evidence that for the quartic model ($$\alpha=0$$), the FPU chain still thermalizes in the thermodynamics limit, on time scales that far exceed the time-scale of metastability. So large equilibrium times do persist in the thermodynamic limit and follow a stretched exponential in $$1/\varepsilon\ .$$

In spite of these hopeful developments, it is until now unclear how the eventual thermalization of FPU can be explained theoretically. What exactly happens on very large time scales can be considered a completely open question.