This is a preprint of: “Asymptotic dynamics of a difference equation with a parabolic equilibrium”, Bartomeu Coll, Armengol Gasull, Rafel Prohens, Qual. Theory Dyn. Syst., vol. 19, 2020. DOI: [10.1007/s12346-020-00406-0] ASYMPTOTIC DYNAMICS OF A DIFFERENCE EQUATION WITH A PARABOLIC EQUILIBRIUM B. COLL† , A. GASULL‡ AND R. PROHENS∗ Abstract. The aim of this work is the study of the asymptotic dynamical behaviour, of solutions that approach parabolic fixed points in difference equations. In one dimensional difference equa- tions, we present the asymptotic development for positive solu- tions tending to the fixed point. For higher dimensions, through the study of two families of difference equations in the two and three dimensional case, we take a look at the asymptotic dynamic behaviour. To show the existence of solutions we rely on the parametrization method. 1. Introduction and main results In studying difference equations, one interesting problem is to know the asymptotic dynamical behaviour, of the positive solutions that ap- proach equilibrium points. This question naturally arises, for instance, when describing applications in biological or economical systems. See [10, 19, 23, 24, 32, 36], for instance. In particular, how populations evolve to equilibrium states and, more concretely, if they tend to rest points in the same or different asymptotically way is of great interest in mathematical modeling. This problem appears not only in modeling real-life processes but also in Celestial Mechanics ([16, 25, 18, 28, 30]), in complex analytical dynamics ([26]) and in many more others fields where the subject is concerned with the study of invariant submanifolds of fixed parabolic points of maps. Related to the determination of the complete asymptotic expansion of solutions tending to fixed points in difference equations are the next three examples appearing in the literature, in which this issue is un- solved. Key words and phrases. Difference equations, parabolic equilibrium, asymptotic development 2010 Mathematics Subject Classification: 39A10, 41A60, 39A30, 37C75 . The second author is partially supported by Spanish MCYT/FEDER grant num- ber MTM2016-77278-P and Gov. Catalunya grant number 2017SGR 1617; the first and third authors are partially supported by Spanish and European Regional De- velopment Funds (ERDF, FEDER) MICINN MTM2017-83568-P grants. The third author is also partially supported by MICINN MTM2014-54275-P Spanish grant. 1 2 B. COLL, A. GASULL AND R. PROHENS The first one was proposed by Berg and Stević in 2002, see [6] and [33] respectively, concerning the difference equation xn−k xn+1 = , (1) 1 + xn + · · · + xn−k+1 for the case k = 1. In [6] the asymptotic expansion of xn is bounded as 2 2 a 2 2 b + 2 ln n + 3 ln2 n ≤ xn ≤ + 2 ln n + 3 ln2 n, n n n n n n where a < 2 < b. As a consequence, a proof of the existence of a particular solution of this equation such that {xn }n → 0 is also given. The second example concerns the difference equation (1) in the case k = 2. In general, i.e. for each k ∈ N, Stević in 2006, see [37], proved the existence of a positive solution, {xn }n , converging to zero, by assuming that the first five terms in the asymptotic expansion of xn have the following form: a b ln n + c d ln2 n + e ln n xn ∼ + + . n n2 n3 for some a, b, c, d, e ∈ R. Finally, the third example concerns the difference equation xn−k xn = , (2) 1 + xn−1 . . . xn−k+1 when k = 3. This equation was introduced by Berg in 2008 and by Berg-Stević in 2011, see [8, 9] respectively. In these works the existence of a solution such that {xn }n → 0 is proved for all k, by assuming that the asymptotic expansion of xn is given by   1 b 1 2 xn ∼ √ a + + 2 (c ln n + d ln n + e) , n n n where the coefficients are fixed in such a way that xn is proved to be a solution of equation (2). In our work, we answer to the problem of determining the complete asymptotic expansion of xn , for one dimensional difference equations having a parabolic equilibrium point. As a consequence, we solve the problem of giving the asymptotic development of the solution, xn , men- tioned in the above first example. This problem is still open in the other two previous examples. It is worth mentioning some previous works that address the problem of obtaining the first terms of asymptotic developments and motivated our work. See, for instance, [5, 7, 20, 31, 35, 38, 39]. These articles, using various tools, manage to obtain the first terms of the development of some families related to those presented in our work. Among all these references, due to their proximity to some of the techniques that we use in our work, we want to highlight the article of Stević, [31], where the first three terms of the asymptotic development are found. ASYMPTOTIC DYNAMICS 3 The study of the asymptotic development in terms of n, for positive solutions, {xn }n∈N , approaching parabolic fixed points, of difference equations is the general subject of our work. Since the topic is classical, some of the results in this paper appear scattered in this or that form, but we have not managed to find concrete references for all of them. So we have decided to make our paper as self-contained as possible. By assuming that the fixed point is at the origin, this problem can be stated as follows. Let us take the difference equation xn+1 = F (xn ), xn ∈ U, for all n ∈ N ∪ {0}, (3) where F : U → Rk1 × Rk2 , F ∈ C ∞ (U ), k1 , k2 ≥ 0, k1 + k2 ≤ 3, is an analytic function defined on an open subset of the origin, U ⊂ Rk1 ×Rk2 , with a parabolic fixed point at x = (0, 0). Additionally, from now on, we will consider positive solutions of k +k equation (3), i.e. only solutions of the form {xn }n∈N , xn ∈ (R+ ) 1 2 will be considered. At this point, we would recall that the origin of Rk1 × Rk2 is a par- abolic fixed point of F if F (0, 0) = (0, 0) and DF (0, 0) = Id; i.e. that map is tangent to the identity. Also, the set of points whose positive iterates converge to the fixed point is invariant by the map and it is called stable invariant set or stable invariant manifold. To decide whether a parabolic fixed point has associated a stable or unstable manifold is still, in general, an open problem. However, it is worth mentioning the existence and uniqueness results in some cases. See [1, 2, 3, 15, 17, 29], for instance. To develop algorithms for the computation of local approximations of invariant manifolds of parabolic fixed points is also of interest. See [4], for instance. Conditions implying a kind of weak hyperbolicity for the fixed point will be specified later. In this work we will prove that there are families of difference equa- tions (3) such that, on the stable manifold of the origin, they present positive solutions, tending to the fix point, with different asymptotic behaviours. We want to comment that, for system (3), we are not only interested in to get the asymptotic development of positive solutions tending to the fixed point located at the origin but also to prove that these are all the positive sequences that tend to the origin. These results are achieved by using the parameterization method [11, 12, 13, 21, 22], in the simultaneous search for an invariant stable manifold, W s , of the origin as an immersion K : V ⊂ Rk1 → Rk1 × Rk2 , and a map R : V → V , where V is a domain which contains 0 on its boundary. This is a novel approach for this type of problems. We note that K is a parametrization of the manifold W s , and R is a representation of the dynamics of F on W s , satisfying the invariance 4 B. COLL, A. GASULL AND R. PROHENS equation F ◦ K = K ◦ R. (4) We use the parametrization method in the proof of Theorem 1.3. We also use the invariance equation in finding numerical evidences of how W s can be, assuming that this invariant stable manifold of the origin exists. The corresponding results are stated in Proposition 3.1. Concerning the one dimensional case of the difference equation (3), i.e. when k1 + k2 = 1, the asymptotic behaviour study of {xn }n∈N satisfying that limn→∞ xn = 0, can be performed by taking X∞ xn+1 = F (xn ) = λxn + ai xnk+i , (5) i=1 for some k ≥ 1 with a1 6= 0, and 0 ≤ |λ| ≤ 1. As in the complex case, by using conformally conjugate functions to F , it is possible reduce this study to more simpler cases according λ. More specifically, the use of the linearization theorem applies when 0 < |λ| < 1, meanwhile the conjugation theorem can be used when λ = 0. In Section 2 we prove an adapted theorem to the dynamics of the real case for these values of λ. In the particular case |λ| = 1, the asymptotic behaviour is known in some cases. In complex dynamics, Resman, in [27, Prop. 3], considers the difference equation defined by the parabolic diffeomorphism X F (z) = z + ai z k+i , z, ai ∈ C, i≥1 for some k ∈ N. In this setting, the asymptotic development of zn is given by k X k+1 k+1 zn = gi n−i/k + gk+1 n− k ln n + o(n− k ln n), i=1 where the coefficients gi = gi (k, A, a2 , . . . , ai ), i = 1, . . . , k + 1, are 1 complex-valued functions, and A = (−ka1 )− k . When λ = 1, in the one dimensional real case, we present in next theorem a key result giving the asymptotic development of the solution of the difference equation (5). Although we have not found explicit references, their results seem to be of common knowledge. In any case, in next section we include a detailed proof. As we will see, Theorem 1.1, together with the parametrization method, give the clue for understanding the appearance of logarithms in the asymptotic developments in some solutions of difference equa- tions of order bigger than one, like for instance equations (1) or (2). Theorem 1.1 (Asymptotic development of xn ). Let us consider the one dimensional difference equation (5), xn+1 = F (xn ), where F (x) = x + a1 xk+1 + a2 xk+2 + . . . , a1 6= 0, (6) ASYMPTOTIC DYNAMICS 5 being k a positive integer number and x,ai ∈ R, for all i ≥ 1. Let x0 be an initial condition belonging to an attracting domain of the origin, i.e. such that limn→∞ xn = 0. Then, for each positive integer m, the asymptotic development of xn is given by m X k p !   X 1 X 1 i,j j xn = p+i/k gp ln n + O m+1+1/k , (7) p=0 i=1 n j=0 n where coefficients gpi,j = gpi,j (k, x0 , a1 , a2 , . . . ) are real valued functions. Remark 1.2. We note that, since we are just involved into the control how the functions lni n/nj/k i = 0, 1, 2, . . . , j = 1, 2, . . . , emerge in the asymptotic development of xn , we are not interested into fix the coefficients gpi,j in the expression (7). We also note that a related result to Theorem 1.1 can be found, for example, in [27] for the complex case and in [31] for the real one. To illustrate Theorem 1.1, its consequences, and the role of the pa- rameter x0 in the asymptotic developments of its statement we consider two simple examples: (a) The sequence satisfying the bilinear difference equation xn xn+1 = = xn − x2n + x3n + O(x4n ), 1 + xn with initial condition x0 near the origin, can be explicitly obtained by introducing the new sequence yn = 1/xn . We get 1   x0 n 1 c c2 c3 1 xn = = c = − 2 + 3 − 4 +O , 1 + nx0 1+ n n n n n n5 where c = 1/x0 . Hence the sequences tending to zero are always pa- rameterized by one parameter c that in this example can be easily related with x0 . Moreover, notice that for this particular simple case, all coefficients associated to logarithmic terms are zero. (b) It is clear that the difference equation associated to the celebrated logistic map xn+1 = µxn (1 − xn ), has a parabolic fixed point at the origin when µ = 1. For this special case the sequences with initial conditions satisfying 0 < x0 < 1 tend to the origin. From the above theorem we can compute their asymptotic development at 0. For instance, we obtain that 1 1 xn = + 2 (c − ln(n)) n n    2  1 2 1 2 ln (n) + 3 c + c + − (1 + 2c) ln(n) + ln (n) + o . n 2 n3 In this case the dependence between c and x0 is not made explicit by the theorem. 6 B. COLL, A. GASULL AND R. PROHENS For the case λ = −1, as a consequence of previous theorem, in Propo- sition 3.1 we obtain a similar result. At this point, we wonder if the previous asymptotic behaviour can be extrapolated to higher dimensions. To give an insight on what hap- pens in dimension higher than one, we present two results concerning difference equations in two and three dimensions, see Theorem 1.3 and Proposition 1.4, respectively. In the two dimensional case we study a given family based on the difference equation (1) when k = 1. In this case, equation (1) is the first component of the two-dimensional scheme iteration F , where F is defined through the shift function G as   2 x F = G := G ◦ G, G(x, y) = y, . 1+y In next theorem we study the difference equation given by F , and we prove that in a open region of the first quadrant, the stable invariant manifold of the origin is given by the graph of an analytic function K. Even more, we also prove that on K the dynamics of the first com- ponent of F is given by expression (7) and, in this sense, the dynamical asymptotic behaviour on K is unique. Before proving next theorem, let us introduce some notation. Con- sider the two dimensional difference equation (3) when k1 = k2 = 1, defined in a neighbourhood, U , of the origin. In this case, by way of notation, we introduce the projectors π 1 (x, y) = x, and π 2 (x, y) = y and, for each r > 0, we define WVs ={(x, y) ∈ U : π 1 F m (x, y) ∈ V, π 2 F m (x, y) > 0, m ≥ 0, F m (x, y) → 0, as m → ∞}, where V = (0, r), as the stable invariant manifold of the origin, for the map F , restricted to the first quadrant. In next theorem we summarize previous results. Theorem 1.3. Let us consider the two dimensional difference equa- tion (3) given by   x y(1 + y) F (x, y) = , , (8) 1+y 1+x+y defined in a neighbourhood of the origin U . Then, the origin, (0, 0), is a parabolic fixed point for F . Moreover, there exist r > 0 and an open region in the first quadrant on which (1) the invariant manifold of the origin, WVs , is the graph of an analytic function, K, where 1 1 9 5 41 K(x) = x − x2 + x3 − x4 + x5 − x6 + O(x7 ), (9) 2 2 16 8 64 ASYMPTOTIC DYNAMICS 7 (2) the dynamics on WVs is given by the analytic function R, 3 5 69 15 R(t) = t − t2 + t3 − t4 + t5 − t6 + O(t7 ), (10) 2 2 16 2 (3) on WVs , the first terms of the asymptotic development of the solutions corresponding to the difference equation tn+1 = R(tn ) are given by     1 1 1 1 2 1 + 4c 1 2 tn = + −c + ln n + 3 c − ln n + ln n n n2 2 n 4 4 1 1,0  + 4 g3 + g31,1 ln n + g31,2 ln2 n + g31,3 ln3 n n (11) 1 1,0 1,1 1,2 2 1,3 3 1,3 4  + 5 g4 + g4 ln n + g4 ln n + g4 ln n + g4 ln n n  4  ln n +o , n5 where 5 3 1 3 5 3 1 g31,0 = −c3 − , g31,1 = + c + c2 , g31,2 = − − c, g33,3 = , 48 16 2 2 16 4 8 83 5 91 3 3 g41,0 = + c + c4 , g41,1 = − − c − c2 − 2c3 , 864 24 288 8 4 17 7 3 13 1 1 g41,2 = + c + c2 , g41,3 = − − c, g41,4 = , 48 8 2 48 2 16 being c a constant parameter that is fixed by initial condition of the orbit. For three dimensional dynamics we present positive orbits, {xn }n∈N , of equation (3), exhibiting different asymptotic behaviours when ap- proaching the origin. We proceed by using the invariance equation, to obtain the first terms of their asymptotic developments and so their numerical approaches. We remark that due to the facilities to implement the parametriza- tion method, almost any computer algebra system (CAS) can be used to obtain the rational coefficients appearing in the asymptotic expan- sions given in the above theorem and in forthcoming Propositions 1.4 and 1.5. In our work we have used CAS Maple c . More concretely, we provide two families of difference equations (3). For the former one, we give two positive orbits, {xn }n∈N , with different asymptotic speeds developments in n when approaching the origin. For the latter one, we obtain another asymptotic speed development expression on n, for the dynamics of an orbit tending to the origin. Both examples agree with the following scheme. Consider the family of “shift” functions G(x, y, z) = (y, z, g(x, y, z)) , (12) 8 B. COLL, A. GASULL AND R. PROHENS where g is an analytic function defined in a neighborhood of the origin, with g(0, 0, 0) = 0. For each fixed function g, function G defines a difference equation. We note that, even though the origin is a fixed point of this difference equation, it is not a parabolic one. One way to have a parabolic fixed point at the origin is to consider the difference equation (3) given by F = G3 := G ◦ G ◦ G, (13) for some suitable function g. In the three dimensional case we study two families of difference equations based on the difference equations (1) and (2), where g is given either by x g(x, y, z) = (14) 1+y+z or x g(x, y, z) = . (15) 1 + yz We remark that the first component of the three-dimensional scheme iteration defined by g as in (14), agrees with the difference equation (1) when k = 2, proposed in [37], while the family defined by g as in (15) is related to the difference equation (2) proposed in [8]. Next proposition shows the results obtained on the asymptotic be- haviour for both families. In Section 3 we present, besides its proof, some details on their dynamics. Proposition 1.4. Let us consider the three dimensional difference equation (3), defined in a neighbourhood of the origin, where F and G are given by (13) and (12) with g as in (14). Then, there exist numerical evidences on the following facts. (1) there exist three invariant manifolds of dimension two, Ki (t, s), i = 1, 2, 3, whose first analytic development terms are given by:  1 1 1 K1 (t, s) = t, t − t2 + ( s + )t3 + O(t4 ), 2 6 4  3 2 2 17 33 3 4 ( − s)t + (−2s + s − )t + O(t ) , (16) 2 2 4  3 13 21 K2 (t, s) = ( − s)t2 + (−2s2 + s − )t3 + O(t4 ), t, 2 2 4  1 2 1 1 3 4 t − t + ( s + )t + O(t ) , (17) 2 6 4  3 11 15 K3 (t, s) = ( − s)t2 + (−2s2 + s − )t3 + O(t4 ), 2 2 4  1 2 1 1 3 4 (t + t + (− s + )t + O(t ), t , (18) 2 6 4 ASYMPTOTIC DYNAMICS 9 such that the curve Ki (t, 3/2), i = 1, 2, 3, is a non-negative solution of equation (3), the asymptotic behaviour of it agrees with expression (11). (2) there exists a positive solution of equation (3), defined in a neighbourhood of the origin, analytically approached by  2 10 58 K(t) = t, t − t2 + t3 − t4 + O(t5 ), 3 9 27  4 2 28 3 224 4 5 t− t + t − t + O(t ) , (19) 3 9 27 (3) The dynamics on K is approached by the analytic function R, 182 4 214 5 20762 6 R(t) = t − 2t2 + 6t3 − t + t − t + O(t7 ), (20) 9 3 81 (4) the first terms of the asymptotic development of R are given by   1 1 1 1 2 2 xn = + (−2c + ln n) + 3 + 4c − (1 + 4c) ln n + ln n 2n 4n2 8n 9 1  + 4 g31,0 + g31,1 ln n + g31,2 ln2 n + ln3 n 16n  4  1 1,0 1,1 1,2 2 1,3 3 4  ln n + 5 g4 + g4 ln n + g4 ln n + g4 ln n + ln n + o , 32n n5 where 41 4 11 5 g31,0 = − − c − 8c3 , g31,1 = + 4c + 12c2 , g31,2 = − − 6c, 36 9 6 2 607 41 4 353 70 g41,0 = + c + c2 + 16c4 , g41,1 = − − c − 12c2 − 32c3 , 243 9 3 54 9 19 13 g41,2 = + 14c + 24c2 , g41,3 = − − 8c, 3 3 (21) being c a constant parameter that is fixed by initial condition of the orbit. Proposition 1.5. Let us consider the three dimensional difference equation (3), defined in a neighbourhood of the origin, where where F and G are given by (13) and (12) with g as in (15). Then, there exist numerical evidences on the following facts. (1) there exist a positive invariant curve, K(t), whose first analytic development terms are given by:   1 3 1 5 7 2 3 5 7 K(t) = t, t − t + t + O(t ), t − t + t + O(t ) , (22) 3 3 3 (2) The dynamics on K is approached by the analytic function R 41 7 32 9 2014 11 R(t) = t − t3 + 2t5 − t + t − t + O(t13 ). (23) 9 3 81 10 B. COLL, A. GASULL AND R. PROHENS (3) the first terms of the asymptotic development of R are given by      1 1 c 1 1 7 xn = √ 1+ − + ln n − + ln n 2n n 2 8 32n2 9   1 55 31 2 + + ln n − ln n 256n3 18 9   3  1 1,0 1,1 1,2 2 1,3 3  ln n + 4 g4 + g4 ln n + g4 ln n + g4 ln n + o , n n9/2 where 847 37 17 1 g41,0 = − , g41,1 = − , g41,2 = , g41,3 = − , 62208 6912 4608 1536 (24) being c a constant parameter that is fixed by initial condition of the orbit. This paper is organized as follows. In Section 2, in the one dimen- sional real case, when λ = 1, we present the proof of Theorem 1.1. In this section, we also study the case λ = −1, where some results on the asymptotic dynamic behaviour of solutions are stated. See Proposi- tion 2.3. In Section 3, we give the proofs of the Theorem 1.3, Proposi- tion 1.4 and, Proposition 1.5. 2. Asymptotic dynamical properties in one dimension In this section we present some properties about asymptotic dynamic behaviour of the solutions, near the fix point at the origin, of the one dimensional difference equation (5). First, in Lemma 2.1, we introduce a result recalling the asymptotic behaviour, except in the rationally neutral cases, i.e. except when |λ| = 1. In Theorem 1.1, the rationally neutral case λ = 1 is studied. Ad- ditionally, in Lema 2.3 we present some results on the asymptotic dy- namic behaviour for the λ = −1 case. For this value of λ, it is worth- while to mention that the asymptotic behaviour of the solutions tending to the origin of the difference equation (5) are studied by considering the iterates of the difference equation xn+1 = F 2 (xn ), taking either x0 or F (x0 ) as the initial condition, as we detail later. Lemma 2.1. Consider the difference equation (5). For each x0 close to zero, the solution xn with initial condition x0 satisfies that: (1) if 0 < |λ| < 1, then xn = ϕ(y0 λn ) = α1 y0 λn + α2 y0 λ2n + O(λ3n ), where α1 6= 0 and ϕ is an analytic diffeomorphism at 0 such that ϕ(0) = 0, and y0 = ϕ−1 (x0 ). ASYMPTOTIC DYNAMICS 11 (2) if λ = 0 and F (x) = ap xp + O(xp+1 ) with ap = 6 0 and p ≥ 2, then n n n n xn = ϕ(y0p ) = α1 y0p + α2 y02p + O y03p , where α1 6= 0 and ϕ is an analytic diffeomorphism at 0 such that ϕ(0) = 0, and y0 = ϕ−1 (x0 ). Proof. In the case 0 < |λ| < 1, i.e. when the origin of equation (5) is an attracting hyperbolic fixed point, following Koenigs Theorem (see [14, Th. 2.1]), there is a unique, up to a scale factor, analytic one- to-one real map y = φ(x), ϕ = φ−1 , of a neighbourhood of the origin onto itself, which conjugates F (x) to the linear function G(y) = λy; what finishes the proof in this case because the sequences generated by yn+1 = λyn are yn = y0 λn . In the λ = 0 case, the origin of equation (5) is a superattracting fixed point and, following Boettcher result, see [14, Th. 4.1] for instance, we get that there is a complex analytic map y = φ(x), with ϕ = φ−1 , of a neighbourhood of the origin onto itself which conjugates F (x) to y p . So, the proof follows because the sequences generated by yn+1 = ynp are n yn = y0p .  Next lemma is an easy technical result which is useful to prove The- orem 1.1. Its proof is included for the sake of completeness. Lemma 2.2. Given n0 ∈ N, let us consider the recurrence relation un+1 − un = f (n) + o (f (n)) , n ≥ n0 , n ∈ N, (25) where f is a real and continuous non-negative function, defined on R. Suppose, additionally, that un0 > 0 and that f is a monotonous function on [n0 , +∞). Consider a function F , F ∈ C 1 (R), such that F 0 (x) = f (x), x ≥ n0 , and suppose that one of the following hypotheses holds: (1) limx→+∞ F (x) = ∞ and limn→+∞ f (n)/f (n + 1) = 1, (2) limx→+∞ F (x) ∈ R Then there exists a function G(x), x ≥ n0 , such that G0 (x) = f (x) and satisfying un = G(n) + o (G(n)) . (26) Proof. To prove equality (26), let us prove the equivalent condition un lim = 1. (27) n→+∞ G(n) By assuming hypothesis (1) and taking G(x) = F (x), we get that limx→+∞ G(x) = ∞. Hence, by applying the Stolz-Cèsaro criteria, to prove (27), we must compute un+1 − un lim . n→+∞ F (n + 1) − F (n) 12 B. COLL, A. GASULL AND R. PROHENS Straightforward calculations show that un+1 − un f (n) + o (f (n)) f (n) = , (28) F (n + 1) − F (n) f (n) f (ξn ) where ξn ∈ (n, n + 1). By assuming that f is a monotonous decreasing function on [n0 , +∞), we have that f (n) f (n) 1≤ ≤ , f (ξn ) f (n + 1) for all n ∈ N. By considering the limiting case of former inequalities sequence, when n tends to infinity, and by using the second hypothesis of case (1), we conclude that the limit of expression (28) exists and, hence, that desired conclusion follows. In the case when f is a monot- onous increasing function on [n0 , +∞), by using analogous arguments, the same conclusion holds. Assume now hypothesis (2). Since limx→+∞ F (x) = a ∈ R then, we have that f is a monotonous decreasing to zero function, by one hand, and that it is not restrictive to assume that a > 0, on the other. Since, for all n ∈ N and for all p ∈ N, p ≥ n + 1, we have that Xp Z p p−1 X f (k) ≤ f (s) ds ≤ f (k), (29) k=n+1 n k=n P then, from the first inequality, then we also have that ∞ k=n+1 f (k) < ∞. Concerning the sequence {um }m∈N , since p p p X X X up+1 − un+1 = f (k) + o (f (k)) = f (k) + o (f (n + 1)) , k=n+1 k=n+1 k=n+1 (30) we get that {um }m∈N is a Cauchy sequence and, then there exists l ∈ R, such that limp→∞ up = l. By considering the limiting case of the inequalities sequence given by (29), when p tends to infinity, using equality (30), we get l − un+1 + o (f (n + 1)) ≤ a − F (n) ≤ l − un + o (f (n)) . From the first and second inequalities we get l − a + F (n − 1) + o (f (n)) ≤ un , un ≤ l − a + F (n) + o (f (n)) , respectively. As a consequence, G(n − 1) + o (f (n)) un G(n) + o (f (n)) ≤ ≤ , (31) G(n) G(n) G(n) where G(n) = F (n) + l − a. From the recurrence relation (25) and using the hypothesis, we have that l > 0. Using this fact, it can be proved that G(n − 1) lim = 1. n→+∞ G(n) ASYMPTOTIC DYNAMICS 13 As a consequence, by considering the limiting case when n tends to infinity in the expression (31), equality (27) holds.  By using the previous lemma, we prove our main result on the as- ymptotic development of the solution xn of the difference equation (5), in the rationally neutral case λ = 1. Theorem 1.1. Let us consider {xn }n∈N , a solution of xn+1 = F (xn ) = xn + a1 xk+1 n + a2 xk+2 n + ..., (32) where a1 6= 0, and satisfying that limn→∞ xn = 0. Let us consider the change of variables 1 1 x = Aω − k , where A = (−ka1 )− k . (33) We remark that, since limn→∞ xn = 0, the case k even and a1 > 0 can not be taken into account. Hence, the change of variables is well defined. Even more, the convergence to zero of xn means that we choose the initial condition x0 ≈ 0 such that either, k is odd and x0 a1 < 0, or k is even and a1 < 0. As a consequence, we remark that the convergence to zero only depends on a1 , k and x0 ; that is, X ai xin → 0, (34) i≥1 as n goes to infinity, for all ai , i ≥ 2. Using previous change of variables applied to the difference equa- tion (32), we obtain the following recurrence relation for ωn X∞ ci ωn+1 = ωn + 1 + i/k , (35) i=1 ω n for some real values ci , depending on the initial coefficients aj . Recur- rence (35) can be re-written as ωn+1 − ωn = 1 + o(1), because, from requirement (34), we have X ∞ X∞ ci ci i = x → 0, i/k i=1 ωn i=1 Ai n as n goes to infinity. By applying Lemma 2.2.1 to the recurrence on ωn , we get that ωn = n + o(n). fn = n. In this way, we introduce the recurrence Let us define ω pn = ωn − ω fn . (36) In accordance with previous definition, we have ∞ ∞ q X ci X ci X di − kq pn+1 − pn = i/k = i/k = i/k + o(n ), i=1 ω n i=1 (n + o(n)) i=1 n 14 B. COLL, A. GASULL AND R. PROHENS for some di real value, and where q is an arbitrary, but fixed, natural number q  k. By applying Lemma 2.2.2 to the previous recurrence on pn , we get that k−1 k(q+1) X ei X ei pn = + ek ln(n) + + o(n−q ). i=1 n(i−k)/k i=k+1 n(i−k)/k From equality (36), we obtain k−1 k(q+1) X ei X ei fn + pn = n + ωn = ω + ek ln(n) + + o(n−q ). i=1 n(i−k)/k i=k+1 n (i−k)/k If we undo the change of variables (33) on previous expression of ωn , we obtain  1/k A 1 xn = 1/k , (37) n 1+x where k(q+1) X ei 1 x= i/k + ek ln n + o(n−q−1 ). i=1 n n i6=k Now, we are going to prove that if we use Taylor development for expression (37), on the variable x, in a neighbourhood of the origin, up to order p on x, then xn writes as p k m ! X X 1 X 1 1,p+1 j p+1 xn,p = m+i/k g i,j m ln n + p+1+1/k g p+1 ln n + Rp (n), m=0 i=1 n j=0 n (38) for some coefficients g i,j =p g i,j p (k, x0 , a1 , a2 , . . . ) which are real valued functions; and where Rp (n) is a function including the Taylor remainder term. Let us prove formula (38) by using mathematical induction on p. For p = 0, performing a Taylor expansion of first order on x in expression (37), we get   A 1 xn = 1/k 1 − x + R0 (n) , n k where R0 (n) is the Taylor remainder term. This expression agrees with the one given by formula (38). We observe that, in this case,   k(q+1) A X ei R0 (n) = − (i+1)/k + o(n−q−1−1/k ) + R0 (n) . k i=k+1 n Hence, we proved the base case of the induction process. Let us prove the induction step, i.e. that if formula (38) holds for p, then it holds for p + 1. ASYMPTOTIC DYNAMICS 15 If we define the functions γi,j (n) = lni n/nj/k , i = 0, 1, 2, . . . j = 1, 2, . . . , we note that, x and each xn,p is a linear combination of a subset of them. For xn,p+1 , let us check that we obtain the functions γi,j that appear in expression (38). To get xn,p+1 we need a Taylor expansion, of order p + 2, of expres- sion (37), including the xp+2 term. Hence, concerning the γi,j functions, those appearing in xn,p+1 are the ones we had up to order xp+1 plus the new ones coming from the xp+2 term. So, we need to check that expression (38) for xn,p+1 includes these new functions, that are {γi,k(p+1)+j (n), i = 0, 1, . . . , p + 1 j = 1, 2, . . . , k} ∪ {γp+2,k(p+2)+1 (n)}. (39) From the induction hypothesis, we can obtain xn,p+1 by multiplying xn,p times x. We observe that the functions γi,j (n) involved in the expression of xn,p , i.e. in the terms of the Taylor development, up to order xp+1 , are {γi,j (n), i = 0, 1, . . . , p, j = 1, 2, . . . , k(p + 1)} ∪ {γp+1,k(p+1)+1 (n)}, (40) at least, and the ones involved in the expression of x are {γ0,j (n), j = 1, 2, . . . , k(q + 1) j 6= k} ∪ {γ1,k (n)}, (41) at least. We note that the functions γi,j (n) given in expression (39) can be obtained from the product of the ones given in expression (40) by the ones given in expression (41), what finishes the proof.  In the next proposition, which is a sequel of Theorem 1.1, we present some results on the asymptotic development of the solution xn of the difference equation (5), in the rationally neutral case λ = −1. We note that, we are only interested in those initial conditions be- longing to the stable manifold of the fixed point at the origin, i.e. we only consider x0 such that limn→∞ F n (x0 ) = 0. As a consequence, in next result, the case k even and a1 < 0 is excluded. Proposition 2.3. Consider the difference equation (5) with λ = −1, that is X∞ xn+1 = F (xn ) = −xn + ai xnk+i . i=1 For each initial condition, x0 , belonging to an attracting domain of the origin, the corresponding solution xn satisfies that: (1) if k is even, then ∞ X k p ! X 1 X xn = p+i/k gpi,j lnj n , (42) p=0 i=1 n j=0 16 B. COLL, A. GASULL AND R. PROHENS (2) if k is odd, then ∞ X k p ! X 1 X j xn = gpi,j ln n , (43) p=0 i=1 np+i/k j=0 for some integer number k > k, where gpi,j = gpi,j (k, x0 , a1 , a2 , . . . ) are real valued functions. Proof. We get the asymptotic behaviour of an orbit, {xn }n∈N , from the study of the two sub-orbits {x2n }n∈N and {x2n+1 }n∈N , i.e. from the study of the behavior of the iterations through F 2 , starting with initial conditions x0 and F (x0 ), respectively. We note that F 2 (x) = x + a1 ((−1)k+1 − 1)xk+1 + a2 ((−1)k+2 − 1)xk+2 + o(xk+2 ). In the case k even, since a1 6= 0, we apply Theorem 1.1 to the iterations of F 2 on the corresponding initial conditions. Hence, the asymptotic behaviours of {x2n }n∈N and {x2n+1 }n∈N are, both, given by expression (42). In the case k odd, we observe that the coefficient of xk+1 in F 2 vanishes. Hence, either there exist a positive integer number j > 1 and aj 6= 0, such that F 2 (x) = x + aj xk+j + o(xk+j ), or F 2 (x) = x. Since x0 belongs to an attracting domain of the origin, the latter case can not be. Whence, we can also apply Theorem 1.1 to the iterations of F 2 . Hence, the asymptotic behaviours of {x2n }n∈N and {x2n+1 }n∈N are, both, given by expression (43). We remark that, since j > 1, the non-linear terms in F 2 are of order o(xk+1 ) and, hence, expression (43) is taken with a certain k > k.  3. Some results in higher dimensions In this section, to give an insight to the asymptotic behaviour in higher dimensions, we prove Theorem 1.3, Proposition 1.4 and, Propo- sition 1.5, concerned to two particular families of difference equations, in two and three dimensions, respectively. Next technical lemma give us some properties of the difference equa- tion (3) by taking F as in (8). Lemma 3.1. For k = 2, let us consider the two dimensional difference equation (3), defined in a neighbourhood of the origin U , given by (8), i.e. with   x y(1 + y) F (x, y) = , . 1+y 1+x+y Then, the following properties follow. ASYMPTOTIC DYNAMICS 17 (1) The origin, (0, 0), is a parabolic fixed point. (2) Points on the coordinate axes are fixed points. (3) V (x, y) = x2 + y 2 is a Lyapunov function. Proof of Theorem 1.3. Let us consider the two dimensional difference equation (3), defined in a neighbourhood of the origin U , given by (8). Let us prove statement (1) of this theorem. First, we remark that from Lemma 3.1, the origin, (0, 0), is a parabolic fixed point. Now, let us prove that the hypotheses of Theorem 4.1 of [1] are fulfilled and, as a consequence, the invariant manifold of the origin, WVs , is the graph of an analytic function. Following the notation in [1], for each r > 0, we write V (r) = (0, r) ⊂ R+ . By way of notation, for each ρ > 0, we define V 1 (ρ) = {ρ}. In a neighbourhood of the origin, we observe that F is an analytic function given by  F (x, y) = x−xy +xy 2 +o(k(x, y)k3 ), y −xy +xy 2 +x2 y +o(k(x, y)k3 ) . We apply the above mentioned theorem to the map F , after performing the analytic change of variables ϕ, given by X∞ ∞ X i j x1 = ai,j x y , y1 = bi,j xi y j . i,j=1 i,j=1 where a0,1 + a1,0 > 0, b0,1 = −b1,0 = 2b0,2 + 2b1,1 + 2b2,0 , (2b0,2 + b1,1 )b1,0 > 0. If we take ρ > 0 such that (2b0,2 + b1,1 )ρ2 < 2, (a0,1 + a1,0 )2 b1,0 then the hypotheses, H1–H4 of [1, Theorem 4.1], are fulfilled for F (x1 , y1 ). Consequently, there exists r > 0 such that WVs is the graph of an ana- lytic function in V . Undoing previous change of variables, we have that the first quadrant of the (x1 , y1 )-plane corresponds to an open region of the (x, y)-plane whose boundary is given by two disjoint curves. One of these curves is concave and tangent to the y = x straight line at the origin. The other one, reaches the origin with slope −a10 /a01 > 1. Locally, both curves bound an open region on which the invariant manifold of the origin, Wϕs−1 (V ) , is the graph of an analytic function, h. Since h is an analytic function, an approximation of its Taylor ex- pansion is obtained as in expression (9). To prove statement (2), i.e. to obtain a representation, R, of the dynamics on the curve K, we use the parametrization method for one- dimensional invariant manifolds of higher dimensional parabolic fixed points introduced in [2]. We check that the hypotheses of Theorem 2.1 18 B. COLL, A. GASULL AND R. PROHENS of [2] are satisfied and, hence, there exist a C ∞ map K which is a parametrization of the one-dimensional invariant manifold of map F , and a polynomial R(t) such that F ◦ K = K ◦ R. (44) As in statement (1), we apply the above mentioned theorem to the map  F (x1 , y1 ) = F 1 (x1 , y1 ), F 2 (x1 , y1 ) . In this setting, we have F (0, 0) = 0, DF (0, 0) = Id, and ∂ 2F 1 ∂ 3F 2 D2 F 2 (0, 0) = 0, (0, 0) < 0, (0, 0) = 0, ∂x21 ∂x21 as we would prove. The computation of R is done by matching powers of t in (44). Statement (3) of this theorem follows from Theorem 1.1, applied to the difference equation tn+1 = R(tn ), where R is given by expression (10).  To give an insight on the three dimensional asymptotic behaviour, next we present two examples of difference equations (3). By assuming analytical behaviour of their solutions, we numerically approach them. These examples exhibit different dynamical asymptotic development in terms of n, for positive solutions, {xn }n∈N , tending to the origin. Next Lemmas 3.2 and 3.3 are technical results useful to prove Propo- sition 1.4 and Proposition 1.5, respectively. Lemma 3.2. Let us consider the difference equation (3) where F and G are given by (13) and (12) with g as in (14). Then the following prop- erties follow. (1) The origin, (0, 0, 0), is a parabolic fixed point and V (x, y, z) = x2 + y 2 + z 2 is a Lyapunov function. (2) The planes x = 0, y = 0 and z = 0 are invariants by the iteration of the difference equation. Furthermore, all points of the coordinate axes are equilibrium points. (3) On each coordinate plane, the dynamical behaviour of F coin- cides with the one given in Theorem 1.3. Lemma 3.3. Let us consider the difference equation (3) where F and G are given by (13) and (12) with g as in (15). Then the following prop- erties follow. (1) The origin, (0, 0, 0), is a parabolic fixed point and V (x, y, z) = x2 + y 2 + z 2 is a Lyapunov function. (2) All points of the coordinate planes: x = 0, y = 0 and z = 0, are fixed points. ASYMPTOTIC DYNAMICS 19 Proof of Proposition 1.4. To prove statements (1) and (2) of this propo- sition, let us take g as given in expression (14), and let us apply the parameterization method for looking an invariant stable manifold of the origin. By imposing the invariant equation (4), we get the parameterization of the first analytic development terms, of invariant stable curves, in terms of an extra parameter, s, such that s ≤ 3/2. Let’s call Ki (t, s), i = 1, 2, 3, to such non-negative curves as they are established in expres- sions (16), (17) and (18) in the statement (1). We point out that, since this method provides us a representation of the dynamics of F on the in- variant stable curves, the asymptotic behaviour of Ki (t, 3/2), i = 1, 2, 3, agrees with expression (11). Furthermore, from Lemma 3.2.(2) and (3), by using Theorem 1.3, we get on each coordinate plane the existence of an analytic solution given by expression (11). Additionally, the parameterization method also provides us the first analytic development terms of a positive invariant stable curve, K(t), on the invariant stable manifold of the origin, as it is given in expres- sion (19) in statement (2). Furthermore, and also as a consequence of the application of the previous parametrization method, the first ana- lytic terms of the dynamics on K, R, are approached by expression (20) in the statement (3). Finally, by applying Theorem 1.1 to the first analytic terms of R, we get expression (21) in the statement (4).  Proof of Proposition 1.5. To prove statements (1) and (2) of this propo- sition, we proceed as in the proof of Proposition 1.4. More concretely, let us take g as it is given in expression (15) and, looking for an invari- ant stable manifold of the origin, we use the parameterization method. As a consequence, we get expression (22) giving the first analytical terms of the invariant curve K(t) and, additionally, we also get the approach of the dynamics on K, given by R(t), as in expression (23). Then, by applying Theorem 1.1 to the first analytic terms of R, we get expression (24).  In Figure 1, we locally depict the invariant surface (three leaves) on which the solutions of the difference equation (3) where F and G are given by (13) and (12) with g as in (14) present two different be- haviours, according to their initial value. More concretely, solutions tending to the origin and solutions going to the fixed points on the axes. On the intersection of the three leafs given by K(t), given by expression (19), there is a solution, {xn }n , having asymptotic expan- sion given by expression (21); that is, this solution goes to the origin as 1/(2n). On each coordinate plane, as it proved in Lemma 3.1, there 20 B. COLL, A. GASULL AND R. PROHENS is a solution {xn }n , whose asymptotic expansion is given by expres- sion (11), that is a solution tending to the origin as 1/n. Finally, solutions in between the two aforementioned ones, and on the invari- ant surface, belong to the stable invariant manifold of the origin, WVs . Figure 1. According to Proposition 1.4, local depict of the stable invariant manifolds, Ki (t, s) and K(t) (see developments (16), (17), (18) and (19), respectively), of the origin, of the difference equation (3), where F and G are given by (13) and (12) with g as in (14). Different perspectives. Two different speeds, 1/(2n) and 1/n, are observed in the asymptotic developments. In Figure 2, we locally depict the stable invariant manifold of the origin, WVs , of the difference equation (3), where F and G are given by (13) and (12) with g as in (15). The manifold WVs is depicted as the invariant curve K(t) given by expression √ (22). On this curve the solutions {xn }n tend to the origin as 1/ 2n. Coordinate planes, as it proved in Lemma 3.3, are fulfilled with equilibrium points. Acknowledgements The second author is partially supported by Spanish MCYT/FEDER grant number MTM2016-77278-P and Gov. Catalunya grant number 2017SGR 1617; the first and third authors are partially supported by Spanish and European Regional Development Funds (ERDF, FEDER) MICINN MTM2017-83568-P and MICINN MTM2014-54275-P grants. The authors want to thank prof. Rafael Ortega for his comments and suggestions during the realization of this work. We also would like to thank the anonymous referees for their careful reading of the paper and their constructive criticism. ASYMPTOTIC DYNAMICS 21 Figure 2. According to Proposition 1.5, local depict of the stable invariant manifold curve, K(t), of the origin, as it is given in expression (22). This is the case of the difference equation (3), where F and G√ are given by (13) and (12) with g as in (15). Speed 1/ 2n is observed in the asymptotic development. References [1] I. Baldomà and E. Fontich. Stable manifolds associated to fixed points with linear part equal to identity. J. Differential Equations, 197(1), (2004),45–72. [2] I. Baldomà, E. Fontich, R. de la Llave, and P. Martı́n. The parameterization method for one-dimensional invariant manifolds of higher dimensional para- bolic fixed points. Discrete Contin. Dyn. Syst., 17(4), (2007), 835–865. [3] I. Baldomà, E. Fontich and P. Martı́n. Invariant manifolds of parabolic fixed points (I). Existence and depence on parameters. arXiv:1603.02533v1 [math.DS]. 2016. [4] I. Baldomà, E. Fontich and P. Martı́n. Invariant manifolds of para- bolic fixed points (II). Approximations by sums of homogeneous functions. arXiv:1603.02535v1 [math.DS]. 2016. [5] L. Berg. “Asymptotische Darstellungen und Entwicklungen, Hoch- schulbücher für Mathematik, Band 66, VEB Deutscher Verlag der Wissenschaften, Berlin, 1968. [6] L. Berg. On the asymptotics of nonlinear difference equations. Zeitschrift für Analysis und ihre Anwendungen. Journal for Analysis and its Applications, 21, no.4, (2002), 1061–1074. [7] L. Berg. Inclusion theorems for non-linear difference equations with applica- tions. J. Difference Equat. Appl. 10, no.4, (2004), 399–408. [8] L. Berg. On the asymptotics of the difference equation xn−3 = xn (1 + xn−1 xn−2 ). J. Difference Equ. Appl. 14, no.1, (2008), 105–108. [9] L. Berg and S. Stević. On the asymptotics of the difference equation yn (1 + yn−1 . . . yn−k+1 ) = yn−k . J. Difference Equ. Appl. 17, no.4, (2011), 577–586. [10] R. J. Beverton and S. J. Holt. On the Dynamics of Exploited Fish Populations. vol. 19, Fish. Invest., London, 1957. [11] X. Cabré, E. Fontich, and R. de la Llave. The parameterization method for in- variant manifolds. I. Manifolds associated to non-resonant subspaces. Indiana Univ. Math. J., 52-2, (2003), 283–328. 22 B. COLL, A. GASULL AND R. PROHENS [12] X. Cabré, E. Fontich, and R. de la Llave. The parameterization method for invariant manifolds. II. Regularity with respect to parameters. Indiana Univ. Math. J., 52-2, (2003), 329–360. [13] X. Cabré, E. Fontich, and R. de la Llave. The parameterization method for invariant manifolds. III. Overview and applications. J. Differential Equations, 218-2, (2005), 444–515. [14] L. Carleson and T. Gamelin. “Complex Dynamics” Universitext: Tracts in Mathematics, 1993, Edit. Springer-Verlag New York, Inc. ISBN-13: 978-0- 387-97942-7. [15] J. Casasayas, E. Fontich and A. Nunes. Invariant manifolds for a class of parabolic points. Nonlinearity, 5, (1992), 1193–1210. [16] R. W. Easton. Parabolic orbits in the planar three-body problem. J. Differential Equations, 52, (1984), 116–134. [17] E. Fontich. Stable curves asymptotic to a degenerate fixed point. Nonlinear Anal. 35, (1999), 711–733. [18] R. McGehee. A stable manifold theorem for degenerate fixed points with appli- cations to celestial mechanics. J. Differential Equations, 14, (1973), 70–88. [19] E. A. Grove, C. M. Kent, G. Ladas, S. Valicenti and R. Levins. Global stabil- ity in some population models. Communications in Difference Equations. Pro- ceedings of the 4th International Conference on Difference Equations (Poznan, 1998), Gordon and Breach, Amsterdam, 2000, pp. 149–176. [20] L. Gutnik and S. Stević. On the Behaviour of the Solutions of a Second-Order Difference Equation. Discrete Dynamics in Nature and Society, (2007), ID 27562, 14 pp. [21] A. Haro, and R. de la Llave. A parameterization method for the computation of invariant tori and their whiskers in quasi-periodic maps: rigorous results. J. Differential Equations, 228-2, (2006), 530–579. [22] A. Haro, and R. de la Llave. A parameterization method for the computation of invariant tori and their whiskers in quasi-periodic maps: explorations and mechanisms for the breakdown of hyperbolicity. SIAM J. Appl. Dyn. Syst., 6-1, (2007), 142–207. [23] H. F. Huo and W. T. Li. Permanence and global stability of positive solutions of a nonautonomous discrete ratio-dependent predator-prey model. Discrete Dy- namics in Nature and Society, 2005 (2005), no. 2, 135–144. [24] Y. Kuang and J. M. Cushing. Global stability in a nonlinear difference-delay equation model of flour beetle population growth. Journal of Difference Equa- tions and Applications, 2 (1996), no.1, 31–37. [25] R. Martı́nez and C. Pinyol. Parabolic orbits in the elliptic restricted three body problem. J. Differential Equations, 111 (1994), 299–339. [26] J. Milnor. “Dynamics in One Complex Variable” Institute for Mathematical Sciences, SUNY, Stony Brook, NY, 1991. [27] M. Resman, ε-Neighborhoods of orbits and formal classification of parabolic diffeomorphisms. Discrete Cont. Dyn-A, 33-8, (2013), 3767–3790. [28] C. Robinson. Homoclinic orbits and oscillation for the planar three-body prob- lem. J. Differential Equations, 52, (1984), 356–377. [29] C. Simó. Stability of degenerate fixed points of analytic area preserving map- pings. “Bifurcation, Ergodic Theory and Applications (Dijon, 1981)” pp. 184– 194, Astérisque, vol. 98, Soc. Math. France, Paris, 1982. [30] D. L. Slotnick, Asymptotic behavior of solutions of canonical systems near a closed, unstable orbit. “Contributions to The Theory of Non-linear Oscilla- tions, Vol. IV” pp. 85–110, Annals of Mathematics Studies, no. 41, Princeton University Press, Princeton, N.J., 1958. ASYMPTOTIC DYNAMICS 23 [31] S. Stević. Asymptotic behaviour of a sequence defined by iteration. Mat. Vesnik, (3-4), 48 (1996), 99–105. [32] S. Stević. Behavior of the positive solutions of the generalized Beddington-Holt equation. Panamerican Mathematical Journal, 10 (2000), no. 4, 77–85. [33] S. Stević. On the recursive sequence xn+1 = xn−1 /g(xn ). Taiwanese J. Math. 6-3, (2002), 405–414. [34] S. Stević. Asymptotic behavior of a sequence defined by iteration with applica- tions. Colloquium Mathematicum, 93 (2002), no. 2, 267–276. [35] S. Stević. Asymptotic behavior of a nonlinear difference equation. Indian Jour- nal of Pure and Applied Mathematics, 34 (2003), no. 12, 1681–1687. [36] S. Stević. Asymptotic behaviour of a class of nonlinear difference equations. Discrete Dynamics in Nature and Society, 2006, ID 47156, 1–10. [37] S. Stević. On positive solutions of a (k+1)th order difference equation. Applied Mathematics Letters, 19-5, (2006), 427–431. [38] S. Stević. On monotone solutions of some classes of difference equations. Dis- crete Dynamics in Nature and Society, 2006, ID 53890, 1–9. [39] S. Stević. On a Discrete Epidemic Model. Discrete Dynamics in Nature and Society, 2007, ID 87519, 10 pp. †∗ Dept. de Matemàtiques i Informàtica, IAC3 Institute of Applied Computing & Community Code, Universitat de les Illes Balears, 07122 Palma de Mallorca, Illes Balears. Spain Email address: † tomeu.coll@uib.cat, ∗ rafel.prohens@uib.cat ‡ Dept. de Matemàtiques, - . Universitat Autònoma de Barcelona, Edifici C, 08193 Bellaterra, Barcelona. Spain Email address: gasull@mat.uab.cat