Wave turbulence description of interacting particles: Klein-Gordon model with a Mexican-hat potential

Wave turbulence description of interacting particles: Klein-Gordon model with a Mexican-hat potential

Basile Gallet Laboratoire SPHYNX, Service de Physique de l’État Condensé, DSM, CEA Saclay, CNRS UMR 3680, 91191 Gif-sur-Yvette, France    Sergey Nazarenko Mathematics Institute, University of Warwick, Coventry, CV4 7AL, UK    Bérengère Dubrulle Laboratoire SPHYNX, Service de Physique de l’État Condensé, DSM, CEA Saclay, CNRS UMR 3680, 91191 Gif-sur-Yvette, France
Abstract

In field theory, particles are waves or excitations that propagate on the fundamental state. In experiments or cosmological models one typically wants to compute the out-of-equilibrium evolution of a given initial distribution of such waves. Wave Turbulence deals with out-of-equilibrium ensembles of weakly nonlinear waves, and is therefore well-suited to address this problem. As an example, we consider the complex Klein-Gordon equation with a Mexican-hat potential. This simple equation displays two kinds of excitations around the fundamental state: massive particles and massless Goldstone bosons. The former are waves with a nonzero frequency for vanishing wavenumber, whereas the latter obey an acoustic dispersion relation. Using wave turbulence theory, we derive wave kinetic equations that govern the coupled evolution of the spectra of massive and massless waves. We first consider the thermodynamic solutions to these equations and study the wave condensation transition, which is the classical equivalent of Bose-Einstein condensation. We then focus on nonlocal interactions in wavenumber space: we study the decay of an ensemble massive particles into massless ones. Under rather general conditions, these massless particles accumulate at low wavenumber. We study the dynamics of waves coexisting with such a strong condensate, and we compute rigorously a nonlocal Kolmogorov-Zakharov solution, where particles are transferred non-locally to the condensate, while energy cascades towards large wave numbers through local interactions. This nonlocal cascading state constitute the intermediate asymptotics between the initial distribution of waves and the thermodynamic state reached in the long-time limit.

I Introduction

Some nonlinear partial differential equations (PDEs) are simple enough and yet they describe deep and nontrivial universal processes in a broad range of physical applications. Well-known examples of such PDEs are the Korteweg-de-Vries and the Gross-Pitaevskii (GP) equations. It would not be an exaggeration to put into such a class a universal model that allows to understand the concept of spontaneous breaking of global symmetry. Consider for instance the following complex nonlinear Klein-Gordon equation,

(1)

where . The linear part of this equation resembles the standard Klein-Gordon equation, where the square mass coefficient is replaced by . Such a negative square mass () was originally believed to lead to superluminal particles. However, it was soon realised that such a field cannot sustain a localised particle-like excitation because the base state is unstable.

Equation (1) conserves the total energy

(2)

where the last two terms correspond to the famous Mexican hat potential ; we therefore refer to equation (1) as the Klein-Gordon Mexican Hat (KGMH) model.

Equation (1) is invariant to a uniform shift in phase of , i.e., , with  const. The invariant corresponding to this symmetry is the “charge”,

(3)

where denotes the complex conjugate of . In this paper we focus mostly on the situation , being briefly described in section III.3.

In high energy physics, the KGMH model is a nonlinear -model gelman (); song () with a Mexican-hat potential. It is a simple relativistic model that describes the interaction between two fields: the first one corresponds to a massive particle, while the second one corresponds to a massless Goldstone boson. As such, the KGMH model provides a simplified description of a system of sigma mesons interacting with pions: the mass of the pions is neglected in this model, as it is much less than the one of the sigma meson.

In this paper we study the KGMH system within the framework of wave turbulence, with direct numerical simulations (DNSs) to check the predictions. Wave turbulence falkovich1992kolmogorov (); nazarenko2011wave () deals with ensembles of weakly-interacting waves. It proceeds through the derivation of a kinetic equation that governs the time evolution of the spectral density of the wave field. Among its solutions are thermal equilibrium and generalized Rayleigh-Jeans distributions. However, the main use of the kinetic equation is to describe out-of-equilibrium situations. Indeed, for scale-invariant systems, one can compute out-of-equilibrium power-law solutions to the kinetic equation: these are the celebrated Kolmogorov-Zakharov spectra. More generally, the kinetic equation governs the evolution of an ensemble of waves (or particles), starting from an arbitrary weakly nonlinear initial condition.

