Numerical dynamics of integrodifference equations: Periodic solutions and invariant manifolds in $C^\alpha(\Omega)$

Integrodifference equations are versatile models in theoretical ecology for the spatial dispersal of species evolving in non-overlapping generations. The dynamics of these infinite-dimensional discrete dynamical systems is often illustrated using computational simulations. This paper studies the effect of Nystr\"om discretization to the local dynamics of periodic integrodifference equations with H\"older continuous functions over a compact domain as state space. We prove persistence and convergence for hyperbolic periodic solutions and their associated stable and unstable manifolds respecting the convergence order of the quadrature/cubature method.


Growth and dispersal in discrete time
Difference equations of the form are commonly used to model the temporal evolution of single species which evolve in nonoverlapping generations, reproduce at specific time intervals or are censused at intervals (metered modells). The growth functions g t : R + → R + capture typical features of a particular population and are of Beverton-Holt, Ricker or Allee type, as well as related forms (see e.g. [22, pp. 11ff, Sect. 2.2]). On larger or inhomogeneous habitats also spatial effects such as dispersal have to be taken into account. This is B Christian Pötzsche christian.poetzsche@aau.at 1 Institut für Mathematik, Universität Klagenfurt, Universitätsstraße 65-67, 9020 Klagenfurt, Austria achieved in terms of a dispersal kernel k t (x, y) ≥ 0 indicating the probability of the species to move from position x to position y in the habitat Ω ⊂ R κ . Standard kernels are of Laplace-or Gauß-type, among others, and are discussed in for instance [22,pp. 17ff,Sect. 2.3]. This yields a recursion of the form u t+1 (x) = Ω k t (x, y)g t (u t (y)) dy for all x ∈ Ω, (1.1) which is denoted as integrodifference equation (IDE for short, cf. [18,22]). The state space variable u t (x) describes the population density in the t-th generation located at x ∈ Ω. The right-hand side of (1.1) maps the density in generation t to the density in the next generation t + 1 in two distinct stages: During the initial sedentary stage individuals grow, reproduce or die, and the local population density u t (x) evolves into g t (u t (x)). During the subsequent dispersal stage, the individuals move according to the dispersal kernel k t (x, y). Of course the above considerations extend to vector-valued growth functions g t and matrix-valued dispersal kernels k t in order to model several interacting or single but structured populations. In conclusion, difference equations such as (1.1) describe growth and dispersal, therefore, can be considered as a discrete time counterpart to reaction-diffusion equations, but with greater flexibility in the choice of the kernel. From a mathematical perspective, IDEs are infinite-dimensional discrete-time dynamical systems. Besides being popular tools in theoretical ecology over the recent years, they canonically arise as time discretizations of integrodifferential equations, as time-1-map of evolutionary partial differential equations or in the iterative solution of (nonlinear) boundary value problems [23, p. 190]. It is understood that IDEs involve an integral operator which in (1.1) is of Hammerstein-, but more general of Urysohn-type. Indeed, for our purposes a sufficiently flexible class are the recursions u t+1 (x) = Ω f t (x, y, u t (y)) dy for all x ∈ Ω, (I 0 ) whose natural state spaces consists of continuous or integrable functions over a compact subset Ω (the habitat in applications from ecology).

