Pre-asymptotic dynamics of the infinite size Neumann (p=2 spherical) model

TL;DR

This paper investigates the pre-asymptotic dynamics of the infinite size Neumann (p=2 spherical) model, providing detailed analysis and exploring the potential for describing the steady state using a Generalised Gibbs Ensemble.

Contribution

It extends previous work by analyzing the pre-asymptotic behavior and discussing the steady state description in the classical disordered p=2 spherical model.

Findings

Detailed characterization of pre-asymptotic dynamics

Discussion on the applicability of Generalised Gibbs Ensemble

Summary of asymptotic results from prior work

Abstract

In this contribution we further study the classical disordered p=2 spherical model with Hamiltonian dynamics, or in integrable systems terms, the Neumann model, in the infinite size limit. We summarise the asymptotic results that some of us presented in a recent publication, and we deepen the analysis of the pre-asymptotic dynamics. We also discuss the possible description of the asymptotic steady state with a Generalised Gibbs Ensemble.

Tables

Table 1. Table 1: This table summarises the asymptotic values of the Lagrange multiplier, the correlation with the initial condition and the infinitely long time delay of the self correlation in the three first columns. The last two columns show the parameter dependence of the two characteristic frequencies ω ± subscript 𝜔 plus-or-minus \omega_{\pm} .

	$z_f$	${\chi }_{\mathrm{st}}$	$q_0$	$q$	$m{\omega }_{-}^2=z_f-2J$	$m{\omega }_{+}^2=z_f+2J$
I	$T^{\prime }+J^2/T^{\prime }$	$1/T^{\prime }$	0	0	${\left(\sqrt{T^{\prime }}-J/\sqrt{T^{\prime }}\right)}^2$	${\left(\sqrt{T^{\prime }}+J/\sqrt{T^{\prime }}\right)}^2$
II	$2J$	$1/J$	0	0	0	$4J$
III	$2J$	$1/J$	$\ne 0$	$\ne 0$	0	$4J$

Equations

$$H_{\mathrm{pot}}[\{s_i\}]$$

$$[J_{ij}]=0,[J_{ij}^2]=\frac{J^2}{N},$$

$$H_{\mathrm{constr}}[\{s_i\},z]=\frac{z}{2}\left(\sum _{i=1}^N{s}_i^2-N\right)$$

$$H_{\mathrm{kin}}[\{{\dot{s}}_i\}]=\frac{m}{2}\sum _{i=1}^N({\dot{s}}_i{)}^2,$$

$$H_{\mathrm{syst}}=H_{\mathrm{kin}}+H_{\mathrm{pot}}+H_{\mathrm{constr}}.$$

$$m{\ddot{s}}_i(t)+z(t)s_i(t)=\sum _{j(\ne i)}{J}_{ij}{s}_j(t),$$

$$\underset{t^{\prime }\to t^{-}}{\lim }m{\partial }_{t^2}C(t,t^{\prime })+z(t)=-2e_{\mathrm{pot}}(t).$$

$$z(t)=-2e_{\mathrm{pot}}(t)+2e_{\mathrm{kin}}(t).$$

$$z(t)=2e_f-4e_{\mathrm{pot}}(t)=-2e_f+4e_{\mathrm{kin}}(t).$$

$$\underset{t\to \infty }{\lim }\underset{N\to \infty }{\lim }.$$

$$C_{ab}(t,t^{\prime })=\frac{1}{N}\sum _{i=1}^N[⟨s_i^a(t)s_i^b(t^{\prime }){⟩}_{i.c.}],$$

$$C_{ab}(t,0)=\frac{1}{N}\sum _{i=1}^N[⟨s_i^a(t)s_i^b(0){⟩}_{i.c.}],$$

$$R_{ab}(t,t^{\prime })={\frac{1}{N}\sum _{i=1}^N\frac{\delta [⟨s_i^a{⟩}_{i.c.}^{(h)}]}{\delta h^b(t^{\prime })}|}_{h=0},$$

$$[m{\partial }_t^2-z(t)]R(t,t_w)=\int d{t}^{\prime }\Sigma (t,t^{\prime })R(t^{\prime },t_w)+\delta (t-t_w),$$

$$[m{\partial }_t^2-z(t)]C(t,t_w)=\int d{t}^{\prime }\left[\Sigma (t,t^{\prime })C(t^{\prime },t_w)+D(t,t^{\prime })R(t_w,t^{\prime })\right]$$

$$+\frac{{\beta }^{\prime }{J}_0}{J}\sum _{a=1}^n{D}_a(t,0)C_a(t_w,0),$$

$$[m{\partial }_t^2-z(t)]C_a(t,0)=\int d{t}^{\prime }\Sigma (t,t^{\prime })C_a(t^{\prime },0)+\frac{{\beta }^{\prime }{J}_0}{J}\sum _{a=1}^n{D}_b(t,0)Q_{ab},$$

$$D(t,t_w)=J^{2}C(t,t_w),$$

$$D_a(t,0)=J^2{C}_a(t,0),$$

$$\Sigma (t,t_w)=J^{2}R(t,t_w).$$

$$z(t)=-m{\partial }_t^{2}C(t,t^{\prime }){|}_{t^{\prime }\to t^{-}}+2J^2\int d{t}^{\prime \prime }R(t,t^{\prime \prime })C(t,t^{\prime \prime })+\frac{J{J}_0}{T^{\prime }}\sum _a(C_a(t,0){)}^2.$$

$$q_{\mathrm{in}}=1-T^{\prime }/J_0.$$

$$C_1(t,0)\equiv C_{a=1b=1}(t,0)\text{mbox}and{C}_{b\ne 1}(t,0)\equiv C_{a=1b\ne 1}(t,0),$$

$$\displaystyle\begin{array}[]{ll}&C_{1}(t,0)=C(t,0)\qquad\mbox{implying}\qquad C_{1}(0,0)=C(0,0)=1\;,\\ &C_{b\neq 1}(0,0)=Q_{1,2}\qquad\mbox{implying}\qquad C_{b\neq 1}(0,0)=q_{\rm in}\;.\end{array}$$

$$C(t,t)=1,R(t,t)=0,$$

$${\partial }_{t}C(t,t^{\prime }){|}_{t^{\prime }\to t^{-}}={\partial }_{t}C(t,t^{\prime }){|}_{t^{\prime }\to t^{+}}=0,$$

$${\partial }_t{C}_{b\ne 1}(t,t^{\prime }){|}_{t^{\prime }\to t^{-}}={\partial }_t{C}_{b\ne 1}(t,t^{\prime }){|}_{t^{\prime }\to t^{+}}=0,$$

$${\partial }_t{C}_{b\ne 1}(t,0){|}_{t\to 0^{+}}=0,{\partial }_{t}R(t,t^{\prime }){|}_{t^{\prime }\to t^{-}}=\frac{1}{m},$$

$$\left(m{\partial }_t^2+z(t)\right)R(t,t^{\prime })$$

$$\left(m{\partial }_t^2+z(t)\right)C(t,t^{\prime })$$

$$\left(m{\partial }_t^2+z(t)\right)C_{b\ne 1}(t,0)$$

$$z(t)$$