Such out-of-equilibrium initial conditions are encountered in particle physics and cosmology: they correspond for instance to heavy-ion collisions produced experimentally, or to the early universe. There has therefore been a significant recent interest in applying the wave turbulence approach in these domains Berges:2010ez (); Gasenzer:2013era ().

From the point of view of wave turbulence theory, the KGMH model is a new object of study, which is interesting and important for several reasons: spontaneous symmetry breaking leads to two kinds of waves, or particles, with distinct dispersion relations. This contrasts with most wave turbulence systems, where one usually deals with a single kind of waves, like e.g. in the GP equation. In addition to transfers of energy from one scale to another, the KGMH system displays “reactions”: decay of massive particles into massless ones, and fusion of massless particles into massive ones. This brings new features to the wave turbulence dynamics, with strongly nonlocal interactions between different kinds of particles, that prohibit standard Kolmogorov-Zakharov cascades.

Nevertheless, the KGMH model shares some similarities with the simpler GP equation: both have two positive invariants. Indeed, in addition to energy, the KGMH equation has an adiabatic invariant which is similar to the number of particles. Here, “adiabatic” means that this invariant is conserved exactly in the asymptotic limit of small nonlinearity only. This adiabatic invariance is sufficient to trigger wave condensation PhysRevLett.95.263901 (), a phenomenon that we study both theoretically and numerically. In the particle physics context, this corresponds to Bose-Einstein condensation of pions zimaniy (); Bunatian:1983dq (); song ().

The decay of the massive mode (-meson) into massless Goldstone bosons (pions) is the leading thread of this paper. In section II we propose a simple mechanical analogue of the KGMH model (1), which provides physical insight into the underlying dynamics. In section III we study the linear and weakly nonlinear dynamics around the minimum of the Mexican-hat potential (the uniform fundamental state). We compute the dispersion relations for the massive and massless modes, together with a reduced Hamiltonian which takes into account the dominant three-wave interaction between these modes. This dominant interaction involves two massless modes and one massive mode. We analyze it in section IV by deriving a simple set of ODEs for a single triad of interacting waves, which highlights the decay instability of a monochromatic massive wave into two massless monochromatic waves.

In section V we develop the wave turbulence description of the KGMH model. We compute kinetic equations governing the evolution of the spectra of massive and massless modes. We identify two invariants of these equations: the energy, and an adiabatic invariant somewhat similar to a number of particles. The kinetic equations admit Rayleigh-Jeans equilibrium solutions, but no local Kolmogorov-Zakharov out-of-equilibrium solutions. We focus on such thermodynamic equilibria in section VI. For low enough energy per particle, the thermodynamic state displays wave condensation, the classical counterpart to the Bose-Einstein condensation (BEC) of pions. These findings are confirmed by DNSs, which show unambiguously that BEC is observed over a very long time, even in systems where the number of particles is only an adiabatic invariant.

The out-of-equilibrium dynamics involve strongly nonlocal processes which we consider in sections VII and VIII. The first section addresses the decay instability for isotropic ensembles of massive waves: on the one hand, the kinetic equations allow to predict the distribution of the massless waves that appear during the linear stage of the instability. On the other hand, the thermodynamics of section VI predict the the reaction yield of the decay instability: starting from massive-waves (-mesons) at a single wavenumber, we predict the fraction of these waves that has decayed into massless waves (pions) in the long-time limit.

Section VIII finally addresses the behavior of ensembles of waves coexisting with a strong condensate of massless particles. We derive rigorously a reduced kinetic equation that describes a nonlocal transfer of particles to the condensate, together with a local Kolmogorov-Zakharov energy cascade towards high wave numbers. We see this cascading state as an intermediate asymptotics that describes the dynamical process of condensation: it makes the link between the initial distribution of particles and the final thermodynamic condensed state.