Numerical dynamics
In applied sciences the long-term behavior of IDEs is willingly illustrated using numerical simulations. For this purpose, [22, pp. 112-113] suggests to replace the integral in (I 0 ) by the trapezoidal or the Simpson rule. Both are special cases of general Nyström methods u t+1 (x) = η∈Ω n w η f t (x, η, u t (η)) for all x ∈ Ω (I n ) based on convergent quadrature or cubature rules with weights w η ≥ 0 and nodes η ∈ Ω n over a finite grid Ω n ⊂ Ω. Here, n ∈ N is related to the number of nodes in Ω n and therefore the accuracy of the approximation, see [8]. We point out that Nyström methods yield full discretizations of (I 0 ) and can be evaluated immediately.
While the numerical analysis of integral equations is a well-established field, e.g. [4,16], this paper enriches it by a dynamical aspect: We study and relate the long-term behavior of the iterates u t generated by an IDE (I 0 ) to those of a Nyström discretization (I n ). This brings us to the area of numerical dynamics [21,31,32] addressing the following questions: -Which dynamical or asymptotic properties of an IDE (I 0 ) as t → ∞ are preserved when passing to its Nyström discretizations (I n )? -What can be said about convergence as n → ∞ when the approximations become increasingly more accurate? In particular, are convergence rates of the integration rules preserved?
For the classical qualitative behavior of autonomous ODEs, such problems originate in [6] and are surveyed in [32]. In between various contributions to continuous-time infinite-dimensional dynamical systems generated by functional differential equations [15] or evolutionary (e.g. parabolic) partial differential equations [21,31] arose, both for spatial, as well as for full discretizations. IDEs merely require spatial discretization, but have in common with these problems that conventional error estimates are unsuitable to describe asymptotic behavior. In fact, bounds for the global discretization error typically grow exponentially in time and therefore establish convergence as n → ∞ only over compact time intervals [25,Thm. 4.1]. Thus, techniques extending those of standard numerical analysis are required to tackle the above problems. Previous contributions to the numerical dynamics of IDEs address basics and error estimates [25], as well as the persistence/convergence of globally asymptotically stable solutions [26]. This paper focusses on an another important aspect, namely the local saddle-point structure near periodic solutions to (I 0 ). Related work, but for autonomous evolutionary differential equations near equilibria, is due to [6,14] (ODEs), [1] (parabolic PDEs) and [12] (retarded FDEs).
In contrast, we study time-periodic IDEs (I 0 ) in the vicinity of periodic solutions. We stress that periodic time-dependence is strongly motivated by applications to incorporate seasonal influences. While [25,26] apply to semi-discretizations of (I 0 ) of collocation-or degenerate kernel-type [4,16], which act between finite-dimensional function spaces, but still contain integrals, we tackle Nyström discretizations (I n ), because they can be evaluated immediately. At this point the question for an ambient state space of (I 0 ) arises. A natural choice are the continuous functions C(Ω) over a compact Ω ⊂ R κ . Here however, already for linear integral operators, the discretization error under Nyström methods converges only in the strong, but not in the uniform topology as n → ∞, see [16,Lemma 4.7.6]. Using the theory of collectively compact operators [2] one can still establish that fixed-points of (I 0 ) (and their stability properties) persist [3,33]. Nonetheless, it is not clear how to establish convergence of the associated stable and unstable manifolds of (I n ) to those of the original problem (I 0 ). For this reason we retreat to the Hölder continuous functions C α (Ω) as state space. This set-up is sufficiently general to capture most relevant applied problems [18,22] and has the advantage that a more conventional perturbation theory (see App. A) applies to realize our goals. It should not be concealed, though, that the prize for this endeavor are more involved assumptions and technical preliminaries on Urysohn operators (well-definedness, complete continuity, differentiability). For the sake of a brief presentation they are outsourced to [27,28].
The structure of our presentation is as follows: In Sect. 2 we introduce a flexible framework for general periodic difference equations in Banach spaces and their linearization. Perturbation results for the Floquet spectrum of linear periodic equations are given in Theorem 2.1, while Theorem 2.2 addresses persistence and convergence of hyperbolic solutions and Theorem 2.3 the associated stable and unstable manifolds -when dealing with periodic equations we speak of fiber bundles. Although tailormade for Nyström discretizations of IDEs, these results also apply to collocationor degenerate kernel-discretizations, as well as when studying time-periodic evolutionary differential equations via their time-h-maps. The concrete case of Urysohn IDEs (I 0 ) is saved for Sect. 3 and illustrates how Thms. 2.1-2.3 apply. One obtains convergence of both hyperbolic solutions, and of the functions parametrizing their invariant fiber bundles with a rate given by the Hölder exponent α ∈ (0, 1] of the kernel functions f t in the first variable. Nonetheless, for smooth f t the higher-order convergence rates inherited from the particular quadrature/cubature rules are established. Section 4 contains numerical simulations confirming our theoretical results. An example with separable kernel (logistic IDE) allows explicit results and a direct comparison between exact with numerically obtained periodic solutions in the Hölder norm. The closing Sect. 5 comments on related approaches and simulation techniques. Finally, an "Appendix" summarizes the technical ingredients required in our analysis.
Notation We write R + := [0, ∞) for the nonnegative reals, S 1 := {z ∈ C : |z| = 1} for the unit circle in C, [·] : R → Z is the integer function and |·| denotes norms on finite-dimensional spaces. On the Cartesian product X × Y of normed spaces X , Y , is the product norm and we proceed accordingly on products of more than two spaces. The open resp. closed balls in X with radius r ≥ 0 and center x ∈ X are The bounded k-linear maps from the Cartesian product X k to Y are denoted by L k (X , Y ), L(X , Y ) := L 1 (X , Y ) and L 0 (X , Y ) := Y . Moreover, we abbreviate L k (X ) := L k (X , X ), L(X ) := L(X , X ), G L(X ) are the invertible maps in L(X ) and I X is the identity on X . Furthermore, N (T ) is the kernel and R(T ) the range of T ∈ L(X , Y ); σ (S) is the spectrum and σ p (S) the point spectrum of S ∈ L(X ).

Difference equations and perturbation
Let (X , · ) denote a Banach space.

