# Global Dynamics of a 3 × 6 System of Difference Equations

Authors: S. M. Qureshi; A. Q. Khan
Journal: Discrete Dynamics in Nature and Society
 (2019)

---

## Abstract

In the proposed work, global dynamics of a3×6 system of rational difference equations has been studied in the interior of R+3. It is proved that system has at least one and at most seven boundary equilibria and a unique +ve equilibrium under certain parametric conditions. By utilizing method of Linearization, local dynamical properties about equilibria have been investigated. It is shown that every +ve solution of the system is bounded, and equilibrium P0 becomes a globally asymptotically stable if α1<α2,α4<α5, α7<α8. It is also shown that every +ve solution of the system converges to P0. Finally theoretical results are verified numerically.


---

## Body

## 1. Introduction

The importance of difference equations cannot be overemphasized. These equations model discrete physical phenomena on one hand and on the other hand these integral parts of numerical schemes used to solve differential equations. This widens the applicability of such equations to many branches of scientific knowledge. Discrete dynamical systems are described by difference equations and have applications in many branches of modern science like population dynamics, ecology, psychology, economics, queuing problems, genetics in biology, statistical problems, electrical networks, number theory, neural networks, quanta in radiation, sociology, physics, engineering, economics, probability theory, combinatorial analysis, stochastic time series, geometry, and resource management [1, 2]. Though the study of the global dynamics of higher-order nonlinear difference equations or system of difference equations is a challenging task but it is rewarding. These results pave path towards the development of basic theory of difference equations of higher-order. In recent years several authors have explored the behavior of solution of such difference equations or system of difference equations by studying equilibrium point, local and global dynamics about equilibria, boundedness and persistence, periodicity nature, prime period 2-solution, semicycle analysis, forbidden set, and many more (see [3–19] and references cited therein). For instance, Gibbons et al. [20] have explored the global dynamical properties of following difference equation:(1)xn+1=α1+α2xn-1α3+xn,n=0,1,…,where αj(j=1,2,3) and x-j(j=1,0) are +ve real numbers. Çinar [21] has explored the +ve solution of following difference equation:(2)xn+1=xn-11+∏i=01xn-i,n=0,1,…,where x-j(j=1,0) are +ve real numbers. Amleh et al. [22, 23] have explored boundedness, periodicity nature, and global dynamical properties of following difference equation:(3)xn+1=α1+α2xnxn-1+α3xn-1α4+α5xnxn-1+α6xn-1,n=0,1,…,where αj(j=1,…,6) and x-j(j=1,0) are +ve real numbers. Amleh et al. [22, 23] have pointed that difference equation depicted in (3) consists of numerous challenging and inspiring special cases having order 2, that can also arise from the following planer system when it is reduced to a difference equation (see [24]):(4)xn+1=α1+α2ynyn,yn+1=α3+α4xn+α5ynα6+α7xn+α8yn,n=0,1,…,where αj(j=1,…,8) and x0,y0 are +ve real numbers. Amleh et al. [22, 23] have also pointed that when one or more αj(j=1,…,6) allowed to zero then (3) contains 49 special cases in which 19 are riccati, linear, trivial, reducible to riccati, or linear difference equations while remaining 30 are posed to conjectures and open problems. Shojaei et al. [25] have extended the work of Çinar [21] to explore the dynamical properties of following difference equation:(5)xn+1=α1xn-2α2+α3xnxn-1xn-2,n=0,1,…,where αj(j=1,2,3) and x-j(j=2,1,0) are +ve real numbers. Later Bajo and Liz [26] have explored the dynamical properties of following difference equation which is special case of (3):(6)xn+1=xn-1α1+α2∏i=01xn-i,n=0,1,…,where αj(j=1,2) and x-j(j=1,0) are +ve real numbers. Zhang et al. [27] have extended the work studied by several authors [21, 25, 26], to explore the dynamical properties of following system of difference equations:(7)xn+1=xn-2α2+∏i=02yn-i,yn+1=yn-1α1+∏i=02xn-i,n=0,1,…,where αj(j=1,2) and xj,y-j(j=2,1,0) are +ve real numbers. Motivated from aforementioned studies, our aim is to extend the work studied by [21, 25–27], to explore the global dynamical properties of the following 3×6 difference equations system:(8)xn+1=α1xn-1α2+α3∏i=01yn-i,yn+1=α4yn-1α5+α6∏i=01zn-i,zn+1=α7zn-1α8+α9∏i=01xn-i,n=0,1,…,where αj(j=1,2,…,9) and x-j,y-j,z-j(j=1,0) are +ve real numbers.The flow pattern of the remaining paper, which is our main contribution and meaningful addition towards discrete dynamical systems described by the3×6 system of difference equations, is as follows: Section 2 deals with the study of existence of equilibria in R+3 whereas linearized form of the under consideration system is given in Section 3. Section 4 deals with the study of local dynamical properties about equilibria. Boundeness of the solution of the system is discussed in Section 5 whereas global dynamics about P0 of the system is studied in Section 6. In Section 7, it is proved that every +ve solution of under consideration system converges to P0. Discussion and simulations are given in the last section.