$$\hat{R}(\omega )=\frac{1}{2J^2}\left[-m{\omega }^2+z_f\pm m\sqrt{({\omega }_{-}^2-{\omega }^2)({\omega }_{+}^2-{\omega }^2)}\right]$$

$$m{\omega }_{\pm }^2=z_f\pm 2J$$

$$q\equiv \underset{t-t^{\prime }\to \infty }{\lim }\underset{t^{\prime }\to \infty }{\lim }C(t,t^{\prime }),$$

$$q_0\equiv \underset{t\to \infty }{\lim }C(t,0).$$

$$\mathrm{env}\left[z(t)-z_f\right]\sim t^{-{\alpha }_z},$$

Full text

1. 1 Introduction

2. 2 The model

3. 2.1 The Hamiltonian

4. 2.2 The dynamic equations

5. 2.3 The initial conditions

6. 2.4 The quench

7. 3 Dynamics of the infinite size system

8. 3.1 Correlations, linear response and Lagrange multiplier

9. 3.2 Dynamic equations

10. 3.3 Initial conditions

11. 3.4 Replica structure in the course of time

12. 3.5 Equal time values

13. 3.6 The final set of equations

14. 3.7 Asymptotic analysis

15. 3.8 Phase diagram

16. 3.9 The various overlaps

17. 3.10 Approach to the asymptotic limit

18. 3.10.1 Pre-asymptotic behaviour of $z(t)$

19. 3.10.2 The susceptibility, and the kinetic and potential energies

20. 3.10.3 Correlation with the initial state

21. 3.10.4 The time-delayed self-correlation

22. 4 Mode analysis

23. 4.1 Exact solution

24. 4.2 Decay of the Lagrange multiplier

25. 5 Conclusions

26. A Correlation with the initial configuration

1 Introduction

Interest in the dynamics of quantum systems in perfect isolation has been re-boosted by the large activity in cold atom experiments [1] and, in parallel, the exact solution of one-dimensional models [2, 3, 4]. A huge amount of work has been performed in the quantum context and a host of results are now available for a variety of low-dimensional, mostly integrable, models. In particular, quantum quench protocols have been intensively used. Typically, these consist in following the dynamics of a state prepared with some prescription and evolved with a Hamiltonian of which it is not an eigenstate.

Questions on the equilibration, or not, of these systems are among the most prominent ones addressed in theoretical and experimental studies. Whether an equilibrium, or an equilibrium-like, description of the asymptotic dynamics of isolated systems exists, is a question that strictly makes sense in the limit of an infinite number of degrees of freedom only.

Very similar issues can be addressed in the classical, and not necessarily low-dimensional, context. Therefore, we very recently initiated a line of research aiming at clarifying the evolution of classical interacting isolated systems in the thermodynamic limit [5, 7]. The latter condition is fundamental for our purposes since we wish to analyse the possible thermal properties of the long-term behaviour of the systems after sudden changes in their parameters, the classical equivalent of a quantum quench

Focusing on such thermodynamic limit, two classes of systems have been distinguished: integrable and non-integrable. While for the former, no canonical equilibrium is expected because of the infinitely many integrals of motion, for the latter, interactions are expected to be sufficiently efficient to lead them to standard equilibration. These issues have been and still are actively studied in the context of quantum isolated systems, see e.g. [2, 3, 4] for reviews, but, surprisingly enough, they have not been equally addressed in the realm of classical mechanics.

In the context of classical disordered models, we have checked that there is equilibration of interacting disordered systems, by studying the quenched dynamics of the isolated $p\geq 3$-spin disordered model [5] (excluding glassy features for quenches to the so-called threshold [6]). In a subsequent publication, we started the study of the $p=2$ case [7], a model that, as we explain below, is integrable and therefore not expected to equilibrate to the conventional Gibbs-Boltzmann measure.

Our intention in [7] was to analyse the dynamics of a classical integrable system that is not simply mappable to an ensemble of independent harmonic oscillators. We were especially interested in the limit of an infinite number of degrees of freedom, $N\to\infty$. We wished to establish whether the long-term dynamics approach a steady state that could be described with a Generalised Gibbs Ensemble, a classical version of the ensembles used to characterise the asymptotic states of integrable quantum systems [2, 3, 4, 8, 9]. With this aim, we needed to fully understand the dynamics of the classical model in the long-times limit.

The model that we choose to study is very well-known in the field of disordered systems and it goes under the name of $p=2$ spherical model [10, 11, 12, 13, 14, 15, 15, 17, 18, 19, 20, 21, 22]. In the Hamiltonian version that we will study here, it is a well-known integrable model, that was proposed by Neumann in the 19th century [23] and studied by several authors [24, 25, 26] more recently. In Ref. [7] we studied its dynamics in great detail. In particular, we elucidated the dynamic phase diagram for quenches from initial configurations drawn from a canonical equilibrium measure at an inverse temperature $T^{\prime}$, after quenches that we specify below. This choice of initial conditions corresponds to taking the system in equilibrium with a bath at the chosen initial temperature, $T^{\prime}$, at times $t<0$, and then suddenly switching off the coupling at $t=0$. If no change in the Hamiltonian is done, the dynamics remains in equilibrium at $T^{\prime}$. If, instead, a change in the Hamiltonian is operated at $t>0$, the system is subject to a sudden quench and the subsequent dynamics are out of equilibrium.

In [7] we did not, though, present a detailed study of the way in which the asymptotic states are reached. We fill this gap here investigating the pre-asymptotic dynamics of the model for $N\to\infty$. We will elucidate the pre-asymptotic behaviour and distinguish it from the one of the dissipative problem that is not able to equilibrate for certain choices of parameters either, but for different reasons.

The article is organised as follows. In Sec. 2 we recall the definition of the model and we give some more details on the quench performed. Section 3 is devoted to the study of the pre-asymptotic dynamics. Finally, in Sec. 5 we present our conclusions and some discussion of lines for further research.

2 The model

In this Section we briefly present the definition of the model, the initial conditions that we choose, and the kind of quench that we perform. A more detailed summary of the behaviour of the model in canonical equilibrium and its relaxational dynamics driven by the coupling to a thermal bath can be found in [7].

2.1 The Hamiltonian

The potential energy of the $p=2$ spin model is the one of a system with two-spin interactions mediated by quenched random couplings $J_{ij}$:

[TABLE]

The exchanges $J_{ij}$ are independent identically distributed random variables drawn from a Gaussian probability distribution function with average and variance

[TABLE]

respectively. All couplings are symmetric under the exchange of their indices, $J_{ij}=J_{ji}$, and the parameter $J$ characterises the width of the Gaussian distribution. We use continuous “spin” variables, $-\infty\leq s_{i}\leq\infty$ with $i=1,\dots,N$, globally forced to satisfy (on average) a spherical constraint, $\sum_{i=1}^{N}s_{i}^{2}=N$, with $N$ the total number of spins [10], that is imposed by a term

[TABLE]

added to the Hamiltonian. $z$ is a Lagrange multiplier that depends on the system’s configuration. The spins thus defined do not have intrinsic dynamics. In statistical physics applications their temporal evolution is given by the coupling to a thermal bath, via a Langevin equation [10, 11, 12, 13, 14, 17, 18, 19, 20, 21, 22] or a Monte Carlo rule [15, 16].