Periodic difference equations
Abstractly, we are interested in a family of nonautonomous difference equations with right-hand sides F n t : U t → X on open sets U t ⊆ X , t ∈ Z, parametrized by n ∈ N 0 . In the following, n ∈ N is a discretization parameter such that F n t are understood as approximations converging to the original problem F 0 t as n → ∞ in a sense to be defined below. A nonautonomous set resp., holds. Given an initial time τ ∈ Z, a forward solution to (Δ n ) is a sequence φ = (φ t ) τ ≤t satisfying φ t ∈ U t and the solution identity φ t+1 = F n t (φ t ) for all τ ≤ t, a backward solution fulfills the solution identity for t < τ and for an entire solution The forward solution starting at τ in the initial state u τ ∈ U τ is uniquely determined as composition and denoted as the general solution to (Δ n ); it is defined as long as the compositions stay in U t . A difference equation (Δ n ) is called θ 0 -periodic, if both F n t+θ 0 = F n t and U t+θ 0 = U t hold for all t ∈ Z with some basic period θ 0 ∈ N; an autonomous equation is 1-periodic, i.e. there exists a F n : U → X with F n t ≡ F n , U t ≡ U on Z. A θ 1 -periodic solution to (Δ n ) is an entire solution satisfying φ t ≡ φ t+θ 1 on Z.
Given a fixed θ ∈ N and a sequence u = (u t ) t∈Z with u t ∈ U t , t ∈ Z, let us introduce the open productÛ := U 0 × . . . × U θ−1 and In order to characterize and compute periodic solutions to (Δ n ), n ∈ N 0 , we introduce the nonlinear operatorŝ and use the norm induced via (1.2) on the Cartesian product X θ . The next two results are immediate: Lemma 2.1 Let n ∈ N 0 , (Δ n ) be θ 0 -periodic and θ be a multiple of θ 0 : This characterization of periodic solutions to (Δ n ) via the mappingF has the numerical advantage to avoid the computation of compositions ϕ n (θ +τ ; τ, ·) : U τ → X , τ ∈ Z, and therefore preserves (numerical) backward stability (see [13]).

Lemma 2.2
Let n ∈ N 0 , m ∈ N, (Δ n ) be θ 0 -periodic and θ be a multiple of θ 0 . If every F n t : U t → X, 0 ≤ t < θ 0 , is m-times continuously (Fréchet) differentiable, thenF n :Û → X θ is of class C m and for everyû ∈Û one has

Linear periodic difference equations
Suppose that K n t ∈ L(X ), t ∈ Z, and consider a family of linear difference equations in X parametrized by n ∈ N 0 . As above we understand (L n ), n ∈ N, as perturbations of an initial problem (L 0 ). The transition operator Φ n (t, τ ) ∈ L(X ) of (L n ) is We are interested in θ -periodic equations (L n ), that is allowing us to introduce the period operator Ξ n θ := Φ n (θ, 0) ∈ L(X ) of (L n ). Its eigenvalues are the Floquet multipliers and σ p (Ξ n θ ) is the Floquet spectrum of (L n ). One says a linear difference equation (L n ) is weakly hyperbolic, if 1 / ∈ σ (Ξ n θ ), and hyperbolic, if S 1 ∩ σ (Ξ n θ ) = ∅ holds. In the hyperbolic situation, the spectrum can be decomposed as σ (Ξ n τ ) = σ +∪ σ − with spectral sets With the spectral projections P n , first for t ≥ 0 and then by θ -periodic continuation on Z. This yields θ -periodic nonautonomous have a constant finite dimension, which is denoted as the Morse index of (L n ) and equals the finite sum λ∈σ − dim N (λI X − Ξ n θ ) ι(λ) of algebraic multiplicities. We begin with a perturbation result for hyperbolic linear systems (L n ) under uniform convergence: Then also the period operator Ξ 0 θ ∈ L(X ) of (L 0 ) is compact and there exists a N ∈ N such that the following holds for all n ≥ N or n = 0: Proof Let 0 ≤ s < θ. Due to (i) the sequence K n s − K 0 s n∈N is bounded and consequently we obtain from K n s ≤ K 0 s + K n s − K 0 s and the periodicity condition We proceed by mathematical induction. Thanks to (2.2), for t = 0 the assertion is trivial and for t = 1 it results from (i). In the induction step t → t + 1 we obtain from the induction hypothesis and the triangle inequality, yielding the claim.
Hence, Ξ 0 θ is the uniform limit of by (ii) compact operators Ξ n θ , n ∈ N, and consequently compact [19, p. 416 Since the closed S and the compact σ (Ξ 0 θ ) are disjoint they have a positive distance. Therefore, there is an ε > 0 so that S ∩ B ε (σ (Ξ 0 θ )) = ∅. By the upper semicontinuity of the spectrum [5, p. 80, Lemma 3] and relation (2.5) there is a The hyperbolicity of (L n ) results as above in (a) with S = S 1 . Furthermore, then [30,p. 44,Prop. 3.13] implies that (L n ) possess an exponential dichotomy on Z as claimed with the θ -periodic invariant projectors P n t satisfying and I X − P n 0 = P n − . In particular, by (2.5) and [5, p. 80, Cor. 1] the spectral projections P n − associated to the unstable spectral parts of Ξ n θ , n ∈ N 0 , fulfill that dim R(P n − ) = dim R(P 0 − ) for large n, say for n ≥ n 2 . Thanks to (2.6) this extends to the dimension of the unstable bundles V n − . Finally, we set N : Lemma 4] yields that the spectral projections satisfy lim n→∞ P n − − P 0 − = 0. Together with claim (I) we obtain for t ∈ Z that from the triangle inequality, whose right-hand side converges to 0 as n → ∞.