Ii An analogous mechanical system

Figure 1: A lattice of magnetic pendulums: small magnets placed under the vertical position of each pendulum create the Mexican-hat potential. The interaction between neighboring pendulums results from dipole-dipole interaction and torsional springs between their fixation points. In the continuous limit this system follows the KGMH equation (1).

One can gain some insight into the physics described by equation (1) by considering a system from classical mechanics that follows the same equation. An example of such system is the square lattice of pendulums sketched in figure 1. The pendulums can move in two angular directions, and . The corresponding kinetic energy of a pendulum is , with the moment of inertia of a pendulum. Each pendulum consists of a magnetized bead attached to a rigid rod, the magnetic dipole of the bead being alined with the rod. Neighboring pendulums interact through two mechanisms: first, two dipoles repulse each other. Second, we assume that the fixation points of the pendulums are connected by torsional springs. Dipole-dipole repulsion efficiently propagates longitudinal perturbations through the lattice, while torsional springs propagate transverse oscillations. The stiffness of the torsional springs is chosen to match the strength of the dipole-dipole interaction, so that the interaction energy between the pendulums reads

(4)

where and are the angular degrees of freedom of the n-th pendulum, and we consider dipole-dipole interactions between nearest neighbors only. We further assume that the relative displacements of the pendulums are small and consider the lowest order description of the dipole-dipole interaction.

Finally, to create the Mexican-hat potential we place a magnet under the rest position of each bead. Gravity drives a given pendulum towards the vertical, but the repulsion between the bead and the magnet displaces the equilibrium position away from the vertical, at a distance from the vertical axis. The combined effects of gravity and of the magnet result in a Mexican-hat-shaped potential . If the length of the rods is much longer than , we can use the approximate relation . In the neighborhood of , we approximate the potential by . The Lagrangian of the system follows from the subtraction of the Mexican-hat potential energy and neighboring interaction term to the kinetic term,

(5)

We consider the continuous limit and introduce the complex variable . After a rescaling of time, space, and this Lagrangian takes the following dimensionless form,

(6)

which corresponds to equation (1) in two spatial dimensions.

In this analogous system, the charge corresponds to the sum of the vertical angular momenta of each pendulum around their fixation points. A nonzero value of means that the initial condition (and the subsequent evolution) favors one direction of azimuthal rotation of the pendulums.

Figure 2: Two types of waves propagate in this system: branch corresponds to radial oscillations of the pendulums, whereas branch corresponds to azimuthal oscillations of the pendulums, along the minimum of the Mexican-hat potential.

Iii Weakly nonlinear dynamics around the minimum of the potential.

A key feature of the KGMH model is that it displays spontaneous symmetry-breaking: the fundamental state, defined as the state of lowest energy, is not the usual . Instead, it is a uniform state , with : is uniform, time-independent, and lies in the minimum of the Mexican-hat potential. We wish to study the dynamics of the system around this fundamental state.

Up to a phase-shift in the definition of , we can choose and consider perturbations of this state,

Substitution into equation (1) yields

(7)

We focus on the weakly-nonlinear dynamics arising from the quadratic terms of this equation, and therefore neglect the cubic nonlinearity. In variables and we have

(8)
(9)

The linear parts of these two equations provide two dispersion relations. For a plane wave proportional to , they give respectively

(10)
(11)

where . We refer to these two dispersion relations as branch and branch . They are plotted in figure 2. Branch has a nonzero frequency for vanishing wavenumber, for . In terms of particle physics, this mode has a nonzero rest energy and thus corresponds to massive particles. By contrast, the frequency of -waves vanishes for : branch corresponds to a massless Goldstone boson. This is a basic model for coexisting sigma meson (branch ) and pion (branch ) fields, in which the pion’s mass is neglected.