## 2. Existence of Equilibria

This section deals with the study of existence of equilibria of the system (8) in R+3. The result about the existence of equilibria can be summarized as follows.Lemma 1.
System (8) has at least one and at most eight equilibria in R+3. More precisely (i)
system (8) has a unique boundary equilibrium P0=(0,0,0)∀αj(j=1,2,…,9);(ii)
system (8) has a boundary equilibrium P1=α7-α8/α9,0,α4-α5/α6 if α1=α2,α4>α5andα7>α8;(iii)
system (8) has a boundary equilibrium P2=α7-α8/α9,α1-α2/α3,0 if α1>α2,α4=α5andα7>α8;(iv)
system (8) has a boundary equilibrium P3=0,α1-α2/α3,α4-α5/α6 if α1>α2,α4>α5andα7=α8;(v)
system (8) has a boundary equilibrium P4=α7-α8/α9,0,0 if α1=α2,α4=α5andα7>α8;(vi)
system (8) has a boundary equilibrium P5=0,0,α4-α5/α6 if α1=α2,α4>α5andα7=α8;(vii)
system (8) has a boundary equilibrium P6=0,α1-α2/α3,0 if α1>α2,α4=α5andα7=α8;(viii)
system (8) has a +ve equilibrium P7=α7-α8/α9,α1-α2/α3,α4-α5/α6 if α1>α2,α4>α5andα7>α8. More specifically, P7=α7-α8/α9,α1-α2/α3,α4-α5/α6 is a unique +ve equilibrium of the system (8) if α1>α2,α4>α5andα7>α8.


## 3. Linearized Form of (8)

We have the following map in order to construct the corresponding linearized form of (8):(9)Ω1,Ω2,Ω3,Ω4,Ω5,Ω6→xn+1,xn,yn+1,yn,zn+1,zn,where(10)Ω1=α1xn-1α2+α3∏i=01yn-i,Ω2=xn,Ω3=α4yn-1α5+α6∏i=01zn-i,Ω4=yn,Ω5=α7zn-1α8+α9∏i=01xn-i,Ω6=zn. Moreover JΛ about Λ under the map depicted in (9) is(11)JΛ=0α1α2+α3y2-α1α3xyα2+α3y22-α1α3xyα2+α3y2200100000000α4α5+α6z2-α4α6yzα5+α6z22-α4α6yzα5+α6z22001000-α7α9xzα8+α9x22-α7α9xzα8+α9x22000α7α8+α9x2000010,where Λ=(x,y,z).


## 4. Local Stability of System (8)

Hereafter we will study the local dynamical properties about equilibria of the system (8) by utilizing Theorem 1.5 of [1].


### 4.1. Local Stability about Boundary Equilibria