The model is endowed with conservative dynamics by changing the “spin” interpretation into a “particle” one [5, 7]. In this way, one defines a kinetic energy [27, 28, 29]

[TABLE]

that is added to the potential and constraint terms to yield the total energy of the Hamiltonian spherical $p=2$-spin model

[TABLE]

The model represents now a particle constrained to move on an $N$-dimensional hyper-sphere with radius $\sqrt{N}$. Its position is given by the $N$-dimensional vector $\vec{s}=(s_{1},\dots,s_{N})$ and its velocity by another $N$-dimensional vector $\dot{\vec{s}}=(\dot{s}_{1},\dots,\dot{s}_{N})$ . The $N$ coordinates $s_{i}$ are globally constrained to lie, as a vector, on a hypersphere with radius $\sqrt{N}$. The velocity vector $\dot{\vec{s}}$ is, on average, perpendicular to $\vec{s}$, due to the spherical constraint. The parameter $m$ is the mass of the particle.

2.2 The dynamic equations

The $N$ equations of motion for the isolated system are

[TABLE]

where the Lagrange multiplier needs to be time-dependent to enforce the spherical constraint all along the evolution.

Equation (6) yields an identity between the energy density and the Lagrange multiplier. Its product with $s_{i}(t_{2})$ and the limit $t_{2}\to t_{1}^{-}$ imply

[TABLE]

The first term can be rewritten as $m\lim_{t^{\prime}\to t^{-}}\partial_{t^{2}}C(t,t^{\prime})=-m\lim_{t^{\prime}\to t^{-}}\sum_{i=1}^{N}\dot{s}_{i}(t)\dot{s}_{i}(t^{\prime})=-m\sum_{i=1}^{N}(\dot{s}_{i}(t))^{2}$, and one has

[TABLE]

The Lagrange multiplier takes the form of an action density, as a difference between kinetic and potential energy densities. Using now the conservation of the total energy, $e_{f}=e_{\rm pot}(t)+e_{\rm kin}(t)$, one obtains

[TABLE]

2.3 The initial conditions

The initial condition will be taken to be $\{s_{i}^{0},{\dot{s}}_{i}^{0}\}\equiv\{s_{i}(0),{\dot{s}}_{i}(0)\}$ and chosen in ways that we specify below. We will be interested in using equilibrium initial states drawn from a Gibbs-Boltzmann measure at different temperatures $T^{\prime}$.

2.4 The quench

The quench is performed in such a way that $J_{ij}^{0}\mapsto J_{ij}=J/J_{0}\,J_{ij}^{0}$ instantaneously. This means that each interaction strength is rescaled by a parameter $J/J_{0}$ so quickly that neither the position not the velocity of the particle are changed. It is not hard to prove that under these conditions energy is either injected ($J/J_{0}<1$) or extracted ($J/J_{0}>1$) from the system.

3 Dynamics of the infinite size system

In this paper we will only discuss the dynamics in the $N\to\infty$ limit. Several methods can be applied to derive, in this limit, exact integro-differential equations for only three observables, the Lagrange multiplier, the self-correlation and the linear response, that fully characterise the collective dynamics of the system. The long times limit will only be taken next. The relevant order of limits is, therefore,

[TABLE]

3.1 Correlations, linear response and Lagrange multiplier

The self correlation and linear response are defined as

[TABLE]

where the averages $\langle\dots\rangle_{i.c.}$ are taken over the probability distribution of the initial conditions. The $a$ and $b$ indices run over replica indices, since we will use this method to enforce the initial conditions. The infinitesimal perturbation $h$ is linearly coupled to the spin, $H\mapsto H-h\sum_{i}s_{i}$. The upper-script ${(h)}$ indicates that the configuration is measured after having applied the field $h$. Because of causality, the linear response is non-zero only for $t>t^{\prime}$. The dynamics is causal and, in consequence, the linear response is proportional to a Heaviside theta function $\theta(t-t^{\prime})$.

In the $N\to\infty$ limit, the choice of the interaction matrix $J_{ij}$ is irrelevant in the sense that any one will be a typical one and the dynamics will be statistically the same. Still, in the calculation that we use to derive the Schwinger-Dyson equations an average over quenched randomness is taken. The square brackets denote here and everywhere in the paper this average.

The Lagrange multiplier $z(t)$ is fixed by the condition $C_{aa}(t,t)=1$ and we have already explained that it equals twice the difference between kinetic and potential energy densities.

3.2 Dynamic equations

In the $N\to\infty$ limit exact causal Schwinger-Dyson equations for the conservative dynamics for initial conditions drawn from the Gibbs-Boltzmann probability measure with Hamiltonian $H_{0}$ and inverse temperatrue $\beta^{\prime}$ are derived. They determine the evolution of the time-delayed self-correlation and linear response

[TABLE]

where the replica indices are $a=1,\dots,n\to 0$ since we have used the replica method to deal with $e^{-\beta^{\prime}H_{0}}$ and fix the overlap replica matrix $Q_{ab}$ [30, 31, 32]. The two kernels are the self-energy $\Sigma$ and the vertex $D$ and for the $p=2$ potential they take the forms

[TABLE]

The causality of the linear response, and hence of the self-energy that is proportional to the linear response, imposed the upper limit of the integrals to be either $t^{\prime}$ or $t$. The lower limit is always the initial time that we call $t=0$. The equation for the Lagrange multiplier can be written as

[TABLE]

3.3 Initial conditions

We have chosen to work with initial conditions drawn from the canonical equilibrium probability distribution built with the Hamiltonian $H_{0}$, with the $p=2$ form, and variance of the coupling strengths proportional to $J^{2}_{0}$. The canonical equilibrium properties of the model [10] are such that there is a phase transition at $T_{c}=J_{0}$ between a disordered paramagnetic phase at high temperature and an ordered phase at low temperature. In the replica formalism the two phases are distinguished by the structure of the replica matrix $Q_{ab}$ [10]:

Disordered initial conditions $\;\;\;T^{\prime}>J_{0}$           $Q_{ab}=\delta_{ab}\;.$

Ordered initial conditions $\;\;\;\;\;\;\;T^{\prime}<J_{0}$           $Q_{ab}=(1-q_{\rm in})\delta_{ab}+q_{\rm in}\;.$

with $q_{\rm in}$ the Edwards-Anderson parameter that is equal to

[TABLE]

This model is particularly simple since it is solved by a replica symmetric Ansatz. The two cases listed above are encoded by the second form if one simply keeps in mind that the $T^{\prime}>T_{c}$ corresponds to setting $q_{\rm in}=0$.

The Lagrange multiplier varies with temperature in the disordered phase while it gets fixed to the highest eigenvalue of the interaction matrix at low temperature:

Disordered initial conditions $\;\;\;T^{\prime}>J_{0}$           $z=T^{\prime}+J_{0}^{2}/T^{\prime}\;.$

Ordered initial conditions $\;\;\;\;\;\;\;T^{\prime}<J_{0}$           $z=2J_{0}\;.$