In the lattice of pendulums sketched in figure 1, the massive mode corresponds to oscillations of the pendulums along their radial direction. In this mode, each pendulum feels two restoring forces, one from the potential well and one from the coupling to its neighbors. By contrast, mode corresponds to oscillations of the pendulums along their azimuthal directions. The restoring force originates only from the coupling with the neighboring pendulums, and if all the pendulums move in phase along their azimuthal direction ( mode with infinite wavelength) then there is no restoring force: mode is massless.

iii.1 Hamiltonian formulation

The system (8-9) follows from the following Hamiltonian,

(12)

where and are the momenta conjugate to the variables and respectively. Hamilton equations are

(13)
(14)

We consider these equations in a -dimensional periodic cube of side and introduce Fourier series,

(15)

where the sums are over . Because these four fields are real-valued, . In Fourier space, Hamilton equations therefore become

(16)
(17)

with

(18)

where

(19)

and

(20)

where is a Kroenecker symbol, i.e., and .

iii.2 Normal variables

Let us introduce the normal variables and that diagonalize ,

(21)
(22)

We refer to and as the respective amplitudes of - and -waves (or modes, or particles).

In these variables the equations of motion are

(23)
(24)

where , with

(25)

and

(26)

where stands for “complex conjugate”, and we introduced the shorthand notations , , , and .

The terms in the Hamiltonian (26) correspond respectively to 3-wave processes of type and . In a weakly nonlinear description of the system, the only relevant nonlinearities are those for which both a wavenumber and a frequency resonance conditions can be satisfied. The other nonlinearities produce non-resonant terms only, which do not contribute to the weakly nonlinear dynamics. We show in appendix A that, out of the processes outlined above, the only one allowed by the frequency resonance condition is . We can therefore discard all the other terms in , as they do not contribute to the weakly nonlinear dynamics, and simply write

(27)

where we have switched to a shorthand notation of Kroenecker deltas, , and introduced the interaction coefficient

(28)

The equations of motion become

(29)
(30)

A similar Hamiltonian system with two types of interacting waves was studied in Zakharov1985285 (). However, in that study the dominant process was of the type , which leads to very different dynamics than the process considered here.

iii.3 Finite-charge

We have considered so far only perturbations around the minimum energy state, which is uniform lying in the minimum of the Mexican-hat potential. This corresponds to a vanishing charge invariant, . We now briefly consider nonzero values of . Writing the field in polar coordinates, , the charge reads

(31)

A nonzero charge therefore indicates a preferred direction of rotation in the complex plane. Uniform states with nonzero charge can be sought in the form , where and are real constants. Substitution into the KGMH equation (1) leads to

(32)

i.e. there is a continuous family of uniform solutions to the KGMH equation, corresponding to different values of the charge invariant .

Once again, the lattice of magnetic pendulums provides a simple interpretation for these states. From equation (31) the charge is the sum of the vertical angular momenta of the pendulums about their fixation points. The uniform solution with corresponds to all the pendulums parallel and at rest in the minimum of the Mexican-hat potential, whereas uniform solutions with correspond to all the pendulums parallel and rotating at the same angular frequency around the vertical. For the latter solutions : the field is not exactly in the minimum of the Mexican-hat potential, but instead it is “bobsleighing” on its outer edge. In the system of pendulums, this “bobsleighing” results from the centrifugal force acting on each pendulum.

To study perturbations around the uniform state with nonzero , we write , and consider infinitesimal perturbations . Substitution into the KGMH equation gives at linear order in

(33)

and after taking the real and imaginary parts, with ,

(34)
(35)

As compared to the situation , the fields and are now coupled at linear order by terms proportional to . Inserting a plane wave structure for these two fields and asking for non-trivial solutions yields the following dispersion relation

(36)

which has two branches of solutions,

(37)

In the limit , we see that branch corresponds to -modes, whereas corresponds to b-modes. For , the mode remains massive, with , and the branch remains massless, with . Note that the speed of the low- -particles is now reduced, , and the -particle’s mass is times greater than in the case. At the linear level, the structure of the problem is nevertheless quite similar to the situation.

Interesting new phenomena appear at the nonlinear level. As for the situation, one can compute normal variables and diagonalizing the quadratic part of the Hamiltonian, (resp. ) corresponding to (resp. ). These normal variables are now superpositions of the and fields, together with the corresponding momenta:

(38)
(39)
(40)
(41)

The weakly nonlinear limit is still governed by three-wave interactions, but the dispersion relations (37) now allow for new processes: and are still impossible, but, in addition to , the processes and are compatible with the dispersion relation. This possibly triggers interesting new dynamics. However, because the normal variables and have rather intricate expressions when , the computation of the kinetic equation is more involved than for . The wave turbulence analysis of the present paper therefore focuses on , and we leave the case for future work.

Iv Decay of the massive mode into two Goldstone waves

The dominant three-wave process of the weakly-nonlinear KGMH dynamics can be highlighted by studying a single triad of interacting modes. Consider the dynamical equations in the weak-nonlinearity limit, (8) and (9), and expand the fields as

(42)
(43)

where the slow time scale is . To order , equations (8) and (9) reduce to

(45)
(46)

In Appendix A, we saw that triads of interacting mode are made of two modes and one mode. We consider a single such triad in the solution to the order equations:

(47)
(48)

To obtain non-trivial dynamics at next order, we further assume that the resonance conditions and are met. To order , the equations become

(49)
(50)

The solvability condition requires the right-hand side of these equations to have no terms that are resonant with the operator on the left-hand side. We insert solutions (47) and (48) into equations (49) and (50), before demanding that the rhs of (49) has no term proportional to , and that the rhs of (50) has no terms proportional to or to . These three constraints lead to the following system of ODEs for the slow evolution of the wave amplitudes

(51)
(52)
(53)

Consider a single -wave: is a time-independent solution to the system of ODEs. We study the stability of this -wave by considering infinitesimal perturbations and . The linearized evolution equation for these perturbations is

(54)

where or . The -wave is stable, in the sense that the perturbations of and oscillate but remain small.

By contrast, if one starts off with a single -wave (, ) and perturbs it with infinitesimal and , the subsequent evolution of these perturbations is governed by

(55)

which admits some exponentially growing solutions: the -wave is unstable and decays spontaneously into two -waves. In terms of particles physics, this instability would describe the disintegration of a massive particle (e.g. the -meson) into two massless Goldstone bosons (two pions).

The growth rate of this instability is a first indication of nonlocal interactions in Fourier space: recalling the resonance conditions, is maximal when (see condition (112)), i.e., when and (or vice-versa). Nonlocal interactions occur in at least two cases: if , then . If , then either or . Such nonlocal interactions are an important feature of wave turbulence in the KGMH model.

Figure 3: Destabilization of a small-scale massive wave into massless ones. The left-hand panels are the fields, the middle panels are the fields, and the right-hand panels are the 1D spectra of (blue) and (red-dashed) as a function of . Here , and from top to bottom: , , , and . The unstable triads involve a small-scale massless mode and a large-scale one.
Figure 4: Destabilization of a large-scale massive wave into massless ones. The left-hand panels are the fields, the middle panels are the fields, and the right-hand panels are the 1D spectra of (blue) and (red-dashed) as a function of . Massless waves appear at the crossover scale . Here , and from top to bottom: , , , and .

To study the subsequent evolution of this instability, we solve equation (1) numerically. We focus on 2D numerical simulations, which corresponds to the dimension of the equivalent lattice of pendulums. We truncate the system at large wavenumbers. Such a truncation naturally arises from the discreteness of the lattice of pendulums. For more general applications of equation (1), the truncation can be seen as a crude modelization of the ultraviolet cut-off arising from Bose statistics.

Note that the domain size is an important parameter of this problem: the dispersion relation is not self-similar, with non-relativistic waves for and relativistic ones for . To accurately describe a mixture of relativistic and non-relativistic waves, we need to choose a domain size that is much larger than .

The numerical code is described in appendix D: it uses a standard pseudo-spectral method, with an exact integration of the stiff linear wave operator and Adam-Bashforth time-stepping to deal with the nonlinear terms. Dealiasing is performed using the one-half rule for this equation with cubic nonlinearity.