This subsection is purely dedicated to the study of local dynamical properties about boundary equilibria:P0,P1,P2,P3,P4,P5,P6 of the system (8), respectively.Theorem 2.
ForP0, the following hold: (i)
P0 of the system (8) is a sink if(12)α1<α2,α4<α5,α7<α8;(ii)
P0 of the system (8) is unstable if(13)α1>α2orα4>α5orα7>α8.Proof.
(i) We have the following linearized system of (8) about P0:(14)Υn+1=BΥn,
where(15)Υn=xnxn-1ynyn-1znzn-1,B=JP0=0α1α20000100000000α4α50000100000000α7α8000010. Let the eigenvalues of B are ν1,ν2,⋯,ν6. Also let(16)D=d1000000d2000000d3000000d4000000d5000000d6,denoting the diagonal matrix where(17)d1=d3=d5=1,d2=d4=d6=1-ϵfor0<ϵ<1,and(18)0<ϵ<min⁡1-α1α2,1-α4α5,1-α7α8.By computing DBD-1 one gets (19)DBD-1=0d1d2-1α1α20000d2d1-100000000d3d4-1α4α50000d4d3-100000000d5d6-1α7α80000d6d5-10. Moreover,(20)d1>d2>0,d3>d4>0,d5>d6>0.Equation (20) then implies that(21)d2d1-1<1,d4d3-1<1,d6d5-1<1.From (17) and (18) one gets(22)d1d2-1α1α2=d2-1α1α2=11-ϵα1α2<1,d3d4-1α4α5=d4-1α4α5=11-ϵα4α5<1,d5d6-1α7α8=d6-1α7α8=11-ϵα7α8<1.Since eigenvalues of B are same as DBD-1, so (23)max1≤m≤6νm=DBD-1=max⁡d2d1-1,d4d3-1,d6d5-1,d1d2-1α1α2,d3d4-1α4α5,d5d6-1α7α8<1. Thus P0 of (8) is a sink.
(ii) Similarly it is easy to prove thatP0 of (8) is unstable if α1>α2 or α4>α5 or α7>α8.Theorem 3.
P 1 of (8) is locally unstable.Proof.
The corresponding linearized form of (8) about P1 is(24)Υn+1=B1Υn, where(25)B1=JP1=0α1α20000100000000100001000-α9α7α4-α5α6α7-α8α9-α9α7α4-α5α6α7-α8α90001000010. Again let the eigenvalues of B1 be ν1,ν2,⋯,ν6 and the diagonal matrix D which is depicted in (16) where (17) holds. On computing DB1D-1 one gets(26)DB1D-1=0d1d2-1α1α20000d2d1-100000000d3d4-10000d4d3-1000-α9α7α4-α5α6α7-α8α9d5d1-1-α9α7α4-α5α6α7-α8α9d5d2-1000d5d6-10000d6d5-10. In view of (20) and (21) one has(27)d1d2-1α1α2=d2-1α1α2=11-ϵα1α2<1, which implies that ϵ<1-α1/α2. Also(28)d3d4-1=d4-1=11-ϵ>1. Hence P1 of the system (8) is unstable.Theorem 4.
P 2 of (8) is locally unstable.Proof.
The corresponding linearized form of (8) about P2 is(29)Υn+1=B2Υn, where(30)B2=JP2=01-α3α1α1-α2α3α7-α8α9-α3α1α1-α2α3α7-α8α900100000000α4α500001000000001000010. Using same arrangement as in the proof of Theorems 2–3 and compute DB2D-1; one gets(31)DB2D-1=0d1d2-1-α3α1α1-α2α3α7-α8α9d1d3-1-α3α1α1-α2α3α7-α8α9d1d4-100d2d1-100000000α4α5d3d4-10000d4d3-100000000d5d6-10000d6d5-10. In view of (17), (20), and (21) one has(32)
d
1
d
2
-
1
+
α
3
α
1
α
1
-
α
2
α
3
α
7
-
α
8
α
9
d
1
d
3
-
1
+
α
3
α
1
α
1
-
α
2
α
3
α
7
-
α
8
α
9
d
1
d
4
-
1
=
α
3
α
1
α
1
-
α
2
α
3
α
7
-
α
8
α
9
+
1
+
α
3
α
1
α
1
-
α
2
α
3
α
7
-
α
8
α
9
1
1
-
ϵ
>
1
. And(33)α4α5d3d4-1=α4α511-ϵ<1, which implies that ϵ<1-α4/α5. Also(34)d5d6-1=d6-1=11-ϵ>1. Hence P2 of (8) is locally unstable.Theorem 5.
P 3 of the system (8) is locally unstable.Proof.
Its proof is similar to Theorems3–4.Theorem 6.
P 4 of (8) is locally unstable.Proof.
The corresponding linearized system is(35)Υn+1=B3Υn, where(36)B3=JP4=0α1α20000100000000α4α500001000000001000010. Again using same arrangements as in the proof of Theorems 2–4 and computing DB3D-1, one gets(37)DB3D-1=0α1α2d1d2-10000d2d1-100000000α4α5d3d4-10000d4d3-100000000d5d6-10000d6d5-10. Further,(38)α1α2d1d2-1=11-ϵ>1,(39)α4α5d3d4-1=11-ϵ>1, and(40)d6d5-1=11-ϵ>1. Hence P4 of the system (8) is locally unstable.Theorem 7.
P 5 and P6 of the system (8) are locally unstable.Proof.
Similar to Theorem6.