Perturbation of hyperbolic solutions and invariant bundles
We next address the robustness of θ 1 -periodic solutions φ 0 to general θ 0 -periodic difference equations (Δ 0 ), as well as their nearby saddle-point structure consisting of stable and unstable bundles (see [17,pp. 143ff,Chap. 6], [24, pp. 256ff, Sect. 4.6]) under perturbation. By imposing a natural hyperbolicity condition on the solution φ 0 it is shown that also the perturbations (Δ n ) have (locally unique) periodic solutions φ n for sufficiently large n, which converge to φ 0 in the limit n → ∞.
Let θ := lcm {θ 0 , θ 1 }. We suppose that the right-hand sides F n t of (Δ n ) are continuously differentiable. Our endeavor is based on the variational equations associated to θ -periodic solutions φ n of (Δ n ), n ∈ N 0 . Since the linear equations (V n ) are θ -periodic, the terminology and results from Sect. 2.2 apply to (V n ) with K n t = DF n t (φ n t ) and the period operator Ξ n θ , n ∈ N 0 . In this context, we understand a solution φ n of (Δ n ) as (weakly) hyperbolic, if (V n ) has the corresponding property.
Based on this result, the (Floquet) spectrum of (V n ) can be computed from the (point) spectrum of the cyclic block operator given in Lemma 2.2. This has the numerical advantage of avoiding to evaluate the compositions (matrix products) Ξ n θ .