Equilibrium below $T_{c}$ corresponds to the condensation of the system in the direction of the largest eigenvector of the random matrix [10]; for this reason we will refer to this case as the condensed one.

3.4 Replica structure in the course of time

Although we change the Hamiltonian from $H_{0}$ to $H$ at the initial time $t=0$, the dynamic equations keep the initial replica structure. There will be two kinds of correlations with the initial condition

[TABLE]

where we singled out the replica labeled one and we shortened the notation to keep only one subscript. The interpretation of the correlations $C_{1}(t,t^{\prime})$ and $C_{b\neq 1}(t,t^{\prime})$ is similar to the one of real replicas. The former is the self-correlation between one replica $\{s_{i}\}(t^{\prime})$ and the same one further evolved until a later time $t$, $\{s_{i}\}(t)$. For this reason, $C_{1}(t,t^{\prime})\mapsto C(t,t^{\prime})$. The latter is the correlation between replica labeled $b$, let us say $\{\sigma_{i}\}(t^{\prime})$, at time $t^{\prime}$ and the singled out replica one evolved until time $t$ and represented by $\{s_{i}\}(t)$. It turns out that only $C_{b}(t,0)$, the correlation with replica $b$ evaluated at the initial time, intervenes in the other equations. This is the only correlation between different replicas that appears in the calculations.

The border conditions on the correlation functions then are

[TABLE]

3.5 Equal time values

The equal-time conditions are

[TABLE]

for all times $t,t^{\prime}$ larger than or equal to $0^{+}$, when the dynamics start.

3.6 The final set of equations

The set of equations that fully determine the evolution of the system from an initial condition in canonical Boltzmann equilibrium at any temperature $T^{\prime}$ are

[TABLE]

High and low temperature initial states are distinguished by $q_{\rm in}=0$ for $T^{\prime}>J_{0}$, and $q_{\rm in}=1-T^{\prime}/J_{0}$ for $T^{\prime}<J_{0}$, respectively. The equation for $C(t,0)$ is just the one for $C(t,t^{\prime})$ evaluated at $t^{\prime}=0$ so we do not write it explicitly.

3.7 Asymptotic analysis

We have shown in [7] that the dynamics approach an asymptotic limit such that $z(t)$ reaches a constant $z(t)\to z_{f}$ and $R(t,t^{\prime})$ a stationary function of the two times, $R(t-t^{\prime})$. In this limit the response equation can be Fourier transformed and Eq. (31) determines its frequency dependence

[TABLE]

with

[TABLE]

(note the unusual choice of sign for the imaginary part that we adopted in [7]). The numerical solution of the full dynamical equations show that the relevant sign is the minus one. $\hat{R}(\omega)$, and also $R(t)$, are independent of the initial temperature for $T^{\prime}<T_{c}=J$ while they depend on temperature through $z_{f}$ for $T^{\prime}>T_{c}=J$.

In terms of the physical parameters, $\hat{R}(\omega)$ is real for $|-m\omega^{2}+z_{f}|>2J$. In two of the dynamical phases that we identified in Ref. [7], $z_{f}=2J$, and this implies that the lower characteristic frequency vanishes $\omega_{-}=0$ and the imaginary part of the linear response is gapless.

3.8 Phase diagram

The phase diagram in Fig. 1 has three phases that are distinguished by the asymptotic value of the Lagrange multiplier, $z_{f}$, the static susceptibility, $\chi_{\rm st}=\hat{R}(\omega=0)$, the long-time limits of the delayed self-overlap,

[TABLE]

and the overlap with the initial condition,

[TABLE]

The first four columns in Table 1 summarise the values that these quantities take in the three phases.

3.9 The various overlaps

In Fig. 2 we show the numerical value of the overlap with the initial condition, $q_{0}\equiv\lim_{t\to\infty}C(t,0)$, and we compare it to $\sqrt{q_{\rm in}q}$, where $q$ is the one of the self-overlap, $q\equiv\lim_{t-t^{\prime}\to\infty}\lim_{t^{\prime}\to\infty}C(t,t^{\prime})$, and $q_{\rm in}=1-T^{\prime}/J_{0}$ is the overlap of the initial condition. The data are plotted as a function of the parameter $J/J_{0}$ that controls the energy injection or extraction. The agreement of the numerical data points is very good for the three initial conditions considered, that are made explicit in the key (see App. A and Ref. [33]).

3.10 Approach to the asymptotic limit

In this Section we present the rest of the original results in this contribution. They concern the way in which the functions reach the asymptotic values $z_{f}$, $q_{0}$ and $q$, recalled in Table 1 and confronted in Fig. 2. We will find that the approach is always algebraic with oscillations, and we will study the exponent as well as the relevant frequencies involved in the periodic behaviour. The dependence of the exponents for $|z(t)-z_{f}|$, $|C(t,0)-q_{0}|$ and $|C(t,t^{\prime})-q|$ on the control parameters $T^{\prime}/J_{0}$ and $J/J_{0}$ is given in Figs. 5 (a) and (b), and Fig. 16, respectively.

3.10.1 Pre-asymptotic behaviour of $z(t)$

In [7] we established that the long-time limit of $z(t)$ is given by the two values in the second column in Table 1 at the two sides of the diagonal $T^{\prime}/J_{0}=J/J_{0}$ in the phase diagram. Let us now study the approach to this value.

We first use the algebraic decay proposal

[TABLE]

where env stands for “the envelope of”, with the values of $z_{f}$ in Table 1, to extract the exponent $\alpha_{z}$. In practice, we will use the upper envelope of the data that, in all cases excluding the critical quenches, coincides with the results obtained from using the lower envelope. The critical quench is special and we will discuss it in detail below.

In Fig. 3 we display four panels with parameters such that $T^{\prime}$ is fixed to $0.5=T^{\prime}<J_{0}=T^{0}_{c}$ and $J/J_{0}$ changes in such a way that: there is injection of energy (upper left), the parameters are on the critical line (upper right), and there is extraction of energy (the two lower panels). A detailed study of these and other cases with $T^{\prime}<J_{0}$ indicates that the exponent $\alpha_{z}$ varies on the left of $T^{\prime}=J$, it takes the value $\alpha_{z}=3/2$ at the dynamic transition at $T^{\prime}=J$, and it remains fixed to $3/2$ on the right of it.

The decay in the critical quench in panel (b) seems to be steeper than what is shown with the dashed line that represents the fit. The reason for this is that the upper ($z(t)-z_{f}>0$) and lower ($z(t)-z_{f}<0$) envelopes of the function to fit do not have the same decay, at least up to the times $t\approx 300$ that we are able to explore with the numerical computation. The dashed line shown in the plot is an “average” between the fits of the upper and lower envelopes over the available times. In fact, a fit to the upper envelope only yields $\alpha_{z}>1.5$, while a fit to the lower envelope yields $\alpha_{z}<1.5$ (the lower envelope is not visible in the plot). Much longer times would be needed to reach a regime in which the envelopes have the same algebraic decay during the computational times and thus make a clean fit in this case.

A summary of the $\alpha_{z}$ values for $T^{\prime}/J_{0}<1$ are shown in Fig. 5 (a) with red points.