### 4.2. Local Stability about the Unique Positive Equilibrium

Hereafter local dynamical properties aboutP7 of (8) are explored. The following result shows that P7 of system (8) is locally unstable.Theorem 8.
P 7 of (8) is unstable.Proof.
The corresponding linearized form of (8) about P7 is(41)Υn+1=B4=JP7Υn, where(42)B4=01-α3α1α1-α2α3α7-α8α9100000001-α9α7α4-α5α6α7-α8α9-α9α7α4-α5α6α7-α8α90000-α3α1α1-α2α3α7-α8α9000001-α6α4α1-α2α3α4-α5α6-α6α4α1-α2α3α4-α5α6000001010. Also(43)DB4D-1=0d1d2-1-α3α1α1-α2α3α7-α8α9d1d3-1d2d1-10000000d4d3-1-α9α7α4-α5α6α7-α8α9d5d1-1-α9α7α4-α5α6α7-α8α9d5d2-10000-α3α1α1-α2α3α7-α8α9d1d4-100000d3d4-1-α6α4α1-α2α3α4-α5α6d3d5-1-α6α4α1-α2α3α4-α5α6d3d6-100000d5d6-10d6d5-10. Moreover(44)d1d2-1+α3α1α1-α2α3α7-α8α9d1d3-1+α3α1α1-α2α3α7-α8α9d1d4-1=α3α1α1-α2α3α7-α8α9+1+α3α1α1-α2α3α7-α8α911-ϵ>1,(45)d3d4-1+α6α4α1-α2α3α4-α5α6d3d5-1+α6α4α1-α2α3α4-α5α6d3d6-1=α6α4α1-α2α3α4-α5α6+1+α6α4α1-α2α3α4-α5α611-ϵ>1, and (46)α9α7α4-α5α6α7-α8α9d5d1-1+α9α7α4-α5α6α7-α8α9d5d2-1+d5d6-1=α9α7α4-α5α6α7-α8α9+1+α9α7α4-α5α6α7-α8α911-ϵ>1. Hence P7 of the system (8) is locally unstable.


## 4.1. Local Stability about Boundary Equilibria