Proof
Keeping n ∈ N 0 fixed, we abbreviate K t = DF n t (φ n t ), t ∈ Z, and observe that the θ th power of DF(φ n ) given in Lemma 2.2 becomes a block diagonal operator Referring to [30, p. 42 . Now the Spectral Mapping Theorem [5, p. 65, Thm. 2] yields the assertion for the spectra. Concerning the point spectrum the claim follows directly from the corresponding eigenvalue-eigenvector relations and the solution identity for (V n ).
Our next result establishes persistence of hyperbolic periodic solutions to (Δ 0 ): , n ∈ N, are uniformly continuous on bounded sets uniformly in n ∈ N, the family DF n t n∈N is equicontinuous for all 0 ≤ t < θ 0 and for every n ∈ N there exists a 0 ≤ t < θ 0 such that DF n t has compact values.
If φ 0 is a weakly hyperbolic θ 1 -periodic solution to (Δ 0 ) and there exists a function then there exist reals ρ 0 > 0 and N 0 ∈ N such that the following hold for all n ≥ N 0 : As the subsequent proof and Lemma 2.3 reveal, the constant K 0 ≥ 0 essentially depends on the distance of the Floquet spectrum of φ 0 to the point 1 ∈ C. The value of K 0 blows up as this distance shrinks to 0, i.e. when (weak) hyperbolicity is lost.
(a) Because the assumptions (i'-iii') of Theorem A.1 hold, we can choose ρ, δ > 0 so small that (A.3) holds for e.g. q := 1 2 . Moreover, there exists a unique fixed point functionφ : We establish that the solutions φ n are weakly hyperbolic. For this purpose, let ε > 0. First, thanks to (2.8) there exists a n 1 ∈ N such that We know from Theorem A.1(c) that lim n→∞ sup t∈Z φ n t − φ 0 t = 0 and since DF n t is equicontinuous by assumption (ii), there exists a n 2 ∈ N such that Combining the last two inequalities readily yields DF 0 t (φ 0 t ) − DF n t (φ n t ) < ε for all t ∈ Z and n ≥ max {n 1 , n 2 }, which establishes the limit relation Second, assumption (ii) implies that the period operators Ξ n θ of (V n ), n ∈ N, contain a compact factor and hence are compact [19,p. 417,Thm. 1.2]. Thus, Theorem 2.1(a) applies to K n t := DF n t (φ n t ), n ∈ N 0 , and shows that φ n are weakly hyperbolic. Finally, given N 0 andφ as in (a) one has which concludes the proof of (a).
(c) In case the solution φ 0 is hyperbolic, then due to (2.10) and the compactness of the period operators Ξ n θ (see above), Theorem 2.1(b) applies to K n t := DF n t (φ n t ), n ∈ N 0 . It follows that the solutions φ n are hyperbolic as well.
The dynamics of (Δ n ) in the vicinity of hyperbolic solutions φ n is determined by a saddle-point structure consisting of local stable and unstable manifolds resp. fiber bundles [24,p. 256ff,Sect. 4.6] (in the periodic case). These sets allow a dynamical characterization and, given some r 0 > 0, we define the local stable fiber bundle there exists a solution (φ t ) t≤τ of (Δ n ) with φ τ = u τ and φ t − φ n t − −−− → t→−∞ 0 associate to φ n . The following result relates the fiber bundles of the perturbed equations (Δ n ), n ∈ N, to that of the initial problem (Δ 0 ). It requires that F n t n∈N is equidifferentiable in each u ∈ U t , that is there exists a DF n t (u) ∈ L(X ) such that holds uniformly in n ∈ N.
We can now show that the saddle-point structure near hyperbolic periodic solutions to (Δ 0 ) is preserved under perturbation (see Fig. 1).

Theorem 2.3 (Perturbed stable and unstable fiber bundles)
Let θ = lcm {θ 0 , θ 1 } and m ∈ N. Suppose that the θ 0 -periodic difference equations (Δ n ), n ∈ N 0 , fulfill: (ii) DF n t : U t → L(X ), n ∈ N, are uniformly continuous on bounded sets uniformly in n ∈ N, the family DF n t n∈N is equicontinuous for all 0 ≤ t < θ 0 and for every n ∈ N there exists a 0 ≤ t < θ 0 such that DF n t has compact values. If φ 0 is a hyperbolic θ 1 -periodic solution to (Δ 0 ) satisfying (2.7), (2.8) and (P t ) t∈Z is the invariant projector onto the stable vector bundle V 0 + of (V 0 ) (cf. Theorem 2.1), then there exist ρ 1 > 0 and integers N 1 ≥ N 0 so that the following holds for n ≥ N 1 or n = 0, and the θ -periodic hyperbolic solutions φ n ensured by Theorem 2.2: (a) The local stable fiber bundle W n + of (Δ n ) allows the representation as graph of a mapping w n + : Z × X → X with w n + (τ + θ, u) = w n + (τ, u) = w n + (τ, P τ u) ∈ N (P τ ) for all τ ∈ Z and u ∈ X . Moreover, w n + (τ, 0) ≡ 0 on Z, the Lipschitz mappings w n + (τ, ·) are of class C m and the stable fiber bundles of (Δ n ) and (Δ 0 ) are related via The local unstable fiber bundle W n − of (Δ n ) allows the representation as graph of a mapping w n − : Z × X → X with and u ∈ X . Moreover, w n − (τ, 0) ≡ 0 on Z, the Lipschitz mappings w n − (τ, ·) are of class C m and the unstable fiber bundles of (Δ n ) and (Δ 0 ) are related via v)), and have the same finite dimension.
In order to achieve convergence as n → ∞ via (2.11) and (2.12) one needs the derivatives DF n t to tend to DF 0 t on bounded sets and, thanks to Theorem 2.2, continuity of the derivative DF 0 t , 0 ≤ t < θ 0 . A concrete illustration follows in Sect. 3. Remark 2.1 (Alternative representation of W n + and W n − ) With someρ 1 > 0 the local stable and unstable fiber bundles of φ n allow the alternative characterization as graphs over the vector bundles V n + resp. V n − of the variational equations (V n ), rather than over the vector bundles V 0 + resp. V 0 − of (V 0 ) (cf. [24, pp. 256ff, Sect. 4.6]) as in Theorem 2.3. In addition, then the associate mappingsw n + (τ, ·),w n − (τ, ·) possess values in N (P n τ ) resp. in R(P n τ ) for all τ ∈ Z. According to Theorem 2.1(c) the corresponding invariant projectors for (V n ) satisfy lim n→∞ P n t − P 0 t = 0 for all t ∈ Z. Therefore, W n − and W 0 − share their finite dimension. Proof Since the existence of W 0 + , W 0 − and their properties are well-established in the literature [24, pp. 187ff], we focus on their persistence and the convergence estimates (2.11) and (2.12). Let φ n = (φ n t ) t∈Z denote the θ -periodic solutions of (Δ n ) guaranteed by Theorem 2.2 for n ≥ N 0 . The associate equations of perturbed motion are θ -periodic and have the trivial solution. The general solutions ϕ n of (Δ n ) andφ n to (Δ n ) are related byφ n (t; τ, u) = ϕ n (t; τ, u + φ n τ ) − φ n t for all τ ≤ t. (a) For each fixed τ ∈ Z the sequence space + τ := (φ t ) τ ≤t : φ t ∈ X and lim t→∞ φ t = 0 , is complete w.r.t. the sup-norm φ ∞ := sup τ ≤t φ t . Forρ > 0 so small that φ t <ρ implies φ t +φ n t ∈ U t for all t ∈ Z and n ≥ N 0 we introduce the operator for all τ ≤ t. Then u τ = P τ u τ + [I X − P τ ]u τ ∈ X is contained in the stable bundle of (Δ n ) if and only if φ :=φ n (·; τ, u τ ) satisfies (cf. [6, proof of Thm. 3.1]) Our approach to (2.13) using the Lipschitz inverse function Theorem A.2 is based on the representation T n + = A + + G n + with for all τ ≤ t. Note that the derivatives DF 0 t : U t → L(X ) exist by assumption (i).
(I) Claim: is θ -periodic and therefore A + is bounded. In order to show that A + is invertible, given v τ ∈ R(P τ ) and a sequence ψ ∈ + τ , we observe that A + φ = (v τ , ψ) has the unique solution holds for all n ≥ N 1 . Due to the limit relation (2.10) in the proof of Theorem 2.2 there is an we finally obtain for all φ,φ ∈ B ρ (0, + τ ) that (III) In this step we apply the Lipschitz inverse function Theorem A.2 to solve the nonlinear equation (2.13) in the Banach spaces X = + τ , Y = R(P τ ) × + τ , points x 0 := 0, y 0 := (P τ u τ , 0), the Lipschitz constant l := 1−β 4K and σ := and there exists a unique solution φ n + (u τ ) ∈ B ρ (0, + τ ) to (2.13). Then the function w n + parametrizing the stable bundle of (Δ n ) is w n Becauseφ solves (Δ n ) and φ solves (Δ 0 ), this simplifies to and it remains to estimate the right-hand side in this inequality. Since U t is assumed to be convex, we apply the Mean Value Theorem [19,p. 341,Thm. 4.2] and arrive at Here, let us point out that for all integers τ ≤ t one has the relation φ t =φ 0 (t; τ, v τ + w 0 The argument is dual to the proof of (a), but now one works in the sequence space − τ := (φ t ) t≤τ : φ t ∈ X and lim t→−∞ φ t = 0 being complete in the supnorm. One applies Theorem A.2 with X = − τ , Y = N (P τ ) × − τ and x 0 := 0, y 0 := (u τ − P τ u τ , 0), l := 1−β 4K , σ := 1−β 2K to the nonlinear operator

Urysohn integrodifference equations
Let us now illustrate the applicability of our abstract perturbation results from Sect. 2, when the initial problem (Δ 0 ) is an integrodifference equation whose right-hand side is an Urysohn operator over a compact nonempty Ω ⊂ R κ . For the sake of having well-defined and smooth mappings F 0 t , t ∈ Z, in an ambient setting, several assumptions on the kernel functions f t are due: fulfill the following assumptions for all 0 ≤ t < θ 0 and 0 ≤ k ≤ m: exists as continuous function, (ii) for all r > 0 there exists a continuous function h r : Ω → R + such that for all x,x, y ∈ Ω, z ∈ Z t ∩B r (0), (iii) for all r > 0 there exists a function c r : R + × Ω → R + satisfying the limit relation lim δ 0 sup y∈Ω c r (δ, y) = 0, such that |z −z| ≤ δ implies ≤ c r (δ, y) |x −x| α for all x,x, y ∈ Ω,z ∈ Z t ∩B r (0).
Since the compact domain Ω is fixed throughout, we conveniently abbreviate and obtain the open sets U t := u ∈ C α d : u(Ω) ⊂ Z t for all t ∈ Z. For our subsequent analysis it is important to note that Hypothesis 3.1 implies the corresponding assumptions made in [28,Sect. 2]. In detail, one has: Proposition 3.1 (Properties of (I 0 )) Let t ∈ Z. If Hypothesis 3.1 holds, then the Urysohn operator F 0 t = F 0 t+θ 0 : U t → C α d is well-defined, completely continuous and of class C m with compact derivative Combined with the solution identity this shows that entire solutions φ to (I 0 ) inherit the smoothness of the kernel function, i.e. φ t ∈ C α d , t ∈ Z. Yet for kernel functions of convolution type a higher smoothness can be expected (cf. [28,Sect. 2

.3]).
Proof Above all, (I 0 ) and (3.1) show that F 0 t is θ 0 -periodic in t. The results from [28] formulated in an abstract measure-theoretical set-up apply to F 0 t with the κ-dimensional Lebesgue measure μ = λ κ . By [28,Thm. 2.6], F 0 t is well-defined and due to [28, Cor. 2.7(i)] also completely continuous. In [28,Thm. 2.12] it is shown that F 0 t is of class C m and [23, p. 89, Prop. 6.5] implies that DF 0 t (u), u ∈ U t , is compact. Corollary 3.1 Let t ∈ Z and 2 ≤ m. If for every r > 0 there exists a continuous function l r : is Lipschitz on C 0 d -bounded sets, that is, for each r > 0 there exists a L r ≥ 0 such that with the Lipschitz constant L r := max sup ξ ∈Ω Ω l r (ξ, y) dy, Ω h r (y) dy .
(I) We derive that y) dy u −ū α v α holds after passing to the least upper bound over all x ∈ Ω.
Along with IDEs (I 0 ) we now consider their Nyström discretizations. They are based on quadrature (κ = 1) or cubature rules (κ > 1), i.e. a family of mappings determined by a grid Ω n ⊂ Ω of finitely many nodes η ∈ Ω n and weights w η ≥ 0; the dependence of w η on n ∈ N is suppressed here. A rule (Q n ) is called (cf. [16]) convergent, if lim n→∞ Q n u = Ω u(y) dy holds for all u ∈ C 0 d , Thanks to [16,p. 20,Thm. 1.4.17], convergence implies stability. In order to evaluate the right-hand side of (I 0 ) approximately, we replace the integral by a convergent integration rule (Q n ), n ∈ N. The resulting Nyström method (see [4,16] for integral equations) yields the family of difference equations Proposition 3.2 (Properties of (I n )) Let t ∈ Z. If Hypothesis 3.1 holds, then the discrete Urysohn operator F n t = F n t+θ 0 : U t → C α d , n ∈ N, is well-defined, completely continuous and of class C m with compact derivative Moreover, if (Q n ) is stable, then F n t n∈N is equidifferentiable, DF n t are uniformly continuous on bounded sets uniformly in n ∈ N and DF n t n∈N is equicontinuous. Proof The grids Ω n , n ∈ N, are a family of compact and discrete subsets of Ω. If we equip them with the measure μ(Ω n ) := η∈Ω n w η , then due to [28, Ex. 2.2 and Rem. 2.5] the abstract measure-theoretical integral from [28] becomes Ω n f t (x, y, u(y)) dμ(y) = η∈Ω n w η f t (x, η, u(η)) for all x ∈ Ω and leads to the discrete integral operators in (I n ). Given this, well-definedness, complete continuity and smoothness of F n t result from [28] as in the proof of Proposition 3.1. From now on, assume that (Q n ) is stable and choose u ∈ U t .