In Fig. 4 we repeat this analysis for $1.25=T^{\prime}>J_{0}=T^{0}_{c}$ and we see that, also for this kind of initial states, the exponent $\alpha_{z}$ depends on the parameters on the left of $T^{\prime}=J_{0}$ while it sticks to the value $3/2$ at $T^{\prime}=J$ and beyond it. The behaviour of $\alpha_{z}$ in between $J=J_{0}$ and $T^{\prime}=J$ is harder to determine from the numerical data. The exponent seems to reach the value $3/2$ at $J=J_{0}$ and it either remains constant beyond this (no quench) value or, if it varies, it does very smoothly, increasing slightly beyond $3/2$, reaching a shallow maximum and then decreasing again towards $3/2$ at $T^{\prime}=J$. A deeper analysis of this intermediate regime is easier to perform using a larger $T^{\prime}/J_{0}$ value that makes the interval $J_{0}<J<T^{\prime}$ wider, and we discuss these cases below.

A summary of the $\alpha_{z}$ values for $T^{\prime}/J_{0}=1.25$ are shown in Fig. 5 (a) with blue points.

With the aim of determining the actual behaviour of $\alpha_{z}$ in the uncertain region, we repeated this analysis for much higher values of the initial temperature, $T^{\prime}/J_{0}=2,\,3$, in such a way that the interval between $J=J_{0}$ and $T^{\prime}=J$ widens. We found that for these $T^{\prime}/J_{0}$, the exponent $\alpha_{z}$ is equal to $3/2$ or very close to it also for $T^{\prime}<J$, that is to say, on the left of the critical line. It is possible that the change in behaviour occurs at $J=J_{0}$ (the line separating quenches with injection of energy from those with extraction) since we detect a very small deviation from $3/2$ for $J<J_{0}$. However, this is hard to establish beyond doubt. Although the error bars that we estimate are smaller for $J<J_{0}$ than for $J>J_{0}$, the actual values of the exponent are pretty close to $3/2$ even for $J<J_{0}$. The summary of our findings for $\alpha_{z}$ at $T^{\prime}/J_{0}=2$ are shown with red points and the ones at $T^{\prime}/J_{0}=3$ are shown with blue points in Fig. 6 (a),

It is quite clear that the algebraic envelope captures part of the approach of $z(t)$ to its asymptotic value $z_{f}$ but not everything: the curves in Figs. 3 and 4 show strong oscillations. We now study the periodicity of the approach to the asymptotic value by focusing on

[TABLE]

with $\alpha_{\small z}$ the exponent extracted from the analysis in Eq. (39) and $z_{f}$ fixed to the values given in Table 1. This equation is the Fourier transform of the pure oscillating part of $z(t)$, in the long time limit. It is clear that our data will be a bit corrupted by the fact that the numerical time interval is forcefully rather short, going from a minimal time $t_{\rm min}$ where we expect to have reached the preasymptotic state to the maximal time in the simulation.

Figure 7 shows the oscillating part of $z(t)$, that is $t^{\alpha_{\small z}}[z(t)-z_{f}]/z_{f}$, with the exponent $\alpha_{\small z}$ defined in Eq. (39) and fixed to the value $\alpha_{\small z}=1.500$ obtained from fitting the envelope. The parameters are $T^{\prime}=0.5\ J_{0}$,  $J=1.5\ J_{0}$, energy is extracted from the sample with this quench, and $z_{f}=2J$. The panels on the left display the time dependence in two time intervals: at relatively short times when the decay is still far from the asymptotic limit, and at sufficiently long times to have reached it. The right panel is the Fourier transform in Eq. (40) with a clear peak at the frequency

[TABLE]

This figure shows that the characteristic frequency of the linear response function, $\omega_{+}$, the only non-vanishing one for this set of parameters, also determines the periodic dependence of the approach to the constant long-time limit of $z_{f}$. Once we have found the characteristic frequency $2\omega_{+}$ from the Fourier analysis, we plotted the form $A\cos(2\omega_{+}t+\phi)$ together with the numerical data in the lower panel on the left. The agreement between the two is perfect within the scale of the figure, see Fig. 7 (a).

In the next figure, Fig. 8, we investigate a quench to the critical line. The oscillating part of $z(t)$ displayed in the left panel proves that this quantity has not yet reached a limit in which all the decaying time-dependence is captured by the power law. More precisely, we determined the exponent $\alpha_{z}$ from the upper envelope and the data show that the lower part of the time-dependence is still varying quite a lot. We believe that this is due to the very slow dynamics in a quench to the critical line and that, basically, for this reason the system has not yet reached the pre-asymptotic regime. Still, the Fourier analysis points towards the special role played by the frequency $2\omega_{+}$ and the appearance of a (small) peak at $\omega_{-}=0$, signalling that this frequency will detach from zero beyond criticality and will play a role in the pre-asymptotic behaviour of $z(t)$ as well. The comparison between the numerical data and the curve $A\cos(2\omega_{+}t+\phi)+B$ is shown in Fig. 8 and the agreement is rather good, apart from the remanent time dependence in the lower part of the curve, that is seen as a deviation between the blue (fit) and red (data).

Finally, in Figs. 9 and 10 we treat cases with large energy injection taking the system from a condensed state or the disordered state into the phase in which $z_{f}$ depends on the two control parameters, that is to say, beyond the critical diagonal to the left of it. On this side of the transition both frequencies $\omega_{-}$ and $\omega_{+}$ are non zero and they influence the pre-asymptotic behaviour. This is not so clear in Fig. 9 where the system is taken out of a condensed state, but it is quite remarkable in Fig. 10. In the case $T^{\prime}/J_{0}=0.5$ and $J/J_{0}=0.25$ the peak at $2\omega_{+}$ is negligible with respect to the one at $2\omega_{-}$ (but should be real, not an effect of numerical error), see Fig. 9. Instead, in the case $T^{\prime}/J_{0}=1.25$ and $J/J_{0}=0.75$ (see Fig. 10) the higher frequency $2\omega_{+}$ is clearly present.

We summarise the results in this section in the equation

[TABLE]

that condenses in a single equation the time dependence of the Lagrange multiplier in the first correction to the asymptotic constant value.

3.10.2 The susceptibility, and the kinetic and potential energies

The dynamic susceptibility also approaches a constant asymptotic value with an algebraic decay. The first two panels in Fig. 11 show examples of this behaviour for $T^{\prime}<J_{0}$. The black and red continuous curves follow the upper and lower envelopes. The dotted black line is a fit to the power law decay. For higher values of $J/J_{0}$ we still find $b=1.5$. The summary of exponent values found is given in panel (c) in the same figure both for $T^{\prime}/J_{0}<1$ and $T^{\prime}/J_{0}>1$.

The “oscillatory” part of $\chi(t)-\chi_{f}$, defined as $\chi_{\rm osc}(t)=(\chi(t)-\chi_{f})\cdot t^{\alpha_{\chi}}$ where $\alpha_{\chi}$ is the exponent obtained by fitting the envelope of $\chi(t)-\chi_{f}$ with a power law, as done in Fig. 11 (a) and (b), follows the same time dependence as the one of $z(t)$. More precisely,