This subsection is purely dedicated to the study of local dynamical properties about boundary equilibria:P0,P1,P2,P3,P4,P5,P6 of the system (8), respectively.Theorem 2.
ForP0, the following hold: (i)
P0 of the system (8) is a sink if(12)α1<α2,α4<α5,α7<α8;(ii)
P0 of the system (8) is unstable if(13)α1>α2orα4>α5orα7>α8.Proof.
(i) We have the following linearized system of (8) about P0:(14)Υn+1=BΥn,
where(15)Υn=xnxn-1ynyn-1znzn-1,B=JP0=0α1α20000100000000α4α50000100000000α7α8000010. Let the eigenvalues of B are ν1,ν2,⋯,ν6. Also let(16)D=d1000000d2000000d3000000d4000000d5000000d6,denoting the diagonal matrix where(17)d1=d3=d5=1,d2=d4=d6=1-ϵfor0<ϵ<1,and(18)0<ϵ<min⁡1-α1α2,1-α4α5,1-α7α8.By computing DBD-1 one gets (19)DBD-1=0d1d2-1α1α20000d2d1-100000000d3d4-1α4α50000d4d3-100000000d5d6-1α7α80000d6d5-10. Moreover,(20)d1>d2>0,d3>d4>0,d5>d6>0.Equation (20) then implies that(21)d2d1-1<1,d4d3-1<1,d6d5-1<1.From (17) and (18) one gets(22)d1d2-1α1α2=d2-1α1α2=11-ϵα1α2<1,d3d4-1α4α5=d4-1α4α5=11-ϵα4α5<1,d5d6-1α7α8=d6-1α7α8=11-ϵα7α8<1.Since eigenvalues of B are same as DBD-1, so (23)max1≤m≤6νm=DBD-1=max⁡d2d1-1,d4d3-1,d6d5-1,d1d2-1α1α2,d3d4-1α4α5,d5d6-1α7α8<1. Thus P0 of (8) is a sink.
(ii) Similarly it is easy to prove thatP0 of (8) is unstable if α1>α2 or α4>α5 or α7>α8.Theorem 3.
P 1 of (8) is locally unstable.Proof.
The corresponding linearized form of (8) about P1 is(24)Υn+1=B1Υn, where(25)B1=JP1=0α1α20000100000000100001000-α9α7α4-α5α6α7-α8α9-α9α7α4-α5α6α7-α8α90001000010. Again let the eigenvalues of B1 be ν1,ν2,⋯,ν6 and the diagonal matrix D which is depicted in (16) where (17) holds. On computing DB1D-1 one gets(26)DB1D-1=0d1d2-1α1α20000d2d1-100000000d3d4-10000d4d3-1000-α9α7α4-α5α6α7-α8α9d5d1-1-α9α7α4-α5α6α7-α8α9d5d2-1000d5d6-10000d6d5-10. In view of (20) and (21) one has(27)d1d2-1α1α2=d2-1α1α2=11-ϵα1α2<1, which implies that ϵ<1-α1/α2. Also(28)d3d4-1=d4-1=11-ϵ>1. Hence P1 of the system (8) is unstable.Theorem 4.
P 2 of (8) is locally unstable.Proof.
The corresponding linearized form of (8) about P2 is(29)Υn+1=B2Υn, where(30)B2=JP2=01-α3α1α1-α2α3α7-α8α9-α3α1α1-α2α3α7-α8α900100000000α4α500001000000001000010. Using same arrangement as in the proof of Theorems 2–3 and compute DB2D-1; one gets(31)DB2D-1=0d1d2-1-α3α1α1-α2α3α7-α8α9d1d3-1-α3α1α1-α2α3α7-α8α9d1d4-100d2d1-100000000α4α5d3d4-10000d4d3-100000000d5d6-10000d6d5-10. In view of (17), (20), and (21) one has(32)
d
1
d
2
-
1
+
α
3
α
1
α
1
-
α
2
α
3
α
7
-
α
8
α
9
d
1
d
3
-
1
+
α
3
α
1
α
1
-
α
2
α
3
α
7
-
α
8
α
9
d
1
d
4
-
1
=
α
3
α
1
α
1
-
α
2
α
3
α
7
-
α
8
α
9
+
1
+
α
3
α
1
α
1
-
α
2
α
3
α
7
-
α
8
α
9
1
1
-
ϵ
>
1
. And(33)α4α5d3d4-1=α4α511-ϵ<1, which implies that ϵ<1-α4/α5. Also(34)d5d6-1=d6-1=11-ϵ>1. Hence P2 of (8) is locally unstable.Theorem 5.
P 3 of the system (8) is locally unstable.Proof.
Its proof is similar to Theorems3–4.Theorem 6.
P 4 of (8) is locally unstable.Proof.
The corresponding linearized system is(35)Υn+1=B3Υn, where(36)B3=JP4=0α1α20000100000000α4α500001000000001000010. Again using same arrangements as in the proof of Theorems 2–4 and computing DB3D-1, one gets(37)DB3D-1=0α1α2d1d2-10000d2d1-100000000α4α5d3d4-10000d4d3-100000000d5d6-10000d6d5-10. Further,(38)α1α2d1d2-1=11-ϵ>1,(39)α4α5d3d4-1=11-ϵ>1, and(40)d6d5-1=11-ϵ>1. Hence P4 of the system (8) is locally unstable.Theorem 7.
P 5 and P6 of the system (8) are locally unstable.Proof.
Similar to Theorem6.


## 4.2. Local Stability about the Unique Positive Equilibrium

Hereafter local dynamical properties aboutP7 of (8) are explored. The following result shows that P7 of system (8) is locally unstable.Theorem 8.
P 7 of (8) is unstable.Proof.
The corresponding linearized form of (8) about P7 is(41)Υn+1=B4=JP7Υn, where(42)B4=01-α3α1α1-α2α3α7-α8α9100000001-α9α7α4-α5α6α7-α8α9-α9α7α4-α5α6α7-α8α90000-α3α1α1-α2α3α7-α8α9000001-α6α4α1-α2α3α4-α5α6-α6α4α1-α2α3α4-α5α6000001010. Also(43)DB4D-1=0d1d2-1-α3α1α1-α2α3α7-α8α9d1d3-1d2d1-10000000d4d3-1-α9α7α4-α5α6α7-α8α9d5d1-1-α9α7α4-α5α6α7-α8α9d5d2-10000-α3α1α1-α2α3α7-α8α9d1d4-100000d3d4-1-α6α4α1-α2α3α4-α5α6d3d5-1-α6α4α1-α2α3α4-α5α6d3d6-100000d5d6-10d6d5-10. Moreover(44)d1d2-1+α3α1α1-α2α3α7-α8α9d1d3-1+α3α1α1-α2α3α7-α8α9d1d4-1=α3α1α1-α2α3α7-α8α9+1+α3α1α1-α2α3α7-α8α911-ϵ>1,(45)d3d4-1+α6α4α1-α2α3α4-α5α6d3d5-1+α6α4α1-α2α3α4-α5α6d3d6-1=α6α4α1-α2α3α4-α5α6+1+α6α4α1-α2α3α4-α5α611-ϵ>1, and (46)α9α7α4-α5α6α7-α8α9d5d1-1+α9α7α4-α5α6α7-α8α9d5d2-1+d5d6-1=α9α7α4-α5α6α7-α8α9+1+α9α7α4-α5α6α7-α8α911-ϵ>1. Hence P7 of the system (8) is locally unstable.