Hölder continuous kernel functions
We say an integration rule (Q n ) has consistency order α ∈ (0, 1] (cf. [16,p. 21,Def. 1.4.19]), if there exists a c 0 ≥ 0 with Also the midpoint rule Q n M u := b−a n n−1 j=0 u(a + ( j + 1 2 ) b−a n ) is convergent and as in [8, p. 52, Theorem] one derives the quadrature error The trapezoidal rule Q n T u := 1 2 (Q n L R u + Q n R R u) is convergent with the same quadrature error as for the rectangular rules. Finally, let n ∈ N be even. Representing the Simpson rule as convex combination Q n S u : The next two results provide sufficient conditions on the kernel functions f t such that the assumptions (2.7) or (2.8) are satisfied for Nyström discretizations (I n ).
for all x, y,ȳ ∈ Ω and z,z ∈ Z t ∩B r (0). If (Q n ) has consistency order α, then for every r > 0 there exists a c 0 r ≥ 0 such that (3.10)

Differentiable kernel functions
Convergence rates improving the consistency order α ∈ (0, 1] obtained in Theorem 3.1 can be expected for integrands in (I 0 ) being differentiable in y ∈ Ω. Here we follow the convention to consider a function on a not necessarily open set Ω ⊂ R κ as differentiable, if it allows a differentiable extension to an open superset of Ω.
Given p-times continuously differentiable functions u : Ω → R d assume that (Q n ) allows a quadrature or cubature error of the form (see [8]) with constants c p ≥ 0. A smooth framework allows the following improvement of Proposition 3.3: Proposition 3.5 (Higher order convergence of F n t ) Let t ∈ Z, p ∈ N and Ω be convex. Suppose the kernel function f t : Ω 2 × Z t → R d fulfills: If (Q n ) satisfies (3.14), then for every r > 0 there exists ac 0 r ≥ 0 such that and p-times continuously differentiable functions u ∈ U t .
Proof Let t ∈ Z and with u ∈ U t of class C p it is convenient to define t (x, y) := D 1 f t (x, y, u(y)).
The estimate (3.15) for the · 0 -norm is an immediate consequence of the error estimate (3.14) and the higher-order chain rule. Let x,x ∈ Ω and the Mean Value Theorem [19,p. 341,Thm. 4.2] gives

