Anomalous correlators, "ghost" waves and nonlinear standing waves in the β -FPUT system
AAnomalous correlators in nonlinear dispersive wave systems
Joseph Zaleski , Miguel Onorato and Yuri V Lvov Department of Mathematical Sciences, Rensselaer Polytechnic Institute, Troy, New York 12180, USA, Dip. di Fisica, Universit`a di Torino and INFN, Sezione di Torino, Via P. Giuria, 1 - Torino, 10125, Italy,
We show that Hamiltonian nonlinear dispersive wave systems with cubic nonlinearity and randominitial data develop, during their evolution, anomalous correlators. These are responsible for theappearance of “ghost” excitations, i.e. those characterized by negative frequencies, in addition tothe positive ones predicted by the linear dispersion relation. We use generalization of the Wick’sdecomposition and the wave turbulence theory to explain theoretically the existence of anomalouscorrelators. We test our theory on the celebrated β -Fermi-Pasta-Ulam-Tsingou chain and show thatnumerically measured values of the anomalous correlators agree, in the weakly nonlinear regime,with our analytical predictions. We also predict that similar phenomena will occur in other nonlinearsystems dominated by nonlinear interactions, including surface gravity waves. Our results pave theroad to study phase correlations in the Fourier space for weakly nonlinear dispersive wave systems. Keywords:
Nonlinear waves , Fermi-Pasta-Ulam-Tsingou chain, Anomalous correlators
I. INTRODUCTION
Wave Turbulence theory has led to successful predic-tions on the wave spectrum in many fields of physics[1, 2]. In this framework the system is represented asa superposition of a large number of weakly interactingwaves with the complex normal variables a k = a ( k, t ).In its essence, the classical Wave Turbulence theory is aperturbation expansion in the amplitude a k of the nonlin-earity, yielding, at the leading order, to a system of quasi-linear waves whose amplitudes are slowly modulated byresonant nonlinear interactions [1–6]. This modulationleads to a redistribution of the spectral energy densityamong length-scales, and is described by a wave kineticequation. One way to derive the wave kinetic equationis to use the random phase and amplitude approach de-veloped in [2, 7, 8]. The initial state of the system canbe always prepared so that the assumption of randomphases and amplitudes is true. Whether the phases re-main random in the evolution of the system has been anissue of intense discussions. In Wave Turbulence theory,the standard object to look at is the second-order corre-lator, (cid:104) a k ( t ) a ∗ l ( t ) (cid:105) , where (cid:104) . . . (cid:105) is an average over an en-semble of initial conditions with different random phasesand amplitudes. As will be clear later on, under thehomogeneity assumption, the second order correlator isrelated to the wave action spectral density function, i.e.the wave spectrum , n k = n ( k, t ). However, one shouldnote that the complex normal variable, as defined in thethe Wave Turbulence theory, is a complex function alsoin physical space. Therefore, the second-order statisticsare not fully determined by the above correlator. Theso called “anomalous correlator”, (cid:104) a k ( t ) a l ( t ) (cid:105) , see [9, 10],needs also to be computed. Under the hypothesis of ho-mogeneity, will be related to the anomalous spectrum , m k = m ( k, t ), to be defined in the next Section. Indeed,if phases are totally random, this quantity would be zero. We show that, in the nonlinear evolution of the system,this is not the case. Far from it, this quantity is stronglynonzero and, in the limit of weak nonlinearity, we predictanalytically and verify numerically its value.Our ideas are based on the extension of the Wave Tur-bulence Theory to include these anomalous correlators.Notably, conventional Wave Turbulence Theory has beensuccessful in the understanding of the spectral energytransfer in complex wave systems such as the ocean [11],optics [12] and Bose-Einstein condensates [13], one di-mensional chains [14], and magnets [15]. Analogously,anomalous correlators first appeared in the well knownBardeen, Cooper, and Schriffer (BCS) theory of super-conductivity [26]. Subsequently, anomalous correlatorshave been studied in S–theory [9, 10].Recently, anomalous correlations were shown to playan important role in explaining numerical observationsof nondecaying oscillations around a steady state in aturbulence–condensate system modeled by the Nonlin-ear Schr¨odinger equation [27–29]. Such oscillations, cor-responding to a fraction of the wave action being periodi-cally converted from the condensate to the turbulent partof the spectrum, were shown to be directly due to phasecoherence [27]. In [23] a system of Coupled NonlinearSchr¨odinger equations has been considered and specificattention was focussed on the phenomena of recurrenceof incoherent waves observed in the early stages of thedynamics. The authors derived a variant of the kineticequation which includes anomalous correlators; the pe-culiarity of such an equation is that it is capable of de-scribing properly the recurrence phenomena observed inthe simulations.One of the main tools used to derive the theory is theWick’s contraction rule that allows one to split higherorder correlators as a sum of products of second ordercorrelators, plus cumulants. To explain analytically theexistence of the anomalous correlators, it is necessaryto use the more general form of the Wick’s decompo-sition, namely the form that allows anomalous correla-tors. We then demonstrate that the anomalous correla-tors are responsible for creating the “ghost waves”, i.e. a r X i v : . [ n li n . C D ] M a r the waves with the frequency equal to the negative ofthe frequency predicted by the linear dispersion relation-ship. These ideas are tested on a simple, but non trivial,system, i.e. the β -Fermi-Pasta-Ulam-Tsingou (FPUT)chain. The chain model was introduced in the fiftiesto study the thermal equipartition in crystals [16]; itconsists of N identical masses, each one connected bya nonlinear string; the elastic force can be expressed asa power series in the displacement from equilibrium.Fermi, Pasta, Ulam and Tsingou integrated numericallythe equations of motion and conjectured that, after manyiterations, the system would exhibit a thermalization, i.e.a state in which the influence of the initial modes disap-pears and the system becomes random, with all modesexcited equally (equipartition of energy) on average. Suc-cessful predictions on the time scale of equipartition havebeen recently obtained in [14, 17, 18] using the WaveTurbulence approach. In this paper we perform exten-sive numerical simulations with initial random data andlook all at the possible excitations, once a thermalizedstate has been reached. This is all done by analyzingthe spatial -temporal ( k − Ω) spectrum, i.e. the square ofthe space-time Fourier transform of the wave amplitudes.Analyses of the effective dispersion relation in the non-linear system is a well known and widely used theoreticaland numerical tool, see for example [19].We give numerical evidence that in addition to the“normal” waves with frequency ω predicted by the lin-ear dispersion relation for wave number k , there arethe “ghost” excitations with the negative frequencies.Our theoretical analysis reveals that the origin of those“ghost” excitations resides on the nonzero values of thesecond-order anomalous correlator. II. THE MODEL
The theory that we develop hereafter applies to anysystem with cubic nonlinearity. Examples of such sys-tems among others, include deep water surface gravitywaves [20], Nonlinear Klein Gordon [17], β -Fermi-Pasta-Ulam-Tsingou chain. In normal variables a k the Hamil-tonian of these systems assumes the canonical form: H = (cid:88) k ω k | a k | + (cid:88) k ,k ,k ,k (cid:2) ( T (1)1234 a ∗ a a a + c.c ) δ + 12 T (2)1234 a ∗ a ∗ a a δ + 14 T (4)1234 ( a ∗ a ∗ a ∗ a ∗ + c.c ) δ (cid:3) , (1)where ω k = ω ( k ) are the positive frequencies associatedto the wave numbers via the dispersion relation, T ( i )1234 arecoefficients that depend on the problem considered andsatisfy specific symmetries for the system to be Hamilto-nian, c.c. implies complex conjugation, a j = a ( k j , t ) arethe complex normal variables, δ lmij = δ ( k i + k j − k l − k m ) isthe Kronecker Delta. We assume that the only resonantinteractions possible are the ones for which the following two relations are satisfied for a set of wave numbers k + k = k + k , ω ( k ) + ω ( k ) = ω ( k ) + ω ( k ) . (2)With the objective of presenting some comparison withnumerical simulations, out of many physical systemsdescribed by the above Hamiltonian, we select a sim-ple one dimensional system, the β -Fermi-Pasta-Ulam-Tsingou chain. Modeling a vibrating string, this prob-lem consists of a system of N identical particles con-nected locally to each other by a nonlinear oscillator. Inthe physical space the displacements with respect to theequilibrium position q j ( t ) and their momenta p j ( t ), theHamiltonian takes the following form: H = H + H (3)with H = N (cid:88) j =1 (cid:18) p j + 12 ( q j − q j +1 ) (cid:19) ,H = β N (cid:88) j =1 ( q j − q j +1 ) . (4) β is the nonlinear spring coefficient (without loss of gen-erality, we have set the masses and the linear spring con-stant equal to 1). The Newton’s law in physical space isgiven by:¨ q j = ( q j +1 + q j − − q j ) + β (cid:2) ( q j +1 − q j ) − ( q j − q j − ) (cid:3) . (5)We assume periodic boundary condition; our approachis developed in Fourier space and the following definitionsof the direct and inverse Discrete Fourier Transforms areadopted: Q k = 1 N N − (cid:88) j =0 q j e − i πkj/N , q j = N/ (cid:88) k = − N/ Q k e i πjk/N , (6)where k are discrete wave numbers and Q k are the Fourieramplitudes. The displacement q j and momentum p j ofthe j particle are linked by canonically conjugated Hamil-ton equations ˙ p j = − ∂H∂q j , ˙ q j = ∂H∂p j . We then perform the Fourier transformation to Fourierimages of position and momenta, and then additionalcanonical transformation to complex amplitude a k givenby a k = 1 √ ω k ( ω k Q k + iP k ) , (7)where ω k = 2 | sin( πk/N ) | > Q k and P k are theFourier amplitudes of q j and p j , respectively. In termsof a k the equation of motion reads, see [21]: i da dt = ω k a + (cid:88) k ,k ,k (cid:0) T (1)1234 a a a δ + T (2)1234 a ∗ a a δ ++ T (3)1234 a ∗ a ∗ a δ + T (4)1234 a ∗ a ∗ a ∗ δ (cid:1) , (8)where all wave numbers k , k and k are summed from0 to N − δ cd..ab.. = δ ( k a + k b + ... − k c − k d − ... )is the generalized Kronecker Delta that accounts for aperiodic Fourier space, i.e. its value is one when theargument is equal to 0 (mod N ). The matrix elements T (1)1234 , T (2)1234 , T (3)1234 , prescribe the strength of interactionsof wave numbers k , k , k and k . Their values (30) aregiven in Appendix A. A. The ( k − Ω) spectrum The main statistical object discussed in this paper isthe wave number-frequency ( k − Ω) spectrum. Startingfrom the complex amplitude a ( k, t ) we take the Fouriertransform in time so that we get a ( k, Ω); Under the hy-pothesis of homogeneous and stationary conditions, thesecond-order ( k − Ω) correlator takes the following form (cid:104) a ( k i , Ω p ) a ( k j , Ω q ) ∗ (cid:105) = N ( k i , Ω p ) δ ( k i − k j ) δ (Ω p − Ω q ) , (9)where (cid:104) .... (cid:105) implies averages over initial conditions withdifferent random phases. N ( k, Ω) is the ( k − Ω) spectrumdefined as follows: N ( a ) ( k, Ω) = 12 π N (cid:90) + ∞−∞ N (cid:88) l =1 R ( l, τ )e − i πkl/N e − i Ω τ dτ, (10)with R ( l, τ ) = (cid:104) a j ( t ) ∗ a j + l ( t + τ ) (cid:105) is the space-time auto-correlation function. The linear ( k − Ω) spectrum - Before diving into the non-linear dynamics, we discuss the predictions in the linearregime. Therefore, we start by neglecting the nonlinear-ity in equation (8) and find the solution in the form a k ( t ) = a k ( t ) e − iω k t . (11)where t is a time at which the solution is known or aninitial condition. We then take Fourier transform in time a ( k, Ω) = a ( k, t ) δ (Ω − ω k ) (12)After multiplication by its complex conjugate and takingaverages over different realizations with the same statis-tics, we get: N ( a ) ( k, Ω) = n ( a ) ( k, t ) δ (Ω − ω k ) , (13)where n ( a ) ( k, t ) is the standard wave spectrum at time t related to the second-order correlator as (cid:104) a ( k i , t ) a ( k i , t ) ∗ (cid:105) = n ( a ) ( k i , t ) δ ( k i − k j ) . (14) and defined via the autocorrelation function as n ( a ) ( k i , t ) = 1 N (cid:88) l (cid:104) a j ( t ) a j + l ( t ) ∗ (cid:105) e − i πkl/N . (15)In the linear regime n ( a ) ( k i , t ) does not evolve in time.Equation (13) implies that in the linear case the ( k − Ω)spectrum is different from zero only for those values of Ωand k for which the dispersion relation is satisfied. Notethat in this formulation ω k is defined as a positive quan-tity; therefore, only the positive branch of the dispersionrelation curve appears in the linear regime. B. Numerical results for the ( k − Ω) spectrum We now test the predictions from equation (13) bothin the linear regime and observe what happens to it inthe nonlinear regime. We perform numerical simulationsof the equations (5) using a symplectic algorithm, see[22]. We use 32 particles in the simulations; such choiceis completely uninfluential for the results presented be-low. In the linear regime, we just prescribe a thermalizedspectrum with some initial random phases of the waveamplitudes a k and evolve the system in time up to a de-sired final time; a Fourier Transform in time is then takento build the ( k − Ω) spectrum. In the nonlinear regime weperform long simulations up to a thermalized spectrum.For a given nonlinearity, 1000 realizations characterizedby different random phases are made and ensemble aver-ages are considered to compute the ( k − Ω) spectrum. Allsimulations have the same initial linear energy and, froman operative point of view, the only difference betweenthem is the value of β . To characterize the strength ofthe nonlinearity, we use the following ratio between non-linear and linear Hamiltonians at the beginning of eachsimulation: (cid:15) = H H ∝ β (16)Results are shown in Figure 1, where, for different valuesof the nonlinear parameter (cid:15) , the spectrum N ( a ) ( k, Ω) isplotted using a colored logarithmic scale. We first focusour attention on the linear regime, (cid:15) = 0: results areshown in Figure 1(a). As well predicted by the the-ory, the plot shows dots in the positive frequency plane,where the frequencies Ω and wave numbers k satisfy thelinear dispersion curve ω k . Increasing the nonlinearity,Figures 1 -(b,c,d), two well known effects appears: thefirst one is a shift of the frequencies, due to nonlinear-ity (this is more evident in Figures 1 -(c,d) where thefrequency scale in the vertical axes has been changed).The second one is the broadening of the frequencies; thisis related to the fact that the amplitude for each wavenumber is not constant in time; therefore, the amplitude-dependent frequencies are not constant in time and theyoscillate around a mean value with some fluctuations.Those results are well understood, at least in the weakly FIG. 1. ( k − Ω) spectrum, N ( a ) ( k, Ω), for different values of (cid:15) : ( a ) (cid:15) = 0, ( b ) (cid:15) = 0 . c ) (cid:15) = 0 .
089 ( d ) (cid:15) = 1 .
12. In thelinear case, ( a ), the N ( a ) ( k, Ω) is different from 0 only when the frequency Ω matches the linear dispersion relation. As thenonlinearity is increased, ( b − d ), a frequency shift, a broadening of the frequencies and a lower branch less intense than theupper one are visible. Waves with negative frequencies are named “ghost” excitations.FIG. 2. Ratio between the number of “ghost” excitations,N ghost , over the total number of waves, N tot ,as a function ofthe nonlinearity. nonlinear regime, and can be predicted using Wave Tur-bulence tools, see [2, 18]. Besides these two effects, start-ing from Figure 1-(b), the presence of a lower branch,whose intensity is much less than the upper one, startsto be visible. The lower curve becomes more importantand, when the nonlinearity is of order one, is of the sameorder of magnitude of the upper one. The total numberof waves in the simulation, N tot , is given by the integral over Ω and the sum over all k of the function N ( a ) ( k, Ω).In the weakly nonlinear regime, N tot is an adiabatic in-variant of the equation of motion (5); the plot highlightsthe existence of waves with negative frequencies, whichwill be named “ghost” excitations. One of the scopesof the present paper is the understanding of the originof such waves. Before entering into the discussion, weshow in Figure 2 the ratio of “ghost” excitations, N ghost ,i.e. N ( a ) ( k, Ω) integrated over negative frequencies andsummed over all wave numbers, divided by the total num-ber of waves, N tot . As can be seen from the plot, there isa monotonic growth of the “ghost” waves that, for verylarge nonlinearity, can reach values up to 25% of the totalnumber.
III. ANOMALOUS CORRELATORS
To explain the presence of “ghost” excitations, we in-troduce the so called second-order anomalous correlator[9, 10, 23]: (cid:104) a k ( t ) a j ( t ) (cid:105) = m k ( t ) δ ( k + j ) , (17)with the anomalous spectrum defined as m ( a ) k ( t ) = 1 N (cid:88) l (cid:104) a j a j + l (cid:105) e − i πkl/N . (18)Similarly, we also introduce the second-order ( k − Ω)anomalous correlator: (cid:104) a ( k i , Ω l ) a ( k j , Ω m ) (cid:105) = M ( a ) ( k i , Ω l ) δ ( k i + k j ) δ (Ω l + Ω m ) , (19)where M ( a ) ( k, Ω) = 12 π N (cid:90) + ∞−∞ N (cid:88) l =1 S ( l, τ )e − i πkl/N e − i Ω τ dτ (20)and S ( l, τ ) = (cid:104) a j ( t ) a j + l ( t + τ ) (cid:105) . The presence in equa-tions (17) and (19) of the Kronecker δ over wave numbersand the Dirac δ over frequency, are related to the hypoth-esis of statistical homogeneity and stationarity, respec-tively. Note that M ( a ) ( k, Ω) is not the Fourier transformin time of m ( a ) k ( t ) and in general both can be complexfunctions. To verify numerically that the anomalous cor-relator is indeed nonzero, we measure numerically thereal part of the second-order correlator (cid:104) a k i ( t ) a k j ( t ) (cid:105) asa function of k and k . Results are plotted in Figure 3where we show the results of two numerical simulationscharacterized by two different values of the nonlinear pa-rameter, (a) (cid:15) = 0 . (cid:15) = 1 .
12. In both cases,a diagonal contribution is visible, pointing out the exis-tence of anomalous correlators in the β -FPUT model. Generalization of the Wick’s decomposition -
Using(17), it is straightforward to extend the Wick’s decompo-sition by taking into account the anomalous correlators,as done in [15]: (cid:104) a ∗ k a ∗ l a p a n (cid:105) = n k n l ( δ kp δ ln + δ kn δ lp ) + m ∗ k m p δ kl δ pn , (cid:104) a ∗ k a l a p a n (cid:105) = n k m p ( δ lk δ pn + δ nk δ lp ) + n k δ pk m l δ nl , (cid:104) a k a l a p a n (cid:105) = m k m l ( δ kp δ ln + δ kl δ pn + δ kl δ pn ) . (21)The above relations will be fundamental for making anatural closure of the moments when calculating analyt-ically the ( k − Ω) spectrum.In Figure 4, we find further evidence justifying this de-composition by plotting the real part of the fourth-ordercorrelator (cid:104) a k a k a ∗ k a ∗ k + k − k (cid:105) with k = 20, computedfrom numerical simulations for (a) (cid:15) = 0 . (cid:15) = 1 .
12. The diagonal lines in both figures, highlightingthe contribution from the second-order anomalous corre-lator, are noticeable. The vertical and horizontal linescorrespond to the trivial resonances in which two wavenumbers are equal (mod N ). A. Theoretical prediction for the anomalouscorrelator in the weakly nonlinear regime
A key step for the development of a theory for theanomalous correlator is the change of variable (near iden-tity transformation) which allows one to remove bound modes, i.e. those modes that are phase locked to free modes and do not obey the linear dispersion relation.The procedure is well known in Hamiltonian mechanicsand well documented for example in [1]. We accomplishthis via the following canonical transformation from vari-able a k ( t ) to b k ( t ) a = b + (cid:88) k ,k ,k (cid:2) B (1)1234 b b b δ + B (3)1234 b ∗ b ∗ b δ ++ B (4)1234 b ∗ b ∗ b ∗ δ (cid:3) , (22)with the coefficients B ( i )1234 selected in such a way to re-move non resonant terms in the original Hamiltonian [24].Their values are given in Appendix A.The transformation is asymptotic in the sense that thesmall amplitude approximation is made and the terms inthe sums on the right hand side are much smaller thanthe leading order term b . The evolution equation forvariable b k ( t ) contains resonant interactions and take thefollowing standard form: i db ∂t = ω b + (cid:88) k ,k ,k T (2)1234 b ∗ b b δ + h.o.t.where higher order terms arising from the transformationhave been neglected.Using the transformation (22) and the generalizedWick’s decomposition (21), we can now build the timeaveraged anomalous spectrum (for details, see AppendixB): (cid:104) m ( a ) k ( t ) (cid:105) t = 2 (cid:16) n ( a ) k + n ( a ) − k (cid:17) (cid:88) j B (3) k, − k,j,j n ( a ) j , (23)where (cid:104) ... (cid:105) t implies averaging over time. For the β –FPUTsystem in thermal equilibrium, where n ( a ) k = T /ω k with T constant, (23) reduces to ω k |(cid:104) m ( a ) k ( t ) (cid:105) t | = 3 N T β . (24)In Figure 5 we compare this prediction for ω k |(cid:104) m ( a ) k ( t ) (cid:105) t | in thermal equilibrium to the values givenby numerical simulations for varying values of nonlin-earity: the results are in good agreement in the weaklynonlinear regime, (cid:15) < . m k ( t ); the subsequent time av-eraging window used was 10 with a sample spacing of∆ t = 0 .
1. For larger nonlinearity, is expected that higherorder terms play a role in the evolution of the anomalouscorrelator.In Figure 6 we show the time evolution of the firstfive modes of |(cid:104) ω k m ( a ) k ( t a ) (cid:105) t a 12. A diagonalcontribution corresponding to k = − k is evident in both figures. As the nonlinearity is increased, the contribution becomeslarger.FIG. 4. Fourth-order correlator | R e [ (cid:104) a k a k a ∗ k a ∗ k + k − k (cid:105) ] | with k =20. (a) (cid:15) = 0 . (cid:15) = 1 . 12. Different horizontal,vertical and diagonal lines are visible. Horizontal and vertical lines corresponds to trivial resonances: k = k , vertical line; k = k , horizontal line; k = − k + N , diagonal line. The latter line corresponds to the presence of an anomalous second-ordercorrelator. The intensity of the lines is larger for larger nonlinearity. is indeed initially zero due to the randomness of phases.Here we use a larger value of nonlinearity (cid:15) = 10 toshow the development of the anomalous correlator in ashorter time window. The amplitudes were initializedso that | a k ( t = 0) | = (cid:113) N − k for k = ± , ± , ± 3, withhigher modes zero. The phases were initially normallydistributed. We observe that the anomalous correlatorgrows with time, reaching a peak in modes 1, 2, 3, be-fore it eventually saturates between all modes equally,with ω k |(cid:104) m ( a ) k ( t a ) (cid:105) t We have now developed all the tools for predicting an-alytically the ( k − Ω) spectrum as defined in the equation(9). Taking the Fourier Transform in time of the canon-ical transformation (see appendix C), using the gener-alized Wick’s decomposition and the hypothesis of sta-tistical stationarity and homogeneity, we get at leadingorder: N ( a ) ( k, Ω) = n ( b ) ( k, t ) δ (Ω − ω k )++ F ( k )Re[ m ( b ) ( k, t )] δ (Ω + ω k ) (25) FIG. 5. Comparison of anomalous spectrum ω k |(cid:104) m ( a ) k ( t ) (cid:105) t | as observed in numerical simulations (dots) with theo-retical predictions given by (24) (dashed lines) for (cid:15) =0 . , . , . , . . FIG. 6. Time evolution of the first five modes of the averagedquantity |(cid:104) ω k m ( a ) k ( t a ) (cid:105) t a 03; for larger non-linearity the theoretical prediction is considerably larger. V. NONLINEAR STANDING WAVES The development of a regime characterized by ananomalous spectrum corresponds to a tendency for thesystem to develop standing waves in the original displace-ment variable q j ( t ). Indeed, the existence of an anoma-lous spectrum implies a correlation between positive andnegative wave numbers. The connection between theanomalous correlator and standing waves can be seen onthe following illustrative example. Consider the restric-tive ensemble of realizations of the linear system whereamplitudes and phases are initiated in Fourier space witha correlation between wave numbers k = 1 and k = − a k ( t = 0) = A e − iφ , if k = 1 A e iφ , if k = − , otherwise , (28)with the random phase φ . In terms of the displacementvariables this would correspond to the system being ini-tially at rest and displaced from equilibrium as a singlewave q j ( t = 0) = 2 A (cid:114) ω cos (cid:16) πjN − φ (cid:17) . Since the system is assumed to be linear, the time evo-lution of complex amplitudes a and a − will be givenby a ± ( t ) (cid:39) A e − i ( ω t ± φ i ) . Averaging over random phase φ , the anomalous corre-lator becomes m = (cid:104) a a − (cid:105) = A e − iω t , analogous to the oscillating term of equation (36) for theanomalous correlator in the nonlinear case with ampli-tudes and phases being initially completely random. Interms of displacement, such initial conditions give q j ( t ) = 2 A (cid:114) ω cos (cid:16) πjN − φ (cid:17) cos( ω t ) , which corresponds to the standing wave pattern. Thus,we see that the phase and amplitude correlations whichresult in a nonzero anomalous correlator are directlylinked to the formation of standing waves in this par-ticular example.This consideration can be generalized for the case ofweakly nonlinear systems and more general initial con-ditions. Indeed, for weakly nonlinear systems the ampli-tudes | a | and | a − | will be changing slowly over manyoscillations, thus maintaining strongly nonzero anoma-lous correlator and standing waves.In Figure 8-(a) we numerically solve the equations ofmotion with initial conditions given by (28). Here weplot a colormap of the displacement q j ( t ) for all massesas a function of time as the system reaches the timescalerequired for statistical thermal equilibrium. The nonlin-earity parameter (cid:15) = 4 . 74, in the regime of strong non-linearity and outside the regime of validity of our theory.Nevertheless, we initially consider this example to displayhow the system behaves when the phase correlations de-velop rapidly. The existence of several regions of standingwave behavior are clearly visible as darker regions in theimage, as the inset Figure 8-(b) shows.It is important to emphasize that Figure 8 showsa single realization of the system, while correlators m k ( t ) , n k ( t ) describe statistical ensemble–averaged quan-tities. Thus the existence of the standing wave patternsis not in violation of the presumed assumption of spatialhomogeneity.Below we give numerical evidence that such coherentstructures can also be observed for smaller values of non-linearity that are within the regime of validity of ourtheory.In Figure 9-(a) we plot the displacement as a functionof time for the system with (cid:15) = 0 . 02, a value of nonlin-earity well within the regime of agreement of our theoryas shown in Figures 7, 5. Here, we prescribe initial condi-tions so that the total energy is initially in the first wavenumber, i.e. a k ( t = 0) = 0 for all k (cid:54) = 1, and plot a singlerealization. This corresponds to a pure traveling wavesolution in the linear system; indeed as seen in Figure 9-(a), the system is initially a traveling wave, representedby series of slanted parallel lines in the colormap of q j ( t ).Conversely, in Figure 9-(b) we show that by the time thesystem has reached the timescale required for statisti-cal thermal equilibrium, a prominent standing wave hasdeveloped, due to the phase correlations between posi-tive and negative lowest wave numbers. Notably, phasecorrelations are not restricted to only the lowest wavenumbers. To emphasize this, we consider the following FIG. 8. Color map of the displacement q j ( t ) for the systemwith (cid:15) = 4 . 74 initialized with particles at rest with initialpositions as a single sine wave. A nonlinear standing wavepattern is visible. spatial frequency filter applied to the displacement˜ q j ( t ) = N/ (cid:88) k = − N/ H k Q k e i πjk/N , (29)where H k = (cid:40) , if k = 5 , , otherwise , is selected to only show thewaves with frequencies corresponding to k = 5 , q j ( t ) in Figures 9-(c) and its inset 9-(d). Here we clearly still observe thesestanding waves in the selected unfiltered wave numbers,meaning that the coherent structures are not limited tothe lowest wave number. Our choice of displaying wavenumbers 5 , (cid:15) = 0 . 54 just outsidethe range of applicability of our theory. We plot theinitial time evolution of the displacements in Figure 10-(a), the time evolution of the displacements in thermalequilibrium in Figure 10-(b), and the displacements af-ter applying a spatial frequency filter to emphasize wavenumber k = 4 in Figure 10-(c). We note the arbitraryfluctuations between the coherent standing waves and be-tween traveling waves in Figure 10-(b,c). VI. CONCLUSION In this paper we have given the numerical evidence thatanomalous correlators develop spontaneously in a classi-cal system. From a theoretical point of view it is possibleto develop a theory for weakly nonlinear dispersive wavesthat accounts for presence of such anomalous correlator.The framework in which the theory has been developedis the Wave Turbulence one. In such theory one usuallyis interested in the second order correlator (cid:104) a k i a ∗ k j (cid:105) whichis strictly related to the wave action spectrum. However,what is clear from numerical simulations of the β -FPUTsystem is that also the correlator (cid:104) a k i a k j (cid:105) can assume FIG. 9. Color map of the displacement q j ( t ) for (cid:15) = 0 . 02: (a) initial traveling wave, 0 < t < < t < + 1000 (c) ˜ q j ( t ), displacement after removing wave numbers k = 1 ... , ...N (d) closerlook at the boxed region in (c). values that are different from zero. This finding has con-sequences on the standard Wave Turbulence theory thatis based on the Wick’s selection rule, i.e. the splitting ofhigher order correlators as a sum of products of secondorder correlators. Following [1, 15], we have generalizedthe Wick’s rule by including the anomalous correlators.We note that we differ from the case described in theS–theory [9, 10] in that there the existence of anomalouscorrelators was connected with coherent pumping in thesystem, with the anomalous correlator being a measureof partial coherence for exiting waves. In our observa-tions and predictions, waves with random initial condi-tions form phase correlations with each other, resulting inan anomalous correlator which is initially zero but thensaturates to a nonzero value as it evolves with time.One of the most striking manifestation of those corre-lators is the appearance of “ghost excitations”, i.e. thosecharacterized by a negative frequencies. A formula forthe content of energy of such excitations as a function ofthe wave spectrum is obtained and compared favorably,in the weakly nonlinear, regime with numerical simula-tions. Moreover, we have shown that the spontaneousemergence of the anomalous correlator is strongly con-nected with the formation of nonlinear standing waves;indeed, the presence of those waves implies a strong cor-relation between the phases of positive and negative wavenumbers.Our approach paves a new road to investigate disper-sive nonlinear systems by taking into account not onlyamplitudes of the waves, as in traditional wave turbu-lence, but also the phases of the waves. We conjecture that the anomalous correlators play an important role inthe theory of extreme events, such as rogue waves, whichform via a mechanism related to phase locking betweendifferent wave numbers [25]. Phase locking also leads tothe existence of solitons in nonlinear media.As was discussed in the introduction, anomalous phasecorrelations have been observed to play a role in caus-ing shifts of wave action from turbulence and conden-sate in the Nonlinear Schrodinger equation [27]. Ourapproach of extending Wave Turbulence Theory to in-clude the anomalous correlator could be generalized toaddress the role these correlations play in the statisti-cal properties of the Nonlinear Schrodinger equation andother integrable systems. On a similar note, recurrencesin an NLS-like model were shown to be directly relatedto the formation of anomalous phase correlations [23];further investigating FPUT recurrences potential ties tothe anomalous correlator is a subject of current work.Finally, we emphasize that the Hamiltonian we con-sidered is of the same family as the one for surface grav-ity waves (after removing by a canonical transformationnonresonant three wave interactions). We predict thatalso the anomalous correlators will play an importantrole in the understanding of statistical properties of oceanwaves. Acknowledgments The authors are grateful to Dr.B. Giulinico for discussions. We are grateful to anony-mous referees, who’s insightful suggestions improved themanuscript considerably. M. O. has been funded by Pro-getto di Ricerca d’Ateneo CSTO160004, by the “Depart-ments of Excellence 2018/2022” Grant awarded by the0 FIG. 10. Color map of the displacement q j ( t ) for (cid:15) = 0 . 54: (a) initial traveling wave, 0 < t < < t < + 1000 (c) ˜ q j ( t ), displacement after removing wave numbers k = 1 ... , ...N (d) closerlook at the boxed region in (c). Italian Ministry of Education, University and Research(MIUR) (L.232/2016) and by Simons Foundation, WaveTurbulence. JZ and YL acknowledge support from NSFOCE grant 1635866. YL acknowledges support fromONR grant N00014-17-1-2852. Appendix A Matrix element in eq. (8) - The matrix elements gov-erning four wave interactions for the variable a k ( t ) are: T (1)1234 = − β e iπ ( − k + k + k + k ) /N (cid:89) i =1 πk i /N ) √ ω i ,T (2)1234 = − T (1)1 − , T (3)1234 = 3 T (1)4231 T (4)1234 = − T (1) − . (30) Matrix elements in the canonical transformation, eq.(22) - The coefficients in eq. (8) suitable for removingnonresonant terms are given by: B (1)1234 = T (1)1234 ω − ω − ω − ω ,B (3)1234 = T (3)1234 ω − ω − ω − ω ,B (4)1234 = T (4)1234 − ω − ω − ω − ω . (31) Appendix B Starting from the transformation in (22) and the gen- eralized Wick’s decomposition in (21), we obtain: m ( a ) k ( t ) = m ( b ) k ( t ) + 2 (cid:16) n ( b ) k ( t ) + n ( b ) − k ( t ) (cid:17) × (cid:88) j B (3) k, − k,j,j n ( b ) j ( t ) , (32)where higher order terms in m k have been neglected. Thenext step consists in building the evolution equation for m ( b ) k ( t ) from equation (23). Interestingly, the evolutionequation for m ( b ) k ( t ) appears as a deterministic dispersivenon homogeneous wave evolution equation [15], i dm ( b ) k dt = 2˜ ω k m ( b ) k + (cid:2) ( n ( b ) k + n ( b ) − k ) (cid:88) j T k − kj − j m ( b ) j (cid:3) , (33)with ˜ ω k = ω k + 2 (cid:80) j T kjkj n ( b )) j . Such equations havebeen derived in the theory of Bose Einstein condensatesand superconductivity.The equation for the spectrum, see [15], is given by dn ( b ) k dt = − m ( b ) k (cid:88) T k, − k,j, − j m ( b ) ∗ j ] . (34)From equations (33) and (34), after some algebra, it ispossible to show that the following interesting relationsholds: d [ | m ( b ) ( k ) | ] dt = d [ n ( b ) k n ( b ) − k ] dtd [ n ( b ) k − n ( b ) − k ] dt = 0 (35)1If n ( b ) k has reached energy equipartition such that n ( b ) k = n ( b ) − k = const/ω k then | m ( b ) k | = n ( b ) k ; therefore, we expectto observe equipartition also for ω k | m ( b ) k | .We now consider the leading order solution of equation(33) m ( b ) k ( t ) = m ( b ) k ( t ) e − i ω k t + higher order terms . and plug it in (32), and assuming that the spectrum n k is in stationary conditions, we get m ( a ) k ( t ) = m ( b ) k ( t ) e − i ω k t + 2 (cid:16) n ( a ) k ( t ) + n ( a ) − k ( t ) (cid:17) × (cid:88) j B (3) k, − k,j,j n ( a ) j ( t ) . (36)Note that we have used the fact that at the leading order n ( b ) k ( t ) (cid:39) n ( a ) k ( t ). Appendix C We consider equation eq. (22) and take the FourierTransform in time, to get: a k i , Ω p = b k i , Ω p ++ (cid:90) (cid:88) j,k,l B (1) ijkl b j, Ω q b k, Ω r b l, Ω s δ klij δ Ω q Ω r Ω s Ω p d Ω qrs ++ (cid:90) (cid:88) j,k,l B (3) ijkl b ∗ j, Ω q b ∗ k, Ω r b l, Ω s δ lijk δ Ω r Ω s Ω p Ω q d Ω qrs ++ (cid:90) (cid:88) j,k,l B (4) ijkl b ∗ j, Ω q b ∗ k, Ω r b ∗ l, Ω s δ ijkl δ Ω p Ω q Ω r Ω s d Ω qrs . (37)The next step is to build the second order correlator (cid:104) a ( k i , Ω l ) a ( k j , Ω m ) ∗ (cid:105) assuming stationarity.We use the generalized Wick’s decomposition in (21),i.e. including the anomalous correlators. The leadingorder result is contained in equation (25).2 [1] G. Falkovich, V. S. Lvov, and V. E. Zakharov, Kol-mogorov spectra of turbulence (Springer, Berlin, 1992).[2] S. Nazarenko, Wave turbulence , Vol. 825 (Springer,2011).[3] J. Benney and A. C. Newell, “Random wave closure,”Studies in Appl. Math. , 1 (1969).[4] A. C. Newell, “The closure problem in a system of ran-dom gravity waves,” Review of Geophysics , 1 (1968).[5] D. J. Benney and P. Saffmann, “Nonlinear interaction ofrandom waves in a dispersive medium,” Proc Royal. Soc , 301–320 (1966).[6] B. B. Kadomtsev, Plasma Turbulence (Academic Press,New York, 1965).[7] Y. Choi, Y. V. Lvov, S. Nazarenko, and B. Pokorni,“Anomalous probability of large amplitudes in wave tur-bulence,” Physics Letters A , 361 (2005).[8] Y. Choi, Y. V. Lvov, and S. Nazarenko, “Probabilitydensities and preservation of randomness in wave turbu-lence,” Physics Letters A , 230 (2004).[9] V. S. Lvov, Wave turbulence under parametric excita-tion: applications to magnets (Springer Science & Busi-ness Media, 2012).[10] V. E. Zakharov, V. S. Lvov, and S. S. Starobinets,“Spin-wave turbulence beyond the parametric excitationthreshold,” Physics-Uspekhi , 896–919 (1975).[11] P. Janssen, The interaction of ocean waves and wind (Cambridge University Press, Cambridge, 2004) p. 379.[12] A. Picozzi, J. Garnier, T. Hansson, P. Suret, S. Ran-doux, G. Millot, and D. N. Christodoulides, “Opticalwave turbulence: Towards a unified nonequilibrium ther-modynamic formulation of statistical nonlinear optics,”Physics Reports (2014).[13] D. Proment, S. Nazarenko, and M. Onorato, “Sus-tained turbulence in the three-dimensional grosspi-taevskii model,” Physica D: Nonlinear Phenomena (2012), 10.1016/j.physd.2011.06.007.[14] D. Proment M. Onorato, L. Vozella and Y. V. Lvov, “Aroute to thermalization in the α -fermi-pasta-ulam sys-tem,” Proceeding of National Academy of Science ,4208–4213 (2015).[15] V.S.Lvov, Wave Turbulence Under Parametric Excita-tions, Applications to Magnets (Springer-Verlag, 1994).[16] E. Fermi, J. Pasta, and S. Ulam, Studies of nonlinearproblems , Tech. Rep. (I, Los Alamos Scientific LaboratoryReport No. LA-1940, 1955).[17] L. Pistone, M. Onorato, and S. Chibbaro, “Thermaliza- tion in the discrete nonlinear klein-gordon chain in thewave-turbulence framework,” EPL (Europhysics Letters) , 44003 (2018).[18] Y. V. Lvov and M. Onorato, “Double scaling in the re-laxation time in the β -fermi-pasta-ulam-tsingou model,”Physical review letters , 144301 (2018).[19] W. Lee, G. Kovacic, and D. Cai, “Generation of disper-sion in nondispersive nonlinear waves in thermal equilib-rium,” Proceedings of the National Academy of Sciencesof the United States of America , 3237–3241 (2013).[20] V Zakharov, “Stability of period waves of finite ampli-tude on surface of a deep fluid,” J. Appl. Mech. Tech.Phys. , 190–194 (1968).[21] Miguel D Bustamante, Kevin Hutchinson, Yuri V Lvov,and Miguel Onorato, “Exact discrete resonances in thefermi-pasta-ulam–tsingou system,” Communications inNonlinear Science and Numerical Simulation , 437–471(2019).[22] H. Yoshida, “Construction of higher order symplectic in-tegrators,” Physics Letters A , 262–268 (1990).[23] M. Guasoni, J. Garnier, B. Rumpf, D. Sugny, J. Fatome,F. Amrani, G. Millot, and A. Picozzi, “Incoherent fermi-pasta-ulam recurrences and unconstrained thermaliza-tion mediated by strong phase correlations,” Physical Re-view X , 011025 (2017).[24] V P Krasitskii, “On reduced equations in the Hamilto-nian theory of weakly nonlinear surface waves,” J. FluidMech. , 1–20 (1994).[25] M. Onorato, S. Residori, U. Bortolozzo, A. Montina, andF. T. Arecchi, “Rogue waves and their generating mech-anisms in different physical contexts,” Physics Reports (2013), 10.1016/j.physrep.2013.03.001.[26] J. Bardeen, L. N. Cooper, J. R. Schrieffer “MicroscopicTheory of Superconductivity,” Physical Review ,162–164 (1957).[27] P. Miller, N. Vladimirova, F. Falkovich “Oscillations ina turbulence-condensate system,” Physical Review E ,065202 (2013).[28] S. Dyachenko, A.C. Newell, A. Pushkarev, V.E. Zakharov“Optical turbulence: weak turbulence, condensates andcollapsing filaments in the nonlinear Schrodinger equa-tion,” Physica D , 96–160 (1992).[29] N. Vladimirova, S. Derevyanko, F. Falkovich “Phasetransitions in wave turbulence,” Physical Review E85