[TABLE]

with $B_{-},B_{+},\varphi_{-},\varphi_{+}$ fitting parameters and $\omega_{\pm}=\sqrt{(z_{f}\pm 2J)/m}$, describes the data very accurately. As beforehand, when $z_{f}=2J$, $\omega_{-}=0$, and thus $\chi(t)$ reduces to a single harmonic. We also note that this form implies

[TABLE]

The kinetic and potential energies, $e_{\rm kin}(t)$ and $e_{\rm pot}(t)$, depend linearly on $z(t)$,

[TABLE]

with $e_{\rm tot}$ the conserved total energy density and, therefore, they have the same time dependence as $z(t)$ once the factor and sign are taken into account.

3.10.3 Correlation with the initial state

A very similar analysis to the one described in the previous Subsection can be applied to the correlation between the instantaneous configuration and the initial one. The envelope to be studied is in this case

[TABLE]

with $q_{0}=\lim_{t\to+\infty}C(t,0)$ being 0 for $T^{\prime}/J_{0}>1$ and different from zero for $T^{\prime}/J_{0}<1$ and $T^{\prime}<J$ (sector III, dashed in green in dynamic phase diagram displayed in Fig. 1). We have to note an important difference with the analysis of the pre-asymptotic behaviour of $z(t)$ performed in the previous section. While $z_{f}$ is continuous, $q_{0}$ is discontinuous at the transition $T^{\prime}=J$. We do not have an analytic expression for $q_{0}$ but we have clear numerical evidence for this claim, that was shown in Fig. 18 (b) in [7]. This implies that the behaviour on the left of the critical line may not be the same as the one on the right of it entailing, as we will see, a discontinuity of the exponent $\alpha_{C}$ itself.

For reference we show in Fig. 12 the behaviour of this correlation for $T^{\prime}/J_{0}=0.5$ and various choices of $J/J_{0}$. To notice here is the case of the critical quench that has been approached from the disordered side, and the asymptotic value $q_{0}=0$ was used. The exponent $\alpha_{C}$ takes a very small value $\alpha_{C}\simeq 0.226$. Right beyond criticality we see, in Fig. 13, that the exponent jumps to $\alpha_{C}\simeq 3/2$, the value it takes in the full trapped dynamics, with $q_{0}\neq 0$. There are also differences between the behaviour of the upper and lower envelopes in close to critical quenches due to the short times available. Very close to criticality the data approach zero from above, and only at the latest times reached numerically we see the oscillations go beyond 0 and the upper and lower envelopes approach similar algebraic behaviours. When we get away from criticality the amplitude of the oscillations around zero is larger at the same time scales and the upper and lower envelope power laws get closer to each other. We think that by exploring longer times we should be able to observe the two envelopes with the same algebraic decay and recover $\alpha_{C}\simeq 1.5$ in all these cases.

Interestingly enough, on the critical line $T^{\prime}=J$, both for $T^{\prime}<J_{0}$ and $T^{\prime}>J_{0}$, the oscillations are suppressed at sufficiently long times, and the correlation approaches $q_{0}$ from above, taking only positive values, see Fig. 12 (b) and the three panels in Fig. 14. Moreover, in all the fits, we took $q_{0}=0$ meaning that we were approaching the critical line from the side $T^{\prime}>J$.

The collection of exponent values $\alpha_{C}$ thus derived are shown in Fig. 5 (b) and Fig. 6 (b). A comparison with the exponent $\alpha_{z}$ shown in the panels (a) of the same figures proves that $\alpha_{C}$ has a similar behaviour, though not identical, to the one of $\alpha_{z}$. Both for $T^{\prime}/J_{0}<1$ and $T^{\prime}/J_{0}>1$ the two exponents are constant and coincide, $\alpha_{z}=\alpha_{C}=3/2$, beyond the transition at $T^{\prime}=J$. Instead, for both kinds of initial states $\alpha_{C}$ weakly depends on the parameters for $J/J_{0}<1$. It is hard to establish beyond any doubt whether $\alpha_{C}=\alpha_{z}$ for these parameters but the numerical data are compatible with this possibility. In the region in between $J/J_{0}=1$ and $T^{\prime}/J=1$ the determination of $\alpha_{C}$ is more difficult. The data do not exclude the possibility of $\alpha_{C}$ being constant and equal to $3/2$ in this region too but a very weak variation could also be hidden within the errorbars.

3.10.4 The time-delayed self-correlation

The time-delayed correlation $C(t_{1},t_{2})$ depends, in principle, on the two times $t_{1}$ and $t_{2}$. In the study in [7] we did not find ageing phenomena in this model. We want to investigate more deeply here the possible two-time dependencies and, in particular, evaluate whether an interrupted ageing scenario, with the two-time correlation decaying as

[TABLE]

and $q_{\infty}(t_{2})=\lim_{\tau\to+\infty}C(t_{2}+\tau,t_{2})$, is possible.

For initial conditions in the ordered phase $T^{\prime}/J_{0}<1$ or quenches from the disordered phase $T^{\prime}/J_{0}$ to parameters on the left of the critical line, that is to say the white zone in the phase diagram in Fig. 1 that we called sector I, we do not expect ageing nor interrupted ageing, that is to say, $a=0$. This is confirmed by our numerical studies (not shown). We can therefore immediately study the stationary dynamics in these cases.

For quenches into the orange sector of the phase diagram, that we called II, one needs to be more careful and check whether there is, or not, interrupted ageing. We can mention at least three reasons for this possibility. Firstly, it was for this kind of quench that we found ageing in the $p=3$ model [5]. Secondly, these kind of energy extraction quenches from disordered initial conditions are the ones that yield interrupted ageing in the isolated $O(N)$ model in the large $N$ limit [34, 35, 36, 37, 38, 39, 40, 41, 42, 43]. Thirdly, a current study that is being carried out on the dynamics of the isolated one dimensional scalar field theory after quenches of different kind displays non trivial two-time dependencies in various correlations after this kind of procedures [44].

In Fig. 15 we plot $C(t_{1},t_{2})$ as a function of $t_{1}-t_{2}$ for various $t_{2}$, say $t_{2}=0,0.01,0.05,0.1,0.5,1,5,10,50$ after quenches into the orange sector of the phase diagram, that is to say, for parameters such that $T^{\prime}/J_{0}>1$, $T^{\prime}<J$. We show the data in linear scale and for short $t_{1}-t_{2}$ only, so as to avoid the many oscillations hiding the behaviour we wish to highlight. In all cases the dynamics soon become stationary, with the $t_{2}$ dependence being lost after, say, $t_{2}=5$. We conclude that there is no relevant interrupted ageing in this problem.

For $t_{2}\gg 1$, that is to say, for times $t_{2}$ such that $C$ has become stationary, we can perform a similar analysis of the algebraic approach to the asymptotic limit $q$

[TABLE]

with $q=\lim_{t_{2}\to\infty}q_{\infty}(t_{2})=\lim_{t_{2}\to\infty}\lim_{\tau\to+\infty}C(t_{2}+\tau,t_{2})$.