## 5. Boundedness

Theorem 9.
Any positive solution{(xn,yn,zn)} of (8) satisfies following inequalities for m≥0:(47)0<xn≤α1α2m+1x-1,ifn=2m+1,0<xn≤α1α2m+1x0,ifn=2m+2.And(48)0<yn≤α4α5m+1y-1,ifn=2m+1,0<yn≤α4α5m+1y0,ifn=2m+2.Finally,(49)0<zn≤α7α8m+1z-1,ifn=2m+1,0<zn≤α7α8m+1z0,ifn=2m+2.Proof.
Obviously (47), (48), and (49) hold for m=0. Now suppose that for m=p≥1 inequalities (47), (48), and (49) are true, i.e.,(50)0<xn≤α1α2p+1x-1,ifn=2p+1,0<xn≤α1α2p+1x0,ifn=2p+2.And(51)0<yn≤α4α5p+1y-1,ifn=2p+1,0<yn≤α4α5p+1y0,ifn=2p+2.Finally,(52)0<zn≤α7α8p+1z-1,ifn=2p+1,0<zn≤α7α8p+1z0,ifn=2p+2.Now, for m=p+1 using (8) one has(53)0<x2p+3=α1x2p+1α2+α3y2p+2y2p+1≤α1α2x2p+1≤α1α2p+2x-1,0<x2p+4=α1x2p+2α2+α3y2p+3y2p+2≤α1α2x2p+2≤α1α2p+2x0.And(54)0<y2p+3=α4y2p+1α5+α6z2p+2z2p+1≤α4α5y2p+1≤α4α5p+2y-1,0<y2p+4=α4y2p+2α5+α6z2p+3z2p+2≤α4α5y2p+2≤α4α5p+2y0.Finally,(55)0<z2p+3=α7z2p+1α8+α9x2p+2x2p+1≤α7α8z2p+1≤α7α8p+2z-1,0<z2p+3=α7z2p+2α8+α9x2p+3x2p+2≤α7α8z2p+2≤α7α8p+2z0.Lemma 10.
If the conditions depicted in (12) hold then {(xn,yn,zn)} of (8) is bounded.Proof.
Direct consequence of Theorem9.


## 6. Global Stability about the EquilibriumP0

Now we will explore the global dynamics of (8) about P0.Theorem 11.
If the conditions depicted in (12) hold then P0 of (8) is globally asymptotically stable.Proof.
SinceP0 of (8) is a sink by Theorem 2 and hence it is enough to prove(56)limn→∞xn,yn,zn=P0. From (8), one gets(57)0<xn+1=α1xn-1α2+α3∏i=01yn-i≤α1α2xn-1,0<yn+1=α4yn-1α5+α6∏i=01zn-i≤α4α5yn-1,0<zn+1=α7zn-1α8+α9∏i=01xn-i≤α7α8zn-1. By induction, one gets(58)x2n-1≤α1α2nx-1,x2n≤α1α2nx0,y2n-1≤α4α5ny-1,y2n≤α4α5ny0,z2n-1≤α7α8nz-1,z2n≤α7α8nz0.Finally if (12) holds then from (58) one gets limn→∞(xn,yn,zn)=P0.


## 7. Rate of Convergence

