Relevance of sampling schemes in light of Ruelles linear response theory

Relev ance of sampling sc hemes in ligh t of Ruelle’s linear resp onse theory V alerio Lucarini ∗ † valerio.lucarini@zmaw.de , T obias Kuna † t.kuna@reading.ac.uk , Jero en W outers ∗ jeroen.wouters@zmaw.de , Da vide F aranda ∗ davide.faranda@zmaw.de April 20, 2022 ∗ Klimac ampus, Universit¨ at Hambur g, Hambur g, Germany † Dep artment of Mathematics and Statistics, University of R e ading, R e ading, UK Abstract W e reconsider the theory of the linear resp onse of non-equilibrium steady states to p erturbations. W e first sho w that by using a general functional decomp osition for space-time dependent forcings, we can define elementary susceptibilities that allow to construct the resp onse of the system to general p erturbations. Starting from the definition of SRB measure, w e then study the consequence of taking different sam- pling sc hemes for analysing the response of the system. W e sho w that only a sp ecific choice of the time horizon for ev aluating the response of the system to a general time-dep enden t p erturbation allows to obtain the form ula first presen ted by Ruelle. W e also discuss the sp ecial case of perio dic p erturbations, sho wing that when they are tak en into con- sideration the sampling can b e fine-tuned to make the definition of the correct time horizon immaterial. Finally , we discuss the implications of our results in terms of strategies for analyzing the outputs of n umer- ical exp erimen ts by providing a critical review of a formula prop osed b y Reick. 1 In tro duction The study of how the prop erties of general non-equilibrium statistical mec hanical systems change when considering a generic perturbation, usually related to v ariations either in the v alue of some in ternal param- eters or in the external forcing, is of great relev ance, b oth in purely 1 mathematical terms and with regards to applications to the natural and so cial sciences. Whereas in quasi-equilibrium statistical mechan- ics it is p ossible to link the resp onse of a system to p erturbations to its unforced fluctuations thanks to the fluctuation-dissipation theorem [9, 28], in the general non-equilibrium case it is not p ossible to frame rigorously an equiv alence b et ween internal fluctuations and forcings. A t a fundamen tal level, this is closely related to the fact that forced and dissipative systems feature a singular inv arian t measure. Whereas natural fluctuations of the system are restricted to the unstable mani- fold, b ecause, by definition, asymptotically there is no dynamics along the stable manifold, p erturbations will induce motions - of exp onen- tially decaying amplitude - out of the attractor with probability one, as discussed in, e.g., [23, 21, 13]. It is worth noting that Lorenz antic- ipated some of these ideas when studying the difference b et ween free and forced v ariability of the climate system [12]. This crucial difficulty inheren t to out-of-equilibrium systems is lifted if the external p erturba- tion is, rather artificially , everywhere tangen t to the unstable manifold, or if the system includes some sto chastic forcing, whic h smo oths out the resulting the inv ariant measure [10]. Recen tly , Ruelle [20, 21, 23] pav ed the wa y to the study of the re- sp onse of general non-equilibrium systems to perturbations b y presen t- ing rigorous results leading to the formulation of a resp onse theory for Axiom A dynamical systems [19], which p ossess a Sinai-Ruelle-Bow en in v ariant measure [26]. Given a measurable observ able of the system, the change in its exp ectation v alue due to an  -p erturbation in the flo w (or in the map, in the case of discrete dynamics) can b e written as a perturbative series of terms proportional to  n , where eac h term of the series can b e written as the exp ectation v alue of some w ell-defined observ able o ver the unp erturbed state. Ruelle’s formula is identical to Kub o’s classical formula [8] when a Hamiltonian system is considered [13]. Whereas Axiom A systems are mathematically non-generic, the applicabilit y of the Ruelle theory to a v ariet y of actual mo dels is sup- p orted b y the so-called chaotic hyp othesis [6], whic h states that systems with man y degrees of freedom b eha v e as if they w ere Axiom A systems when macroscopic statistical prop erties are considered. The chaotic h yp othesis has b een interpreted as the natural extension of the classic ergo dic hypothesis to non-Hamiltonian systems [5]. In the last decade great efforts hav e b een directed at extending and clarifying the degree of applicability of the resp onse theory for non-equilibrium systems along five main lines: • extension of the theory for more general classes of dynamical systems [4, 2]; • in tro duction of effective algorithms for computing the resp onse in dynamical systems with many degrees of freedom [1], in order to supp ort the numerical analyses pioneered b y [18, 3]; • in vestigation of the frequency-dependent resp onse - the suscep- tibilit y - for the linear and nonlinear cases, with the ensuing in- tro duction of a new theory of Kramers-Kronig relations and sum 2 rules for non-equilibrium systems [13, 25] supported b y n umerical exp erimen ts [14]; • study of the response to external p erturbations of non-equilib- rium systems undergoing sto c hastic dynamics [17, 27]; • use of Ruelle’s response theory to study the impact of adding sto c hastic forcing to otherwise deterministic systems [15]. In particular, the resp onse theory seems esp ecially promising for tack- ling notoriously complex problems suc h as those related to studying the resp onse of geophysical systems to p erturbations, whic h include the inv estigation of climate change; see discussions in [1, 14, 16]. In particular, in [16], it is discussed that by deriving from the linear sus- ceptibilit y the time-dep enden t Green function, it is p ossible to devise a strategy to compute climate change for a general observ able and for a general time-dependent pattern of forcing. Recently , resp onse theory is b ecoming of great interest also in so cial sciences such as economics [7]. When developing a response theory , there are tw o p ossible w ays to frame the temp oral impact of the additional p erturbation to the dynamics. Either one considers the impact at a given time t of a p erturbation affecting the system since a very distant past, or one considers the impact in the distant future of p erturbations starting at the present time. When deriving the response form ula, Ruelle takes the first approach and delivers the c orr e ct formula [20]. T aking a differen t p oin t of view and considering the sp ecific case of perio dic p erturbations – whic h, anyw a y , tell us the whole story about the resp onse by linearity –, Reick [18] derives a formula that is w ell suited for analyzing the output of numerical exp erimen ts [14, 16]. Giv en the great relev ance and increasing p opularit y in applications of the resp onse theory introduced by Ruelle, in this pap er, we recon- sider the theory of the linear resp onse of non-equilibrium steady states to p erturbations and try to bridge the theoretical deriv ations and the strategies for designing numerical exp erimen ts and analyzing efficien tly their outputs. In Section 2, w e study the relev ance of the c hoice of the time horizon for ev aluating the impact of the p erturbation and we demonstrate b y direct calculation that the Ruelle approach is the correct one. W e clarify some of the assumptions implicitly considered in his deriv ations. W e then discuss the sp ecial case of perio dic p erturbations, showing that using them as basis for a resp onse theory greatly simplifies the form ulas and the conditions under whic h the formulas are derived. In Section 3, we discuss the implications of our results in terms of strategies for impro ving the quality of numerical sim ulations and of the analysis of their output signals and reconsider Reic k’s form ula [18]. In Section 4 w e present our conclusions and p erspectives for future w ork. 3 2 Linear Resp onse Theory , revised 2.1 Separable p erturbations W e study the linear resp onse of a discrete dynamical system to general time-dep enden t perturbations. All calculations are formal, in the sense that we neglect all higher orders in the p erturbation without deriving an estimate for these terms and we assume that all sums con verge in all senses necessary . The unp erturbed dynamical system is giv en by x t +1 = f ( x t ) , with t ∈ Z , x t ∈ M , M b eing a smo oth manifold and f : M → M a differen tiable map. F or simplicit y we consider a time-indep enden t unp erturbed dynamics, although the following can b e extended to a time-dep enden t case in a s traigh tforw ard manner. Moreo ver, the anal- ysis of the case of a contin uous time flow ˙ x = f ( x ) is p erfectly analo- gous to what is presented in the following and the main corresp onding results will b e mentioned in App endix A. The dynamical system is perturb ed by a time-dep enden t forcing X ( t, x ) as follows: ˜ x t +1 = ˜ f t +1 ( ˜ x t ) := f ( ˜ x t ) + X ( t + 1 , f ( ˜ x t )) . (1) The effect of the perturbation on individual tra jectories is in general difficult to describ e. More can ho wev er be said ab out the statistical prop erties of the system. One can look at the exp ectation v alues of observ ables under inv arian t states of the dynamical system: ρ ( A ) := Z ρ ( dx ) A ( x ) , where ρ ( dx ) is an inv ariant measure of the unp erturbed dynamics i.e. ρ ( A ◦ f ) = ρ ( A ) . for an y observ able A. In general, a dynamical system can possess many in v ariant measures. The physically relev an t measure for dynamical systems is the SRB measure [26]. This measure is physical in the sense that for a set of initial conditions of full Lebesgue measure the time a verages lim t →∞ 1 t P t k =1 A ( f k ( x )) con verge to the exp ectation v alue under ρ . In other words, for an y measure l ( dx ) that is absolutely con tinuous w.r.t. Lebesgue, w e hav e that ρ ( A ) = lim t →∞ 1 t t X k =1 Z l ( dx ) A ( f k ( x )) . (2) W e wan t to determine the linear resp onse of exp ectation v alues under the SRB measure to p erturbations of the dynamical system as in Eq. 1. W e denote b y δ T ρ the difference in the expectation v alue b et w een the p erturb ed and unperturb ed system at time T . In [23], 4 Ruelle presen ts a formula for the linear resp onse due to p erturbations that are separable in time and space: X ( t, x ) = φ ( t ) χ ( x ) . The leading order term of the expansion of δ T ρ ( A ) in X is given by δ T ρ ( A ) ≈ X j ∈ Z G A ( j ) φ ( T − j ) , (3) with G A ( j ) = θ ( j ) Z ρ ( dx ) χ ( x ) D ( A ◦ f j )( x ) , (4) where θ is the Heaviside function. Since δ T ρ ( A ) is expressed as a con- v olution pro duct of G A and φ , the F ourier transform of the resp onse δ ω ρ ( A ) = P T ∈ Z e iT ω δ T ρ ( A ) is giv en b y a product of the F ourier trans- form ˆ φ ( ω ) of the time factor φ ( t ) and a susceptibility function ˆ κ A ( ω ): δ ω ρ ( A ) ≈ ˆ κ A ( ω ) ˆ φ ( ω ) . (5) where ˆ φ ( ω ) = X j ∈ Z e ij ω φ ( j ) , ˆ κ A ( ω ) = X j ∈ Z e ij ω G A ( j ) = X j ≥ 0 e ij ω Z ρ ( dx ) χ ( x ) D ( A ◦ f j )( x ) . (6) Due to the causality of the resp onse function G A ( j ) (i.e. G A ( j ) = 0, j < 0)), the susceptibilit y ˆ κ A ( ω ) is analytic in the upper complex plane and satisfies Kramers-Kronig relations [21, 13]. 2.2 General p erturbations If the p erturbation is of a more general nature (i.e. not separable), we can deduce a linear resp onse formula from Eq. 5, solely based on lin- earit y in the following wa y . Let φ r ( t ) b e a Schauder basis [11] of time- dep enden t functions and ψ s ( x ) a Sc hauder basis of space-dep enden t functions. One can take for example the F ourier basis in time and a wa velet basis in space, or whatever basis ma y be suitable for the system at hand. The pro duct functions φ r ( t ) ψ s ( x ) then form a basis of the time and space dependent functions as a tensor product [24]. More concretely , w e may for an appropriate sense of conv ergence as- sume that an y function X ( t, x ) can b e decomp osed in the pro duct basis φ r ( t ) ψ s ( x ) with co efficien ts a r,s : X ( t, x ) = X r,s ≥ 0 a r,s φ r ( t ) ψ s ( x ) . 5 Since each of the factors in this sum is separable, we can use Eq. 5 and the linearit y of the resp onse to get that the resp onse is given by δ ω ρ ( A ) ≈ X r,s ≥ 0 a r,s ˆ φ r ( ω ) ˆ κ s,A ( ω ) , (7) where ˆ κ s,a is the susceptibilit y function of observ able A , corresp onding to the forcing pattern given by ψ s ( x ). Since the v ectors ψ s ( x ) consti- tute a basis, the functions ˆ κ s,a are elementary linear susceptibilities that allow to construct the response of the system to any pattern of forcing. By inserting the expression of ˆ κ s,A ( ω ) from Eq. 6 in to Eq. 7, it is p ossible to deduce the frequency-dep enden t resp onse of the system. It is expressed as an ensemble a verage of a dot pro duct of F ourier transforms, namely the transforms of the p erturbation term and of the linear tangen t of the observ able, G A ( ω , x ): δ ω ρ ( A ) ≈ X j ≥ 0 e ij ω Z ρ ( dx ) ˆ X ( ω , x ) D ( A ◦ f j )( x ) = Z ρ ( dx ) ˆ X ( ω , x ) G A ( ω , x ) , (8) with G A ( ω , x ) = X j ≥ 0 e ij ω D ( A ◦ f j )( x ) ˆ X ( ω , x ) = X T ∈ Z X ( T , x ) e iω T = X r,s ≥ 0 a r,s ˆ φ r ( ω ) ψ s ( x ) . (9) Instead from Eqs. 3-4 in the time domain δ T ρ ( A ) ≈ Z ρ ( dx ) X j ≥ 0 X ( T − j, x ) D ( A ◦ f j )( x ) . (10) F or the case of a p erio dic p erturbation X ( t + τ , x ) = X ( t, x ), where τ ∈ N , we get as linear resp onse δ T ρ ( A ) ≈ Z ρ ( dx ) τ X n =1 ∞ X m =0 X ( T − n − mτ , x ) D ( A ◦ f n + mτ )( x ) = Z ρ ( dx ) τ X n =1 X ( T − n, x ) G A,n ( x ) , (11) with G A,n ( x ) = ∞ X m =0 D ( A ◦ f n + mτ )( x ) . 6 In order to elucidate some crucial asp ects of the Ruelle’s resp onse theory , we no w propose a direct deriv ation of the linear resp onse to the p erturbation X ( t, x ) b y considering the history of the perturb ed and unp erturb ed tra jectory of the system and verify under which con- ditions we find agreement with Eqs. 8-10. Our goal is to deriv e the leading order term of the expansion of δ T ρ ( A ) with resp ect to X from first principle, i.e. without resorting to the Schauder decomp osition as ab o v e. Such a deriv ation should of course arrive at the same results as those in Eqs. 8-10. 2.2.1 Resp onse at a mo ving time horizon W e describ e the p erturbed measure ˜ ρ T ( A ) such that the system is initialized at time T in an initial condition according to the measure l . W e mo v e the time horizon at whic h w e observe forward and a v erage the time-ev olved measuremen ts. The system is prepared and then observ ed while it is ev olving ov er a sufficiently long time. The measure ˜ ρ T is time-dep enden t as the dynamics ˜ f is also time-dep enden t. F ormally w e take ˜ ρ T to b e the ergo dic mean of the exp ectation v alues of A , starting at time T : ˜ ρ T ( A ) = lim t →∞ 1 t t X k =1 Z l ( dx ) A ( ˜ f k T ( x )) . (12) Here l ( dx ) is an initial measure that is absolutely contin uous with resp ect to Leb esgue and ˜ f k T represen ts k iterations of the p erturb ed dynamics from time T to T + k : ˜ f k T ( x ) = ˜ f T + k ◦ . . . ◦ ˜ f T +1 ( x ) . (13) The difference in exp ectation v alues δ T ρ is the given by δ T ρ ( A ) = ˜ ρ T ( A ) − ρ ( A ) . (14) F ollowing the computation presented in [23] for the separable case, w e can expand the p erturbed dynamics ˜ f around the unp erturb ed dynamics f . W e then try to rewrite the resp onse of the p erturbed system in terms of the SRB measure of the unp erturbed system b y finding an expression for A ( ˜ f k T ( x )) in terms of A ( f k ( x )). W e can approximate up to first order in X the tw o time step future ev olution by expanding around the unp erturbed dynamics f 2 ( x ): ˜ x T +2 = ˜ f T +2 ◦ ˜ f T +1 ( ˜ x T ) ≈ f 2 ( ˜ x T ) + X ( T + 1 , f ( ˜ x T )) .D f ( f ( ˜ x T )) + X ( T + 2 , f 2 ( ˜ x T )) . F or k time steps w e similarly get: ˜ x T + k = ˜ f T + k ◦ . . . ◦ ˜ f T +1 ( ˜ x T ) ≈ f k ( ˜ x T ) + k X j =1 X ( T + j, f j ( ˜ x T )) . ( D f k − j )( f j ( ˜ x T )) . 7 Th us, we can approximate A ( ˜ f T + k ◦ . . . ◦ ˜ f T +1 ( x )) to first order in X as follows: A ( ˜ f T + k ◦ . . . ◦ ˜ f T +1 ( x )) ≈ A ( f k ( x )) + A 0 ( f k ( x ))   k X j =1 X ( T + j, f j ( x )) . ( D f k − j )( f j ( x ))   = A ( f k ( x )) + k X j =1 X ( T + j, f j ( x )) D ( A ◦ f k − j )( f j ( x )) . (15) The linear response of A is obtained by substituting Eq. 15 in to Eq. 14, through Eq. 2 and Eq. 12: δ T ρ ( A ) ≈ lim t →∞ 1 t t X k =1 Z l ( dx ) k X j =1 X ( T + j, f j ( x )) D ( A ◦ f k − j )( f j ( x )) = lim t →∞ 1 t t − 1 X i =0 t − i X j =1 Z l ( dx ) X ( T + j, f j ( x )) D ( A ◦ f i )( f j ( x )) . Using that for i ≥ t the expression is zero, we ha ve δ T ρ ( A ) ≈ X i ≥ 0 Z   lim t →∞ 1 t t − i X j =1 ( f ∗ ) j l ( dx ) X ( T + j, x )   D ( A ◦ f i )( x ) . (16) Note that it is not p ossible to rewrite the sum in j as the ergo dic time mean of l due to the time dep endence of the p erturbation X ( T + j, x ). Therefore, surprisingly , Eq. 16 do es not in general agree with Eq. 10. In particular, by taking the limit on the righ t hand side, w e obtain that the T -dep endence disappears. Say we shift T to T − T 0 in the limit app earing in the ab ov e equation: lim t →∞ 1 t t − i X j =1 ( f ∗ ) j l ( dx ) X ( T − T 0 + j, x ) = lim t →∞ 1 t t − i − T 0 X j 0 =1 − T 0 ( f ∗ ) j 0 + T 0 l ( dx ) X ( T + j 0 , x ) = lim t →∞ 1 t t − i X j 0 =1 ( f ∗ ) j 0  ( f ∗ ) T 0 l ( dx )  X ( T + j 0 , x ) . T aking the reasonable assumption that the result in Eq. 16 do es not dep end on the initial m easure l ( dx ) (this cannot b e obtained from the uniqueness of the SRB measure), the obtained resp onse of the system is time-indep enden t even if the forcing is time-dep endent. 8 Let us compare the result con tained in Eq. 16 with Eq. 3 in the sp ecial case of a time-indep enden t p erturbation X ( t, x ) = χ ( x ). Now Eq. 16 and Eq. 3 agree since Eq. 16 simplifies to: δ T ρ ( A ) ≈ X i ≥ 0 Z ρ ( dx ) χ ( x ) D ( A ◦ f i )( x ) , b ecause lim t →∞ 1 /t P t − i j =1 ( f ∗ ) j l ( dx )) = ρ ( dx ), by the definition of the SRB measure. The form ula giv en by Ruelle [21] is recov ered, as can b e seen by substituting φ ( t ) = 1 into Eq. 3. Ho wev er, already in the case of a time-p eriodic p erturbation X ( t, x ) = X ( t + τ , x ) there is no agreement b et w een Eq. 16 and Eq. 10. In this case the sum o ver j app earing in Eq. 16 can b e written as a double sum, one o ver k p eriods, indexed by m , and one o ver the τ phases in each p erio d, indexed b y n : lim t →∞ 1 t t X j =1 ( f ∗ ) j l ( dx ) X ( T + j, x ) = lim k →∞ 1 k τ k X m =1 τ X n =1 ( f ∗ ) mτ ( f ∗ ) n l ( dx ) X ( T + n, x ) = 1 τ τ X n =1 ρ n ( dx ) X ( t + n, x ) , (17) with ρ n ( dx ) = lim k →∞ 1 k k X m =1 ( f ∗ ) mτ ( f ∗ ) n l ( dx ) . Under the assumption that ρ n = ρ for all n ∈ { 1 , . . . , τ } , i.e. sub- sampling do es not impact the unp erturbed inv ariant measure, the re- sp onse giv es a similar result as the Ruelle formula, but with an a v- eraged p erturbation. Substituting Eq. 17 into Eq. 16, w e obtain a form ula of the form of Eq. 10, with the difference that instead of the true forcing X ( t, x ) the av eraged forcing 1 τ τ X n =1 X ( t + n, x ) app ears. The disagreemen t is apparent, e.g. when one considers a p er- turbation of the form X ( t, x ) = sin( 2 π l τ t ) χ ( x ), whic h obviously results in a zero resp onse. This effect has a clear intuitiv e interpretation. The resp onse at a given time dep ends mostly on the immediate past, hence if one do es not k eep fixed the horizon, one risks to a verage out the v ariabilit y . The previous form ula reflects this intuiti on. One wa y to obtain agreement with F ormula 11 is to choose a specific sampling pro cedure. W e sample with the same p eriodicity τ of the 9 forcing, thus altering the definition of the response. W e define the measures for the p erturbed and unp erturbed system as ˜ ρ 0 T ,p ( dx ) := lim N →∞ 1 N N X k =0 ( ˜ f kτ + p T ) ∗ l ( dx ) ρ 0 ( dx ) := lim N →∞ 1 N N X k =0 ( f kτ + p ) ∗ l ( dx ) . (18) With this definition we obtain using Eq. 15: δ ρ 0 T ,p ( A ) = ˜ ρ 0 T ,p ( A ) − ρ 0 ( A ) ≈ lim N →∞ 1 N N X k =0 N − 1 X m = − 1 N − m X i =1 θ ( mτ + n + p ) ! Z  f mτ + n + p  ∗ l ( dx ) X ( T + mτ + n + p, x ) D ( A ◦ f kτ − mτ − n )( x ) . Using the p eriodicity of X one can obtain δ ρ 0 T ,p ( A ) ≈ τ X n =1 X i ≥ 1 Z lim N →∞ 1 N N − i X m = − 1 ( f n + mτ + p ) ∗ l ( dx ) θ ( mτ + p + n ) ! X ( T + n + p, x )  D ( A ◦ f iτ − n )( x )  = τ − 1 X n =0 Z ρ ( dx ) X ( T + p − n, x ) G A,n ( x ) = δ T + p ρ ( A ) . Hence b y choosing the initial phase p at which we start sampling, we can obtain the resp onse at this phase. This means that we only need to start one long simulation of f and ˜ f and do summations of the differences ( A ◦ ˜ f − A ◦ f )( x ) according to Eq. 18 at all phases p in one p erio d to obtain the en tire response to the p erio dic forcing. By applying a forcing that contains several frequencies, suc h as a block w av e, we can extract the susceptibilit y at all present frequencies in one run b y taking the F ourier transform of the resp onse. Note that if we sample the signal with a p erio dicit y η which is prime with resp ect to the perio d τ of the forcing, we will obtain no p -dep endence (with p , in this case, ranging from 0 to η − 1) in the resp onse. F or all v alues of p we will obtain as a result the resp onse to the time-av eraged forcing. Therefore, the case of sampling at all time steps discussed ab o ve is just the sp ecial case given b y η = 1, where we are basically considering the case of the Nyquist frequency . Instead, if τ and η are not prime with resp ect to each other, the sampling pro cedure will be able to ascertain the p -dep endence of the resp onse of the system at the p erio dicit y given by the common harmonic terms. If the p eriodicity of the forcing is not known, the ab o ve discussion tells us that by doing a sampling at larger and larger perio ds η and c hecking for each of those the phase-dep endence of the resp onse, it is p ossible to deduce the fundamental perio d of the forcing. If the pro ce- dure do es not con verge, we are facing a quasi-p eriodic or con tinuous- sp ectrum forcing for which this approach fails. 10 Therefore, this situation is unsatisfactory . Why do w e only get the correct result for p erio dic p erturbations and fine-tuning the sampling or b y taking constan t p erturbations? 2.2.2 Resp onse at a fixed time horizon This paradox can b e resolved by defining the time-dep enden t SRB measure in Eq. 12 using a different metho d of sampling. W e now consider the time evolution ˜ f k T in this definition to go from time T − k in the past up to the fixed time horizon T , so instead of Eq. 13, w e ha ve: ˜ f k T = ˜ f T ◦ . . . ◦ ˜ f T − k . (19) Note that this approach do es not use the reversed time dynamics but rather a different time p ersp ectiv e in which the final time is fixed as the curren t time and the p erturbation starts in the remote past. The expansion to first order in X around the dynamics of f now b ecomes: ˜ x T = ˜ f T ◦ . . . ◦ ˜ f T − k +1 ( ˜ x T − k ) ≈ f k ( ˜ x T − k ) + k − 1 X j =0 X ( T − j, f k − j ( ˜ x T − k )) . ( D f j )( f k − j ( x T − k )) . (20) Hence, the linear resp onse of ρ ( A ) at time T is given by δ T ρ ( A ) ≈ lim t →∞ 1 t t X k =1 Z l ( dx ) k − 1 X j =0 X ( T − j, f k − j ( x )) D ( A ◦ f j )( f k − j ( x )) = lim t →∞ 1 t X j ≥ 0 Z t − j X i =1 ( f ∗ ) i l ( dx ) X ( T − j, x ) D ( A ◦ f j )( x ) . Note that in contrast to Eq. 16 the indices are such that the time av er- age of the measure and the p erturbation are decoupled. This crucially dep ends on the choice of the sampling. This allows us to use the def- inition of the SRB measure in Eq. 2 and replace the time av erage in the limit by ρ : δ T ρ ( A ) ≈ X j ≥ 0 Z ρ ( dx ) X ( T − j, x ) D ( A ◦ f j )( x ) . (21) This agrees with Eq. 10. Here we do get the anticipated result. Note that this expression gives also a non-zero resp onse for a p erturbation whic h is non-zero only for a finite time, as opp osed to Eq. 16. This sampling is the natural one fore deducing the general linear resp onse theory . Doing the calculation for constan t forcing do es not elucidate the relev ance of the choice of sampling. This sampling cor- resp onds to a Gedankenexperiment where the system is prepared in the distan t past and we observe the difference of the p erturb ed and unp erturbed evolution up to a given instan t T . 11 3 Numerics 3.1 Reic k’s form ula F or p erturbations that are separable ( X ( t, x ) = φ ( t ) χ ( x )) and hav e a single driving frequency Ω ( φ ( t ) = cos (Ω t )), the follo wing sampling sc heme for computing the susceptibilit y for a given observ able A has b een prop osed by Reic k [18]: ˆ κ A (Ω) = lim  → 0 lim N →∞ 1 N  N X t =1 e i Ω t Z ρ ( dx )  A ( ˜ f t 0 ( x )) − A ( f t ( x ))  . (22) This formula has b een later adopted to analyze the output of a simple climate model [16] and a generalization has been proposed to study the nonlinear susceptibilities describing harmonic generation [14]. Appli- cabilit y of this form ula dep ends on p erforming numerical exp eriments where the initial samples approximate the unp erturb ed SRB measure ρ . Using our previous calculations, w e wan t to circumstan tiate the v a- lidit y of the formula. W e apply Ruelle’s resp onse theory to obtain a p erturbativ e expression of Eq. 22 in terms of quantities of the unp er- turb ed dynamics. lim N →∞ 1 N N X t =1 e i Ω t Z l ( dx )  A ( ˜ f t 0 ( x )) − A ( f t ( x ))  ≈ lim N →∞ 1 N N X t =1 e i Ω t Z t X j =1 ( f j ) ∗ l ( dx ) X ( j, f j ( x )) D ( A ◦ f t − j )( x ) . Here w e encounter the same problem as in Eq. 16, namely the coupling of the a verages of the measure and the p erturbation. Indeed, using Reic k’s form ula sampling from an initial measure differen t from the unp erturbed SRB measure, one do es not get a reasonable resp onse, as rep orted in [16]. By sampling according to the unp erturbed SRB measure ρ instead of l , the ab o ve equation b ecomes lim N →∞ 1 N N X j =1 N − j X k =0 e i Ω k e i Ω j Z ρ ( dx ) X ( j, x ) D ( A ◦ f k )( x ) = lim N →∞ 1 N Z ρ ( dx ) G A ( ω , x ) N X j =1 e i Ω j X ( j, x ) . (23) W e insert the inv erse discrete time F ourier transform X ( j, x ) = 1 2 π Z π − π ˆ X ( ω , x ) e − iω j dω 12 in to Eq. 23: Z ρ ( dx ) G A ( ω , x ) lim N →∞ 1 N N X j =1 e i Ω j X ( j, x ) = Z ρ ( dx ) G A ( ω , x ) lim N →∞ 1 N N X j =1 e i Ω j 1 2 π Z π − π ˆ X ( ω , x ) e − iω j dω = Z ρ ( dx ) G A ( ω , x ) lim N →∞ 1 2 π Z π − π ˆ X ( ω , x ) u N (Ω − ω ) dω = lim N →∞ 1 2 π Z π − π dω u N (Ω − ω ) ˆ κ A ( ω ) . (24) where u N (Ω − ω ) = 1 N N X j =1 e i (Ω − ω ) j This can b e rewritten by making use of x + . . . + x N = x 1 − x N 1 − x . as u N (Ω − ω ) = 1 N e i (Ω − ω ) 1 − e iN (Ω − ω ) 1 − e i (Ω − ω ) , whic h conv erges to 0 as N go es to infinit y , except for Ω = ω . A t Ω = ω the sum ov er j gives N . Hence, if ˆ X ( ω , x ) is integrable, we can take the limit in Eq. 24 inside the integral 1 2 π Z π − π ˆ X ( ω , x ) 1 { Ω } ( ω ) dω = 0 , where 1 { Ω } is the indicator function on { Ω } . W e deduce that in the case of a general p erturbation with a contin uous F ourier sp ectrum, Reick’s n umerical approach cannot b e applied. Note also that for finite time steps N (as is alwa ys the case for numerical exp eriments), there is an additional broadening of the signal of order 1 / N , as is apparent from Eq. 24. If on the other hand the F ourier transform is singular, for example X ( ω , x ) = δ (Ω − ω ) χ ( x ) , whic h corresp onds to the mono c hromatic signal X ( j, x ) = 1 2 π e − i Ω j χ ( x ) , w e hav e that lim N →∞ 1 N N X j =1 e i Ω j X ( j, x ) = χ ( x ) 2 π . 13 Therefore, Eq. 23 b ecomes Z ρ ( dx ) G A ( ω , x ) χ ( x ) = ˆ κ A (Ω) , as predicted by Reick. The ab o ve calculation shows how the explicit expansion of Reick’s resp onse formula allows us to in terpret its finite time b eha viour. Eq. 24 sho ws how this sampling scheme amoun ts to filtering the susceptibility ˆ κ A with the function u N . Note that in the case of a sev eral frequencies con tributing to the forcing, w e are again in the case of general p erio dic forcing. The sus- ceptibilit y can in this case b e computed in tw o w ays. Either one uses Reic k’s form ula at every frequency presen t in the signal, whic h amounts to doing sp ectroscop y . On the other hand, one can also compute the full resp onse δ T ρ at all phases ov er on p eriod and apply a F ourier trans- form to this time-dep endent function. The resp onse at any one sp ecific phase can b e efficiently computed with the p eriodic sampling strategy prop osed in Eq. 18. In this approach each v alue for the difference of A b et w een p erturb ed and unp erturb ed is pro cessed only once, compared to the summation b eing done for every frequency with Reick’s form ula. 3.2 Sampling contin uous sp ectra The discussion in the previous subsection demonstrates how sampling according to F orm ula 22 can only giv e a correct result in cases where the F ourier sp ectrum of the perturbation is discrete. In this subsection w e explore how the discussion on the expansion at a fixed time hori- zon can help us find a sampling for the case of a con tinuous F ourier sp ectrum. One p ossibilit y is to sample the resp onse directly from the full for- m ula of the p erturbation of the SRB measure: δ T ρ ( A ) = lim t →∞ 1 t t X k =1 Z l ( dx ) A ( ˜ f T ◦ . . . ◦ ˜ f T − k ( x )) − A ( f k ( x )) This sampling how ever entails some practical difficulties. As can b e seen from the formula, an ergo dic av erage is taken o ver the length of the n umerical run k . Increasing k to k + 1 is equiv alen t to altering the initial conditions from x to ˜ f T − k − 1 ( x ). This tra jectory cannot b e reco v ered from the previously calculated tra jectories of length k . Hence one needs to redo the calculations of the tra jectories for every v alue of k . As w e will see, less costly sampling methods can be devised. T o study the b eha viour of different sampling metho ds, let us define the follo wing quan tity: δ ( k,n ) ρ T ( A ) = Z l ( dx ) A ( ˜ f T ◦ . . . ◦ ˜ f T − k ◦ f n ( x )) − A ( f k + n ( x )) By changing k , w e con trol the length of time ov er whic h we observe the difference b et ween the p erturb ed and unp erturbed dynamics. The 14 initial measure l is furthermore transformed by n applications of the unp erturbed dynamics. By increasing n , the initial measure l conv erges to the unp erturbed SRB measure ρ . Making use of Equation 20, we can expand δ ( k,n ) ρ T to get a b etter idea of the behaviour of this quan tity under differen t limits and ergodic a verages: δ ( k,n ) ρ T ( A ) = k − 1 X j =0 Z f ( k − j + n ) ∗ l ( dx ) X ( T − j, x ) D ( A ◦ f j )( x ) (25) The aim when constructing a sampling sc heme is to take limits and ergo dic av erages ov er k and n in such a wa y that the resp onse given b y Equation 10 is obtained. As we ha ve seen with Reick’s form ula, this conv ergence can dep end on the p erturbation X . F urthermore, as exemplified by the discussion in this section, numerical cost should b e considered. The following ergo dic mean: lim t →∞ 1 t t X k =1 δ ( k,n ) ρ T ( A ) (26) con verges to the resp onse δ T ρ ( A ) given by Eq. 25 for any v alue of n . This can b e shown follo wing the discussion provided in Section 2.2.2, where w e discuss the case n = 0. Another p ossibilit y to obtain the SRB measure ρ in Eq. 25 is to take the limit of n going to infinity for a fixed k . W e obtain: lim n →∞ δ ( k,n ) ρ T ( A ) = k − 1 X j =0 Z ρ ( dx ) X ( T − j, x ) D ( A ◦ f j )( x ) . (27) where we ha v e assumed that lim m →∞ f m ∗ l = ρ . This expression tends to δ T ρ ( A ) in the limit of k → ∞ . F rom a theoretical p oin t of view, an increase in the v alue of n simply translates into a change in the initial measure l . Numerically , though, doing a long initial unp erturb ed run will evolv e the initial measure tow ards the inv ariant measure ρ , hence impro ving con vergence when Eqs. (25)-(27) are considered. In fact, there are a n umber of different options when attempting to reac h a go od numerical conv ergence. These include • increasing the length of unp erturb ed and perturb ed tra jectories ( n and k ) • enlarging the num b er of initial conditions (chosen according to l ) • deciding whether or not ergodic av eraging is performed ov er n and k . In the limit of infinitely long perturb ed runs, these approac hes give the same result. How ever, for finite time they will p erform differently . 4 Conclusions In this pap er we hav e reconsidered Ruelle’s linear resp onse theory b y analyzing the impact of choosing different metho ds of sampling in re- lation to differen t classes of forcings. Explicitly doing an expansion of 15 the perturb ed dynamics around the unperturb ed measures allo ws us to explore whic h sampling methods conv erge and under whic h conditions. The general resp onse form ula is obtained by c ho osing a sp ecific sampling where the system is prepared in the distant past and we ob- serv e the difference of the p erturbed and the unp erturbed dynamics up to a giv en time T . By prop osing a general decomposition of space-time dep enden t forcings using a Schauder decomposition, w e hav e elucidated that it is possible to define elementary linear susceptibilities that allow to construct the resp onse of the system to any pattern of forcing. The other p ossible sampling strategy , where the time horizon is not fixed, does not giv e rise to a natural response theory except for constan t p erturbations. In the case of p eriodic forcings one can obtain a meaningful formula by redefining appropriately the resp onse, finely tuned to the forcing under inv estigation. One needs to subsample the signal with the same p erio d of the forcing and explore all the initial phases. By taking this approach, it is in principle p ossible to discov er the fundamental p eriod of the external p erturbations b y v arying the sampling p eriod. Thanks to our approach w e get a deep er understand- ing of the range of applicabilit y of Reick’s formula, which has been used as a signal processing tool to study the linear resp onse of n umerical mo dels. Nonetheless, this approac h fails if the forcing is not p eriodic, in whic h case we must resort to the fixed-time horizon framework to get a meaningful answer. In fact, our findings explain why considering the fixed-time horizon it is p ossible to analyze a resp onse to forcings that hav e a contin uous F ourier sp ectrum. The clarifications presented in this pap er may b e of relev ance for devising the data pro cessing for actual lab oratory exp erimen ts on nonlinear systems. W e also clarify that it is crucial in practical terms to use an ensem- ble approac h where the initial conditions sample appro ximately the unp erturbed SRB measure. Our calculation is explicitly p erformed for discrete time, but the analogous results for contin uous time are presen ted in App endix A. Moreov er our considerations seem to b e ap- propriate also for the case of nonlinear resp onse [22, 14]. T o summarize, we hav e shown the following • Sampling a general resp onse from an initial time up to a moving time horizon do es not lead to a well-defined sampling metho d. • Starting the simulation at times in the distant past and a veraging the resp onse at a fixed time horizon alwa ys results in the full resp onse of the system at the fixed p oint in time. This approach can b e computationally inefficient. • Sampling a p erio dic resp onse with a moving time horizon results in a resp onse of the system as if it were forced with an av eraged forcing. • In the p eriodic case, the full resp onse can b e computed by sam- pling with a horizon mo ving forw ard in time with steps of one p eriod. This resp onse dep ends on the initial phase. The sus- ceptibilit y can b e computed through a F ourier transform of the 16 resp onse. • A constant forcing can b e considered as a p erio dic forcing with p eriod 1 and can thus b e sampled with a horizon moving with time steps of 1. • F or p eriodic forcings, Reick’s spectroscopic formula also allows to discern the resp onse at different frequencies, i.e. the suscepti- bilit y . It giv es a zero susceptibility for forcings with a contin uous sp ectrum. W e b eliev e that the results presented in this article can be of interest to researchers interested in studying the resp onse of complex systems to modulations of their internal parameters or to external perturba- tions. F or v arious reasons, climate science is an esp ecially promising field of application. First of all, we clarify crucial differences b et ween sampling p erio dic and ap eriodic forcings. This is a crucial issue if one w ants to apply linear resp onse theory to study different scenarios such as the resp onse of the system to monotonically increasing CO 2 lev els (see a forthcoming pap er by the authors) v ersus its response to perio dic forcings suc h as those due to astronomical and astrophysical phenom- ena. In particular, we hav e prop osed a parsimonious but effective w ay for analysing p eriodic - but non-mono c hromatic - forcings. Again, this setting is applicable to climate science due to the presence of cycles with different time scale, such as the daily , yearly and solar cycles. F uture work will address the inv estigation of the resp onse of a non- equilibrium system to a general random field. Moreo ver, we will ana- lyze the impact of the v arious sampling sc hemes described in this paper when studying the output of numerical mo dels. Ac knowledgemen ts VL, JW, DF ac knowledge the financial support of the EU-ER C pro ject NAMASTE-Thermo dynamics of the climate system. TK ackno wledges F. Bonetto for fruitful discussions. A Con tin uous time resp onse form ulas Here we give the formulas for contin uous time systems corresp onding to the ones presented in the main text. The time evolution is in this setting given by a differential equation dx dt = F ( x ) , resulting in a flow x ( t + s ) = f s ( x ( t )). The SRB measure is given b y ρ ( A ) = lim t →∞ 1 t Z t 0 ds Z l ( dx ) A ( f s ( x )) . F or a separable p erturbation dx dt = F ( x ) + χ ( x ) φ ( t ) 17 the susceptibility is given by ˆ κ A ( ω ) = Z ∞ 0 dte iω t Z ρ ( dx ) χ ( x ) D ( A ◦ f t )( x ) . In case of a general p erturbation F ( x ) → F ( x ) + X ( t, x ) the linear resp onse b ecomes: δ T ρ ( A ) ≈ Z ∞ 0 dτ Z ρ ( dx ) X ( T − τ , x ) D ( A ◦ f τ )( x ) The SRB measure with a moving time horizon is: ˜ ρ T = lim t →∞ 1 t Z t 0 ds Z l ( dx ) A ( ˜ f T + t T ( x )) and the SRB measure with a fixed time horizon: ˜ ρ T = lim t →∞ 1 t Z t 0 ds Z l ( dx ) A ( ˜ f T T − t ( x )) where ˜ f t 2 t 1 ( x ) is a tra jectory of the p erturb ed system, starting at time t 1 in x and evolving up to time t 2 . Reic k’s formula now b ecomes: ˆ κ A ( ω ) = lim  → 0 lim ν →∞ 1 ν  Z ν 0 dte i Ω t Z ρ ( dx )  A ( ˜ f t 0 ( x )) − A ( f t 0 ( x ))  18 References [1] R. Abramo v and A.J. Ma jda. Blended response algorithms for lin- ear fluctuation-dissipation for complex nonlinear dynamical sys- tems. Nonline arity , 20:2793–2821, 2007. [2] V. Baladi. On the susceptibility function of piecewise expanding in terv al maps. Comm. Math. Phys. , 275:839–859, 2007. [3] B. Cessac and J.-A. Sepulc hre. Linear response, susceptibility and resonances in chaotic toy mo dels. Physic a D , 225:13–28, 2007. [4] D. Dolgopy at. On differen tiability of SRB states for partially h yp erb olic systems. Invent. Math. , 155:389–449, 2004. [5] G. Galla votti. In J.P . F rancoise, G.L. Nab er, and T. S. Tsun, editors, None quilibrium Statistic al Me chanics (Stationary): Overview , v olume 3 of Encyclop e dia of Mathematic al Physics , pages 530–539. Elsevier Amsterdam, 2006. [6] G. Gallav otti and E.G.D. Cohen. Dynamical ensembles in sta- tionary states. J. Stat. Phys. , 80:931–970, 1995. [7] R.S. Hawkins and M. Aoki. Macro economic relaxation: Adjust- men t pro cesses of hierarchical economic structures. Ec onomics: The Op en-A c c ess, Op en-Assessment E-Journal , 3:1–21, 2009. [8] R. Kub o. Statistical-mec hanical theory of irreversible pro cesses. i. J. Phys. So c. Jpn. , 12:570–586, 1957. [9] R. Kub o. The fluctuation dissipation theorem. R ep. Pr o g. Phys. , 29:255–284, 1966. [10] G. Lacorata and A. V ulpiani. Fluctuation-resp onse relation and mo deling in systems with fast and slow dynamics. Nonlin. Pr o- c esses Ge ophys. , 14:681–694, 2007. [11] J. Lindenstrauss and L. Tzafriri. Classic al Banach sp ac es . Springer, 1996. [12] E.N. Lorenz. F orced and free v ariations of weather and climate. J. Atmos. Sci. , 36:1367–1376, 1979. [13] V. Lucarini. Resp onse theory for equilibrium and non-equilibrium statistical mechanics: Causalit y and generalized Kramers-Kronig relations. J. Stat. Phys. , 131:543–558, 2008. [14] V. Lucarini. Evidence of disp ersion relations for the nonlinear resp onse of Lorenz 63 system. J. Stat. Phys. , 134:381–400, 2009. [15] V. Lucarini. Sto c hastic p erturbations to dynamical systems: a resp onse theory approach. Physic al R eview E, submitte d , 2011. [16] V. Lucarini and S. Sarno. A statistical mechanical approach for the computation of the climatic resp onse to general forcings. Non- lin. Pr o c esses Ge ophys. , 18:7–28, 2011. [17] A.J. Ma jda and X. W ang. Linear resp onse theory for statistical ensem bles in complex systems with time-p eriodic forcing. Comm. Math. Sci. , 8:142–172, 2010. 19 [18] C. H. Reick. Linear resp onse of the Lorenz system. Physic al R eview E , 66(3):036103, 2002. [19] D. Ruelle. Chaotic evolution and str ange attr actors . Cambridge Univ ersity Press Cambridge, 1989. [20] D. Ruelle. Differentiation of SRB states. Communic ations in Mathematic al Physics , 187(1):227–241, 1997. [21] D. Ruelle. General linear resp onse form ula in statistical mec han- ics, and the fluctuation-dissipation theorem far from equilibrium. Phys. L etters A , 245:220–224, 1998. [22] D. Ruelle. Nonequilibrium statistical mec hanics near equilibrium: computing higher-order terms. Nonline arity , 11:5–18, 1998. [23] D. Ruelle. A review of linear resp onse theory for general differen- tiable dynamical systems. Nonline arity , 22(4):855–870, 2009. [24] R. A. Ryan. Intr o duction to tensor pr o ducts of Banach sp ac es . Springer, 2002. [25] A. Shimizu. Universal prop erties of nonlinear resp onse functions of nonequilibrium steady states . Journal of the Physic al So ciety of Jap an , 79:113001, 2010. [26] L.-S. Y oung. What are SRB measures, and which dynamical sys- tems hav e them? Journal of Statistic al Physics , 108:733–754, 2002. 10.1023/A:1019762724717. [27] T. Y uge. Sum rule for resp onse function in nonequilibrium Langevin systems. Phys. R ev. E , 82:051130, 2010. [28] D.N. Zubarev. None quilibrium Statistic al Thermo dynamics . Con- sultan t Bureau New Y ork, 1974. 20 Figure 1: The conv ergence of u N to the indicator function 1 { 0 } . The dashed blue line shows u 10 , the full red line u 50 . 21