Figure 16 summarises the results of the study of this envelope and its fits with a power law. Once again we find that for $T^{\prime}/J_{0}<1$ and $J>T^{\prime}$ (on the right of the critical line) the exponent fixes to $3/2$ while for $T^{\prime}/J_{0}>1$ the data suggest that the behaviour of $\alpha_{C}^{\prime\prime}$ is similar to the one of $\alpha_{C}$. That is, for $T^{\prime}/J_{0}>1$, $\alpha_{C}^{\prime\prime}$ is equal to $3/2$ for $J>J_{0}$. For $J<J_{0}$, $\alpha_{C}^{\prime\prime}$ seems to decrease but it still remains very close to $3/2$.

4 Mode analysis

The dynamics can also be addressed by writing the Langevin equations in the rotated basis constructed with the eigenvectors $\vec{v}_{\mu}$, $\mu=1,...,N$, of the post-quench matrix $J_{ij}$. Indeed, defining

[TABLE]

the $N$ rotated equations of motion read

[TABLE]

where $\lambda_{\mu}$ are the eigenvalues associated to the eigenvectors $\vec{v}_{\mu}$.

This set of equations has to be complemented with the initial conditions $s_{\mu}(0)$ and $\dot{s}_{\mu}(0)$. They are close to the equations for a parametric oscillator, with the difference that here the time-dependent frequency depends on the variables via the Lagrange multiplier. In this formulation the relation with Neumann’s integrable classical system [23] is explicit.

The uniform interaction quench corresponds to

[TABLE]

4.1 Exact solution

Equation (50) can be implicitly solved with the amplitude-phase Ansatz [45, 46, 47, 48] in which the projection of the spin configuration is written as

[TABLE]

and the one of the spin’s projection velocity as

[TABLE]

complemented with the equation ruling the evolution of $\Omega_{\mu}(t)$,

[TABLE]

with the initial condition

[TABLE]

Equations (52) and (53) are reminiscent of the general solution of the harmonic oscillator problem, here with a time-dependent “frequency” $\Omega_{\mu}(t)$. In the numerical calculations we choose the initial condition

[TABLE]

that ensures that $\Omega_{\mu}(t)$ is real at all times.

The spherical constraint fixes the time-dependence of the Lagrange multiplier that can be expressed in different ways. For example,

[TABLE]

In Ref. [7] we gave more details on the numerical solution of the problem set in these terms. For our purposes here, that is finding the preasymptotic approach towards the limit value of $z(t)$, the finite $N$ solution is not well adapted, since we observe that even for relatively large number of modes, say $N=1024$, the power law is blurred at times that are quite short, say, $t=50$. In the next Subsection we propose an approximation that allows us to take the large $N$ limit within this formalism.

4.2 Decay of the Lagrange multiplier

The formalism in the previous subsection is exact but still quite complicated and, as we explained, not so well adapted to the numerical analysis of preasymptotic effects. We can make a series of assumptions, that we will check are verified, to go farther into a quasi-analytic description of the problem.

We therefore propose that if $z(t)$ is sufficiently close to its stationary limit, i.e., if $z(t)\simeq z_{f}$, then

[TABLE]

where $\omega^{(f)}_{\mu}=\sqrt{(z_{f}-\lambda_{\mu})/m}$ and $\phi_{\mu}$ is an unspecified phase. This form for $p_{\mu}(t)$ is consistent with the one obtained differentiating $s_{\mu}(t)$ with respect to time, and this implies that the two phases $\phi_{\mu}$ should be the same. This Ansatz assumes that when $z(t)$ is close enough to its stationary value, all the modes decouple and behave like independent harmonic oscillators, something that is verified in our calculations and that we discussed in [7]. The amplitude of the $s_{\mu}$ oscillation is given by $\sqrt{2T_{\mu}/(z_{f}-\lambda_{\mu})}$, where $T_{\mu}$ is the $\mu$th mode temperature. Inserting this Ansatz in the equation for $z(t)$, Eq. (57), we obtain

[TABLE]

If we are only interested in the asymptotic behaviour and we assume that $\phi_{\mu}$ is approximately the same for all modes we can take $\phi_{\mu}=0$. In the long time limit the factors $\sin^{2}(\omega^{(f)}_{\mu}t)$ and $\cos^{2}(\omega^{(f)}_{\mu}t)$ average to $1/2$ and a constant limit is reached due to dephasing among the different mode contributions,

[TABLE]

which is just another way of expressing the spherical constraint, now very close to the one in canonical equilibrium in which $T_{\mu}=T$ for all $\mu$. This equation fixes the constant asymptotic value of the Lagrange multiplier $z_{f}$. Moreover, we showed in Sec. 5.7 in Ref. [7] that in the limit in which the problem can be thought of as one of independent harmonic oscillators, the mode temperatures are given by

[TABLE]

with $z_{\rm in}$ the value of the Lagrange multiplier in the initial condition. After some careful handling of the sum in Eq. (61) for finite $N$, that is described in Sec. 5.7 of this reference, one rewrites the constraint equation as

[TABLE]

that is to say, the canonical equilibrium equation at temperature $T^{\prime}$ for a model with variance of the disorder distribution equal to $J^{2}$. This result is verified numerically, see Table 1.

In search of the finite time corrections to this asymptotic limit, we then split Eq. (60) in two terms, one equal to the asymptotic limit $z_{f}$ and the other one equal to the finite time correction:

[TABLE]

The issue now is to compute the time dependence of the sum in the second term.

One can evaluate the sums in (63) and (64) in the $N\rightarrow\infty$ thermodynamic limit, introducing the density of eigenvalues,

[TABLE]

Concerning the sum in (63),

[TABLE]

the correct asymptotic limit in Sectors I and IV. We already know that this expression cannot be continued on the other side of the transition at $T^{\prime}=J$ and that the Lagrange multiplier freezes at $z_{f}=2J$ for $T^{\prime}<J$. In order to capture this behaviour one has to go beyond the approximation (62) for the mode temperatures and consider that the edge-mode one scales with $N$, similarly to what happens in the canonical equilibrium calculation with the projection of the spin on the edge mode eigenvector.

As regards the sum in the finite time correction $f(t)$ it becomes

[TABLE]

where $\omega^{(f)}(\lambda)=\sqrt{(z_{f}-\lambda)/m}$ and

[TABLE]

We have not found an analytical solution to this integral but we computed it numerically. We follow this route below.

In Fig. 17 (a) we show the decay of the function $f(t)$ evaluated from Eq. (67) for $T^{\prime}=1.25\,J_{0}$ and $J=0.75\,J_{0}$, a point which belongs to sector I of the phase diagram. The function $f(t)$ decays as a power law with an exponent $\alpha_{z}=1.43$, which is in agreement with the value $\alpha_{z}=1.45$ found using the mean field equations (see Fig. 4). We note that in this case there is no problem in taking the $N\to\infty$ limit in the sum in (64) since $z_{\rm in}>2J_{0}$ and $z_{f}>2J$. In panel (b) in the same figure we show the decay of $f(t)$ for $T^{\prime}=1.25\,J_{0}$ and $J=1.5\,J_{0}$, in Sector II (orange) in the phase diagram. $f(t)$ now decays as a power law with exponent $3/2$, that coincides with the result from the $N\to\infty$ equations that predicts $\alpha_{z}=3/2$ in this sector of the phase diagram. In this case, the continuum approximation of the sum is not fully justified as $z_{f}=2J$. Still the agreement between the algebraic decay of the numerical evaluation of Eq. (67) and the slope $-3/2$ found with the Schwinger-Dyson approach is excellent.

