Explicit Solutions and Bifurcations for a System of Rational Difference Equations

In this paper, we consider the explicit solution of the following system of nonlinear rational difference equations: xn+1 = xn−1/ (xn−1 + r), yn+1 = xn−1yn/ (xn−1yn + r), with initial conditions x−1, x0 and y0, which are arbitrary positive real numbers. By doing this, we encounter the hypergeometric function. We also investigate global dynamics of this system. The global dynamics of this system consists of two kind of bifurcations.


Introduction and Preliminaries
Difference equations play an important role in many disciplines including biology, ecology, physics, economics, and many more [1,2].These equations appear naturally as discrete analogs of differential equations, and as numerical solutions of differential and delay differential equations.Although most of these equations have a simple form, it is extremely difficult to understand their local and global dynamics including the existence, the boundedness, the asymptotic local stability, and the periodicity of solutions.Recently, there has been an increased interest in the quantitative as well as qualitative analysis of difference equations.One type of such difference equations that have been considered and examined by researchers is the nonlinear rational difference equations [3][4][5].The study of rational difference equations is of crucial importance, since we know so little about such equations.
Amleh et al. [6] considered the nonlinear difference equation x n = α + x n−1 /x n , where the parameter α and the initial condition x 0 are arbitrary positive real number.They proved that, if x n is a nontrivial solution of the equation so that there is a n 0 such that x n ≥ α + 1 for n ≥ n 0 , then x n is monotonically convergent to zero.In [7], Cinar studied the periodicity and the positive solution of the rational difference equations system x n+1 = 1/y n , y n+1 = y n /(x n−1 y n−1 ).The global character of solutions of the the system x n+1 = (h + x n )/(a + y n ), y n+1 = y n /(b + x n ) with positive parameters is studied in [8].Kurbanli et al. [9] investigated the positive solutions of the system of difference equations x n+1 = x n−1 /(y n x n−1 + 1), y n+1 = y n−1 /(x n y n−1 + 1).Yalcinkaya et al. [10] investigated the sufficient conditions for the global asymptotic stability of the following system x n+1 = (x n + y n−1 )/(x n y n−1 − 1), y n+1 = (y n + x n−1 /(y n x n−1 − 1).Khyat [11] studied many types of bifurcations for systems of several systems of rational difference equations.Many papers deal with the difference equations system.We refer in particular to [11][12][13][14] where further references can be found.
Our goal in this paper is to investigate the explicit solution of the following system of nonlinear rational difference equations systems where the parameter r and the initial conditions x −1 , x 0 , y 0 are arbitrary positive real numbers.Moreover, we investigate the global dynamics that occur in Equation (1).We show that the system in Equation (1) can have three fixed points and, using linear stability analysis, all transitions are determined under which fixed points lose stability.The bifurcations that occur in Equation ( 1) are numerically studied with a one-parameter bifurcation analysis.

The General Solution
Assume that the initial values are x −1 = b, x 0 = a, y 0 = c, then substituting n = 0, 1, 2, 3, . . .into Equation (1) gives the following expressions for x n and y n : where k = 0, 1, 2, 3, . .., W 0 (a) = 1, W k (a) = a + rW k−1 (a), and the P i represents the denominator for y i after ith iteration where P 0 = 1 and P 1 = bc + r.Formula for the sequence {P i } is a function of initial conditions.We use the above expressions to find On the other hand, by using Equations ( 1) and (2), we obtain from Comparing Equations ( 3) and ( 4), we find This means that We summarize this in the following Lemma.
Lemma 1.The following relations hold for k ≥ 0.
Now, we consider a special case x −1 = x 0 = b, y 0 = c and r = 1 in order to simplify the computations.In this case, we are able to calculate the quantities in Lemma 1.By using the properties of gamma function, it is easy to prove the following Lemma.Lemma 2. For every positive real p and q, the following relation holds: According to Lemma 2, we write , Similarly, we write where F 2 1 is the hypergeometric function.We introduce the notation Hence, we conclude We use the notation v(l) as the greatest integer function.We reach the following proposition.
Proposition 1.For x −1 = x 0 = b, y 0 = c, r = 1 and n ≥ 4, the general form of the solution of Equation ( 1) is Proof.By using Lemma 1, we have on the other hand, by using the expressions for I k + J k and Q m l , we obtain and Now, we obtain Similarly, we obtain and Hence, From Equations ( 5) and ( 6), we are done since Example 1.If we take a = b = 1 3 , c = 3 and r = 1, then, by direct computations, we obtain It is easy to see that P 1 = 2. Now, we have ,

Bifurcation Analysis
Consider the nested system in Equation (1).We introduce a new variable z, which is defined as z n = x n−1 (z n+1 = x n ), and then the system in Equation ( 1) can be expressed by the following three-dimensional system of difference equations: where n = 0, 1, 2, . .., X n = (x n , y n , z n ), F = z n z n +r , z n y n z n y n +r , x n , r ∈ R, and the initial point X 0 is positive.The fixed points of Equation ( 7) can be found by solving the system F(X * , r) − X * = 0. Thus, the system in Equation ( 7) can have three fixed points It is known that a fixed point X * is (asymptotically) stable if all eigenvalues of the Jacobian matrix of Equation ( 7) evaluated at X * have absolute values less than one.Moreover, if at least one of the eigenvalues has absolute value greater than one, then the fixed point X * is unstable [5,15].The Jacobian matrix of Equation ( 7) evaluated at the fixed point X i , having the eigenvalues µ 1 = The eigenvalues evaluated at the fixed points are summarized in Table 1.
To study the type of bifurcations that occur in Equation ( 7), we need to look for parameter values for which Equation ( 7) has a fixed point with eigenvalues satisfying one of the bifurcation conditions (for more details, see [16,17]).Proposition 3. (i) The fixed points X 1 and X 3 lose stability via a 1:4 resonance Neimark-Sacker bifurcation point (NS) when r = −1 for which µ 2,3 = e ±iθ 0 , θ 0 = π 2 .(ii) The fixed points X 2 and X 3 lose stability via a branching point (BP) when r = 1 2 for which µ 1 = 1.The BP corresponds to the appearance where the curves of the fixed points X 2 and X 3 are intersect.This indicates the existence of other cycles with higher periods.(iii) The fixed points X 1 and X 2 lose stability via a fold-flip period-doubling bifurcation point (PD) when r = 1 for which µ 2 = 1 and µ 3 = −1.(iv) We can always find a neutral-saddle point (nS) points for which µ i µ j = 1 along the curves of the fixed points (a) X 2 when r = {−1, 0.3819660112, 2.618033989}.(b) X 3 when r = 0.5698402917.
The bifurcation diagram can be constructed by repeating the system in Equation ( 7), e.g., k-times, for different r and using the fixed points as initial point.Figure 1a plots the bifurcation diagram of the system in Equation ( 7) with r ∈ [−1.5, 3].For each r, the initial points were reset to We performed 10 3 iterations and then we discarded the first 500 iterations to avoid transit effects.This led to the bifurcation diagram in Figure 1a.This diagram is refined by a continuation method (by using the MATLAB package MatContM, see [18][19][20]).The result of the continuation method is shown in Figure 1b.In Figure 2a,b, we present the values for the eigenvalues during fixed point continuations using the three fixed points of Equation (7).This confirms the theoretical results derived above.7) in (r, y)-plane: (a) the nested system in Equation ( 7) is repeated 10 5 times for different r and using the fixed points as initial point; and (b) the nested system in Equation ( 7) is repeated 10 5 times for different r and using continuation method (i.e., MatContM).The bifurcation points labeled in green correspond to the fixed point curve computed along X 1 , the red points correspond to the fixed point curve along X 2 , while the blue ones correspond to the fixed point curve along X 3 .7): (a) the eigenvalues corresponding to X 1 ; and (b) the eigenvalues corresponding to X 2 (dashed) and to X 3 (solid).

Figure 1 .
Figure 1.The bifurcation diagram for the system in Equation (7) in (r, y)-plane: (a) the nested system in Equation (7) is repeated 10 5 times for different r and using the fixed points as initial point; and (b) the nested system in Equation (7) is repeated 10 5 times for different r and using continuation method (i.e., MatContM).The bifurcation points labeled in green correspond to the fixed point curve computed along X 1 , the red points correspond to the fixed point curve along X 2 , while the blue ones correspond to the fixed point curve along X 3 . r

Figure 2 .
Figure 2.The eigenvalues values during fixed point continuation using the three fixed points of Equation (7): (a) the eigenvalues corresponding to X 1 ; and (b) the eigenvalues corresponding to X 2 (dashed) and to X 3 (solid).