Theorem 12.
If (12) holds then error vector(59)en=en1en-11en2en-12en3en-13=xn-xxn-1-xyn-yyn-1-yzn-zzn-1-z, of every +ve solution {(xn,yn,zn)} of (8) satisfies both of the following asymptotic relations:(60)limn→∞enn=ν1,…,6JP0,limn→∞en+1en=ν1,…,6JP0, where |ν1,…,6JP0| is equal to the modulus of one of the eigenvalues of the Jacobian matrix evaluated at P0.Proof.
Let{(xn,yn,zn)} be an arbitrary solution of (8), s.t., limn→∞xn=x,limn→∞yn=y and limn→∞zn=z. To find error term one has from the system (8) (61)xn+1-x=α1xn-1α2+α3∏i=01yn-i-α1xα2+α3y2,=α1α2+α3∏i=01yn-ixn-1-x-α1α3xyα2+α3∏i=01yn-iα2+α3y2yn-y-α1α3xynα2+α3∏i=01yn-iα2+α3y2yn-1-y,yn+1-y=α4yn-1α5+α6∏i=01zn-i-α4yα5+α6z2,=α4α5+α6∏i=01zn-iyn-1-y-α4α6yzα5+α6∏i=01zn-iα5+α6z2zn-z-α4α6yznα5+α6∏i=01zn-iα5+α6z2zn-1-z, and (62)zn+1-z=α7zn-1α8+α9∏i=01xn-i-α7yα8+α9x2,=-α7α9xzα8+α9∏i=01xn-iα8+α9x2xn-x-α7α9xnzα8+α9∏i=01xn-iα8+α9x2xn-1-x+α7α8+α9∏i=01xn-izn-1-z. So(63)xn+1-x=α1α2+α3∏i=01yn-ixn-1-x-α1α3xyα2+α3∏i=01yn-iα2+α3y2yn-y-α1α3xynα2+α3∏i=01yn-iα2+α3y2yn-1-y,yn+1-y=α4α5+α6∏i=01zn-iyn-1-y-α4α6yzα5+α6∏i=01zn-iα5+α6z2zn-z-α4α6yznα5+α6∏i=01zn-iα5+α6z2zn-1-z,zn+1-z=-α7α9xzα8+α9∏i=01xn-iα8+α9x2xn-x-α7α9xnzα8+α9∏i=01xn-iα8+α9x2xn-1-x+α7α8+α9∏i=01xn-izn-1-z.Set(64)en1=xn-x,en2=yn-y,en3=zn-z.In view of (64), (63) becomes(65)en+11=Λn1en-11+Λn2en2+Λn3en-12,en+12=Λn4en-12+Λn5en3+Λn6en-13,en+13=Λn7en1+Λn8en-11+Λn9en-13, where (66)Λn1=α1α2+α3∏i=01yn-i,Λn2=-α1α3xyα2+α3∏i=01yn-iα2+α3y2,Λn3=-α1α3xynα2+α3∏i=01yn-iα2+α3y2,Λn4=α4α5+α6∏i=01zn-i,Λn5=-α4α6yzα5+α6∏i=01zn-iα5+α6z2,Λn6=-α4α6yznα5+α6∏i=01zn-iα5+α6z2,Λn7=-α7α9xzα8+α9∏i=01xn-iα8+α9x2,Λn8=-α7α9xnzα8+α9∏i=01xn-iα8+α9x2,Λn9=α7α8+α9∏i=01xn-i. Moreover (67)limn→∞Λn1=α1α2+α3y2,limn→∞Λn2=-α1α3xyα2+α3y22,limn→∞Λn3=-α1α3xyα2+α3y22,limn→∞Λn4=α4α5+α6z2,limn→∞Λn5=-α4α6yzα5+α6z22,limn→∞Λn6=-α4α6yzα5+α6z22,limn→∞Λn7=-α7α9xzα8+α9x22,limn→∞Λn8=-α7α9xzα8+α9x22,limn→∞Λn9=α7α8+α9x2, that is, (68)Λn1=α1α2+α3y2+ϱ1,Λn2=-α1α3xyα2+α3y22+ϱ2,Λn3=-α1α3xyα2+α3y22+ϱ3,Λn4=α4α5+α6z2+ϱ4,Λn5=-α4α6yzα5+α6z22+ϱ5,Λn6=-α4α6yzα5+α6z22+ϱ6,Λn7=-α7α9xzα8+α9x22+ϱ7,Λn8=-α7α9xzα8+α9x22+ϱ8,Λn9=α7α8+α9x2+ϱ9,
whereϱ1,…,ϱ9→0 as n→∞. Thus we have the system [28]:(69)en+1=A+Bnen,where(70)A=0α1α2+α3y2-α1α3xyα2+α3y22-α1α3xyα2+α3y2200100000000α4α5+α6z2-α4α6yzα5+α6z22-α4α6yzα5+α6z22001000-α7α9xzα8+α9x22-α7α9xzα8+α9x22000α7α8+α9x2000010.(71)Bn=0ϱ1ϱ2ϱ300100000000ϱ4ϱ5ϱ6001000ϱ7ϱ8000ϱ9000010, and(72)Bn→0,n→∞. Therefore limiting system of error terms takes the following form:(73)en+11en1en+12en2en+13en3=JΛen1en-11en2en-12en3en-13.Hence (73) is same as linearized system of (8) about (x,y,z). Particularly about P0 it becomes(74)en+11en1en+12en2en+13en3=0α1α20000100000000α4α50000100000000α7α8000010en1en-11en2en-12en3en-13, which is same as linearized system of (8) about P0.