In Fig. 18 (a) we show the decay of the function $f(t)$ for $T^{\prime}=0.5\,J_{0}$ and $J=0.8\,J_{0}$, in sector III of the phase diagram. We can clearly observe that $f(t)$ decays as a power law with the exponent $3/2$, which is in agreement with the results from the mean-field equations that predict $\alpha_{z}=3/2$ in this sector of the phase diagram. This is another case in which the continuum limit should not hold since both $z_{0}=2J_{0}$ and $z_{f}=2J$.

In panel (b) we show $f(t)$ for $T^{\prime}=0.5\,J_{0}$ and $J=0.25\,J_{0}$, in sector IV of the phase diagram. The behaviour is not as clear in this case. One can try a power law with an exponent $\alpha_{z}=0.815$, the same value that was found for $N\rightarrow\infty$ (see Fig. 3). The power law fits well the quasi-analytic results up to time $t=100$, approximately. For times larger than $t\simeq 100$ there are deviations from the power law behaviour. Since in this case $z_{\mathrm{in}}=2J_{0}$, there is a singularity in $T(\lambda)$ at $\lambda=2J$. The deviation from the power law behaviour could be due to large numerical errors caused by this singularity.

Let us end this section with a series of comments on the behaviour of $f(t)$ in the various sectors of the phase diagram.

The case in which there is no good agreement between the power law behaviour predicted by the results from the Schwinger-Dyson $N\to\infty$ equations and the quasi-analytic results for finite $N$ is sector IV. In this sector $z_{f}>2J$ but $z_{in}=2J_{0}$ in the $N\to\infty$ limit.

We can end with a short discussion of the behaviour of $f(t)$, more precisely, of each factor in the integrand that contribute in different ways:

– The factor $\sqrt{(2J)^{2}-\lambda^{2}}$, originating in the density of eigenvalues, tends to cure the singularities independently of the value of $z_{f}$ and $z_{in}$.

– The factor $T(\lambda)$ diverges if $z_{in}=2J_{0}$ and $z_{f}>2J$ which is the case in the problematic sector IV. If $z_{in}=2J_{0}$ and $z_{f}=2J$ it seems from the results in sector III that this factor does not cause trouble, the zero in the numerator being "cancelled" with the zero in the denominator.

– The factor $(2\lambda-z_{f})/(z_{f}-\lambda)$ diverges in sectors II and III, but the behaviour of the integral is smooth. Seemingly, the divergence is cured by the factor $\sqrt{(2J)^{2}-\lambda^{2}}$ coming from the eigenvalue density.

5 Conclusions

This paper complements the study of the conservative dynamics of the $p=2$ disordered spherical model or, in other words, the Neumann integrable model, that we started in [7]. We focused here on the pre-asymptotic dynamics of the model in the infinite system limit; more precisely, on the way in which the approach to the constant asymptotic values of the Lagrange multiplier, and correlation functions, is achieved. We showed that in all phases the dynamics approach the constant values algebraically with functions that decay as power laws but also oscillate in time. The exponents take the same constant value $\alpha=3/2$ for $z(t)$, $\chi(t)$, $R(t)$, $C(t,0)$ and $C(t-t^{\prime})$ in all quenches such that $J>T^{\prime}$, that is to say, to the right of the critical diagonal $T^{\prime}=J$ in the phase diagram. Instead, they depend on the parameters on the other side of the diagonal $T^{\prime}=J$. Apart from the algebraic decay, the functions are harmonic, with the frequencies given by $\omega_{+}$ and $\omega_{-}$.

The relaxational dynamics of the spherical $p=2$ model and the O(N) field theory share many points in common. In the former, the out of equilibrium dynamics after subcritical quenches correspond to the progressive alignment of the $N$ dimensional vector that collects all the spins, $\vec{s}=(s_{1},\dots,s_{N})$ on the direction of the eigenvector associated to the largest eigenvalue of the interaction matrix $J_{ij}$ [13, 14]. In the latter, the mechanism is very similar, and the progressive condensation takes place on the vanishing wave vector (tendency to local order) [50, 49]. Both models coarsen after subcritical quenches and the growing length scales have the same time-dependence $t^{1/2}$. Moreover, the Lagrange multiplier imposing the spherical constraint, or the equivalent quantity imposing the self-consistent condition on the field modulus, approach their asymptotic values with the same $t^{-1}$ algebraic law.

One can then naturally ask whether the two models, endowed now with conservative dynamics, also share the same scaling properties. However, this is not the case. To start with, interrupted ageing properties after quenches to the critical dynamic phase transition of the quantum O(N) model [34, 35, 36, 37, 38, 39, 40, 41, 42, 43] were derived in [38]. In the model here considered, interrupted ageing in critical quenches, if present, does not last sufficiently long to be measurable and is therefore irrelevant.

Appendix A Correlation with the initial configuration

We recall here a way to derive the relation $q_{0}=\sqrt{q_{\rm in}q}$, obtained for models with a complex structure of metastable states, such as the $p$-spin spherical model with $p\geq 3$. This relation, as we showed numerically in the main text, also holds for the simpler $p=2$ model.

The analysis of $q_{0}$ enables us to draw a link with the TAP states. In fact, the long-term limit of the correlation function with the initial condition, Eq. (37) can also be written as

[TABLE]

where the set $\{m_{i}(0)\}$ and $\{m_{i}(\infty)\}$ are the local magnetisations of the TAP state in which the system is initialised and the one reached asymptotically, at $t\rightarrow+\infty$. Recalling that

[TABLE]

with

[TABLE]

it follows that

[TABLE]

If we finally assume that the angular variables of the initial and final TAP states are identical, and that the unique effect of the quench has been to rescale their “width” the previous formula boils down to

[TABLE]

This relation is checked, numerically, in Fig. 2 and the agreement is perfect, within our numerical accuracy.

As an example, the dissipative spherical $p$-spin model verifies this property for a quench with $T_{\rm in}\in[T_{s},T_{d}]$. In fact, the dynamical equations yield, in this case, the two limits

[TABLE]

where in the first equation we defined $q^{*}$ as the factor between square brackets and the

[TABLE]

Thus, injecting the second equation in (74) into the first one it follows that

[TABLE]

and then

[TABLE]

This last equation is exactly the one that determines the order parameter when the system is at equilibrium with $T=T_{\rm in}$. Therefore, $q^{*}=q_{\rm in}$ and $q_{0}=\sqrt{q_{\rm in}q}$. If after the quench the system stays in the same TAP state, only the order parameter $q$ changes.

Acknowledgements

We acknowledge financial support from ECOS-Sud A14E01, PICS 506691 (CNRS-CONICET Argentina) and NSF under Grant No. PHY11-25915. LFC thanks the KITP Santa Barbara for hospitality during part of the preparation of this work. She is a member of Institut Universitaire de France.