In this section the domain size is with a resolution of . After dealiasing we keep 128 modes in each direction: . In figure 3 we show snapshots of the and fields, together with their one-dimensional spectra, i.e., their spectra integrated over the direction of . The initial condition is a small-scale progressive -wave with , plus some weak background noise. The instability rapidly sets in and -waves appear. As can be seen in equation (55), the maximum growth rate of the instability is achieved when the product of the two wavenumbers is small. Minimizing this quantity while satisfying the resonance conditions suggest that and . Indeed, in figure 3 we observe the emergence a combination of large-scale -waves, together with -waves around the wavenumber of the initial -wave.

By contrast, in figure 4 we show similar snapshots for a simulation initiated with a large-scale -wave, with . Once again, -waves rapidly appear as a result of the decay instability. However, for such a large-scale -wave, the frequency resonance condition imposes that one of the -waves has a frequency (and therefore a wavenumber) of order unity. But then the resonance of the three wave vectors imposes that the second -wave also has a wave number of order unity and almost antiparallel to that of the first -wave. As a consequence, in figure 4 the instability sets in with -waves having wave numbers of order unity.

In the following we extend the description of these nonlocal interactions to random ensembles of waves, using the framework of wave turbulence.

V Wave turbulence

We now consider ensembles of weakly-nonlinear waves with homogeneous statistics in an infinite domain. For such systems wave turbulence can be used to derive wave kinetic equations that govern the evolution of the wave spectra.

v.1 Wave kinetic equations

Let us define the wave action spectra of the modes and :

here the angular brackets denote an ensemble average over many realizations.

In appendix B, we apply standard wave turbulence techniques to equations (29) and (30), which leads to the following kinetic equations for the evolution of the two spectra,

(56)
(57)

v.2 Conservation laws for the kinetic equations

Kinetic equations (56) and (57) conserve the wave energy,

(58)

as well as the following particle invariant,

(59)

This second invariant is nontrivial, and it can be guessed only once the dominant three-wave process is identified. Indeed, the process is similar to a chemical reaction where is a “molecule” made of two “atoms”, and would be the total number of such “atoms”: 1 per -wave, and 2 per -wave. is also the number of particles that would be in the system if the reaction were complete. For brevity, in the following we simply refer to as the number of particles.

An important comment to make here is that is only an adiabatic invariant: it is conserved by the weakly nonlinear dynamics, but there is no similar invariant for the strongly nonlinear KGMH equation. By contrast, is the quadratic part of the global invariant , the latter being conserved even in the strongly nonlinear regime.

v.3 Thermodynamic equilibrium

One can easily check that the kinetic equations admit solutions corresponding to energy equipartition in -space,

(60)

for any temperature  const 111This definition of is common in wave turbulence. In the energy equipartition state, it corresponds to an energy per oscillatory degree of freedom, where is a dimensionless Boltzmann constant, in dimensions., and solutions corresponding to an equipartition of the particle invariant,

(61)

More generally, the equilibrium Rayleigh-Jeans (or thermodynamic) solution is given by

(62)

where  const is a chemical potential.

Finally, we remark that any initial condition with a vanishing spectrum and an arbitrary spectrum is a stationary solution to the kinetic equations (56) and (57). One can also show that spectra with -waves at low wavenumber only are also steady solutions to the kinetic equations. Such distributions with only a single type of particles are described in details in sections VII and VIII.

v.4 Isotropic systems

Consider the class of isotropic spectra, and , where , and integrate out the angular variables describing the directions of and in the kinetic equations (56) and (57). Because and are independent of the directions of and , all that needs to be integrated over the angular variables is . This leads to

(63)
(64)

where

(65)

and

(66)

The integration domains and in equations (63) and (64) are determined by the triangle inequalities for and and the additional restriction of type (112):

(67)
(68)

They are sketched respectively in figures 5 and 6.

Figure 5: Integration domain for (left-hand panel) and (right-hand panel).
Figure 6: Integration domain for (left-hand panel) and (right-hand panel).

Substituting and we have

(69)
(70)

Substitution of for leads to an even simpler system,

(71)