Numerical illustrations
In order to illustrate our theoretical results, we consider an autonomous logistic IDE with separable kernel in 1d. It allows to analyze the behavior of periodic solutions under approximation (largely) explicitly. We demonstrate this by means of Nyström discretizations based on several quadrature rules (taken from [11, pp. 361ff] with the linearly independent functions e 1 , e 2 : being L 2 -orthonormal, i.e. depends on a growth parameter a > 0. The ansatz u t := v t e 1 + w t e 2 with coefficients v t , w t ∈ R leads to the autonomous planar difference equation which fully describes the dynamics of (4.1). It has the invariant set R × {0} and when restricted to this v-axis the behavior is determined by the scalar difference equation If we fix L = π , then (4.3) possesses a pair of 5-periodic solutions emanating from a supercritical fold bifurcation at a ≈ 5.8164. These solutions and their hyperbolicity are illustrated in Fig. 2(the graphics were obtained using a simple continuation scheme with Brent's method [29, p. 256, Sect. 6.2.3] as corrector). In particular, for a = 6 we  of (4.2) resulting in the 5-periodic solution φ 0 t := v * t e 1 to the IDE (4.1). This provides us with a precisely known periodic reference solution for the logistic IDE (4.1) in order to illustrate Theorem 2.2 when applied to its Nyström discretizations (I n ).
In order to determine the 5-periodic solutions φ n of (I n ) we apply a Newton solver to the fixed point problemF n (û) =û havingφ 0 as initial value. For various common summed integration methods we display the development of the They illustrate that the solutions φ n approximate φ 0 preserving the order of the particular quadrature rule. This confirms Theorem 2.2 and specifically the estimate (3.16) of Theorem 3.2. An exception is the fast convergence of the Clenshaw-Curtis method (see Fig. 4) due to the fact that the functions to be approximated are cosine shaped.

Summary and perspectives
This paper studied the local dynamics of abstract periodic difference equations of the form (Δ 0 ) near periodic solutions under discretizations. For sufficiently accurate approximations, in Theorem 2.1 it is established that hyperbolicity of a periodic linear equation persists. Applied to the variational equations (V n ), Theorem 2.2 guarantees that also hyperbolic periodic solutions persist as hyperbolic solutions to the discretizations (Δ n ). The same holds true for their associated stable and unstable manifolds (fiber bundles), as shown in Theorem 2.3.
These abstract perturbation results might apply to various spatial discretizations of evolutionary equations, such as projection methods. We nevertheless illustrate them in terms of integrodifference equations (I 0 ) and their Nyström approximations (I n ), which can be immediately implemented in simulations. Indeed, for these full discretizations the integral is replaced by a convergent quadrature/cubature rule. Its convergence rate regarding the error to the discretized objects (hyperbolic periodic solutions, graphs of their stable and unstable fiber bundles) is preserved throughout. This rate in turn depends on the smoothness of the kernel functions f t . In case they are Hölder continuous, then the convergence rate coincides with the corresponding Hölder exponent (cf. Theorem 3.1), while for higher-order smoothness one obtains polynomial convergence as shown in Theorem 3.2.
Our contribution concentrates on Nyström methods. Nevertheless, rather than using quadrature/cubature formulas, an alternative tool to evaluate the right-hand sides of IDEs is the Fast Fourier Transformation (for short FFT, [22, pp. 106ff, Sect. 8.2]). Although this approach is fairly popular in the literature, it is restricted to Hammerstein IDEs of convolution type u t+1 (x) = Ω k t (x − y)g t (u t (y)) dy, while we deal with general Urysohn equations (I 0 ). In comparison, as pointed out in [22], in order to arrive at a similar accuracy, Nyström methods require less nodes but tend to be slower than FFT methods. Up to the author's knowledge, questions concerning the Numerical Dynamics for FFT discretizations were not tackled so far.
We finally point out that global dynamical features of integrodifference equations were addressed in [9] (for polynomial growth functions) and [10] (for general smooth growth functions). Using rigorous numerics they establish the existence of periodic solutions, connecting orbits, and ultimately chaotic dynamics. The methods involve set oriented numerics and topological tools such as the Conley index, which are fundamentally different from ours. Regarding connecting orbits, see also [20].