## 8. Discussion and Numerical Simulations

This proposed work is about the global dynamical properties of a3×6 system of difference equations, which is our key finding towards discrete dynamical system. We investigated that discrete-time system (8) has 8 equilibria: P0,P1,P2,P3,P4,P5,P6 and the unique +ve equilibrium point: P7 under restriction to the parameters αj(j=1,2,…,9). By method of linearization, we studied the local dynamical properties about each equilibria and conclusion are presented in Table 1. Further we proved that {(xn,yn,zn)} of (8) is bounded and P0 is globally asymptotically stable if (12) holds. Finally, we proved that every +ve solution {(xn,yn,zn)} of the system converges to P0. The above obtained theoretical results can also be verified from following numerical simulations. Thus in the remaining section, theoretical results are verified numerically.Table 1
Number of Equilibria along their qualitative behavior of (8).
E.P Corresponding behavior P 0 locally asymptotically stable if (12) hold; unstable if (13) hold. P 1 unstable. P 2 unstable. P 3 unstable. P 4 unstable. P 5 unstable. P 6 unstable. P 7 unstable.Example 1.
Ifα1=24,α2=27,α3=23,α4=14,α5=17,α6=13,α7=14,α8=16,α9=5 then (8) with xj,yj,zj(j=-1,0), respectively, are 3.7,9,0.9,0.4,0.9,0.4, the following form is taken:(75)xn+1=24xn-127+23∏i=01yn-i,yn+1=14yn-117+13∏i=01zn-i,zn+1=14zn-116+5∏i=01xn-i.The graphs of nvsxn(respectivelyynandzn) of (75) are shown in Figures 1(a), 1(b), and 1(c) while its attractor is shown in Figure 1(d); i.e., if α1=24<α2=27,α4=14<α5=17,α7=14<α8=16, then all solutions are attracted to P0. This coincides with the proof of Theorem 11.Figure 1
Phase portrait of (75).
(a) (b) (c) (d)Example 2.
Ifα1=27,α2=30,α3=0.3,α4=24,α5=37,α6=1.3,α7=4,α8=6,α9=0.5 then (8) with xj,yj,zj(j=-1,0), respectively, are 0.7,0.9,0.9,3.4,0.99999,1.4, the following form is taken:(76)xn+1=27xn-130+0.3∏i=01yn-i,yn+1=24yn-137+1.3∏i=01zn-i,zn+1=4zn-16+0.5∏i=01xn-i.The graphs of nvsxn(respectivelyynandzn) of (76) are shown in Figures 2(a), 2(b), and 2(c) while its attractor is shown in Figure 2(d).Figure 2
Phase portrait of (76).
(a) (b) (c) (d)Example 3.
Ifα1=28,α2=14,α3=3,α4=19,α5=7,α6=13,α7=4,α8=2,α9=5 then (8) with xj,yj,zj(j=-1,0), respectively, are 0.7,0.9,10.9,=0.4,0.99,1.4, the following form is taken:(77)xn+1=28xn-114+3∏i=01yn-i,yn+1=19yn-17+13∏i=01zn-i,zn+1=4zn-12+5∏i=01xn-i.The graphs of nvsxn(respectivelyynandzn) of (75) are shown in Figures 3(a), 3(b), and 3(c). In this case P0 of (77) is unstable.Figure 3
Phase portrait of (77).
(a) (b) (c)



---

*Source: 9797242-2019-07-01.xml*