Processing math: 100%
Research article Topical Sections

Assessment of offshore wind power potential in the Aegean and Ionian Seas based on high-resolution hindcast model results

  • In this study long-term wind data obtained from high-resolution hindcast simulations is used to analytically assess offshore wind power potential in the Aegean and Ionian Seas and provide wind climate and wind power potential characteristics at selected locations, where offshore wind farms are at the concept/planning phase. After ensuring the good model performance through detailed validation against buoy measurements, offshore wind speed and wind direction at 10 m above sea level are statistically analyzed on the annual and seasonal time scale. The spatial distribution of the mean wind speed and wind direction are provided in the appropriate time scales, along with the mean annual and the inter-annual variability; these statistical quantities are useful in the offshore wind energy sector as regards the preliminary identification of favorable sites for exploitation of offshore wind energy. Moreover, the offshore wind power potential and its variability are also estimated at 80 m height above sea level. The obtained results reveal that there are specific areas in the central and the eastern Aegean Sea that combine intense annual winds with low variability; the annual offshore wind power potential in these areas reach values close to 900 W/m2, suggesting that a detailed assessment of offshore wind energy would be worth noticing and could lead in attractive investments. Furthermore, as a rough estimate of the availability factor, the equiprobable contours of the event [4 m/s ≤ wind speed ≤ 25 m/s] are also estimated and presented. The selected lower and upper bounds of wind speed correspond to typical cut-in and cut-out wind speed thresholds, respectively, for commercial offshore wind turbines. Finally, for seven offshore wind farms that are at the concept/planning phase the main wind climate and wind power density characteristics are also provided.

    Citation: Takvor Soukissian, Anastasios Papadopoulos, Panagiotis Skrimizeas, Flora Karathanasi, Panagiotis Axaopoulos, Evripides Avgoustoglou, Hara Kyriakidou, Christos Tsalis, Antigoni Voudouri, Flora Gofa, Petros Katsafados. Assessment of offshore wind power potential in the Aegean and Ionian Seas based on high-resolution hindcast model results[J]. AIMS Energy, 2017, 5(2): 268-289. doi: 10.3934/energy.2017.2.268

    Related Papers:

    [1] Ya-Lei Li, Da-Bin Wang, Jin-Long Zhang . Sign-changing solutions for a class of p-Laplacian Kirchhoff-type problem with logarithmic nonlinearity. AIMS Mathematics, 2020, 5(3): 2100-2112. doi: 10.3934/math.2020139
    [2] Xia Li, Wen Guan, Da-Bin Wang . Least energy sign-changing solutions of Kirchhoff equation on bounded domains. AIMS Mathematics, 2022, 7(5): 8879-8890. doi: 10.3934/math.2022495
    [3] Jinfu Yang, Wenmin Li, Wei Guo, Jiafeng Zhang . Existence of infinitely many normalized radial solutions for a class of quasilinear Schrödinger-Poisson equations in R3. AIMS Mathematics, 2022, 7(10): 19292-19305. doi: 10.3934/math.20221059
    [4] Qing Yang, Chuanzhi Bai . Sign-changing solutions for a class of fractional Kirchhoff-type problem with logarithmic nonlinearity. AIMS Mathematics, 2021, 6(1): 868-881. doi: 10.3934/math.2021051
    [5] Ziqing Yuan, Jing Zhao . Solutions for gauged nonlinear Schrödinger equations on R2 involving sign-changing potentials. AIMS Mathematics, 2024, 9(8): 21337-21355. doi: 10.3934/math.20241036
    [6] Ye Xue, Zhiqing Han . Existence and multiplicity of solutions for Schrödinger equations with sublinear nonlinearities. AIMS Mathematics, 2021, 6(6): 5479-5492. doi: 10.3934/math.2021324
    [7] Jie Yang, Lintao Liu, Haibo Chen . Multiple positive solutions to the fractional Kirchhoff-type problems involving sign-changing weight functions. AIMS Mathematics, 2024, 9(4): 8353-8370. doi: 10.3934/math.2024406
    [8] Zhongxiang Wang, Gao Jia . Existence of solutions for modified Kirchhoff-type equation without the Ambrosetti-Rabinowitz condition. AIMS Mathematics, 2021, 6(5): 4614-4637. doi: 10.3934/math.2021272
    [9] Zhongxiang Wang . Existence and asymptotic behavior of normalized solutions for the modified Kirchhoff equations in R3. AIMS Mathematics, 2022, 7(5): 8774-8801. doi: 10.3934/math.2022490
    [10] Yuan Shan, Baoqing Liu . Existence and multiplicity of solutions for generalized asymptotically linear Schrödinger-Kirchhoff equations. AIMS Mathematics, 2021, 6(6): 6160-6170. doi: 10.3934/math.2021361
  • In this study long-term wind data obtained from high-resolution hindcast simulations is used to analytically assess offshore wind power potential in the Aegean and Ionian Seas and provide wind climate and wind power potential characteristics at selected locations, where offshore wind farms are at the concept/planning phase. After ensuring the good model performance through detailed validation against buoy measurements, offshore wind speed and wind direction at 10 m above sea level are statistically analyzed on the annual and seasonal time scale. The spatial distribution of the mean wind speed and wind direction are provided in the appropriate time scales, along with the mean annual and the inter-annual variability; these statistical quantities are useful in the offshore wind energy sector as regards the preliminary identification of favorable sites for exploitation of offshore wind energy. Moreover, the offshore wind power potential and its variability are also estimated at 80 m height above sea level. The obtained results reveal that there are specific areas in the central and the eastern Aegean Sea that combine intense annual winds with low variability; the annual offshore wind power potential in these areas reach values close to 900 W/m2, suggesting that a detailed assessment of offshore wind energy would be worth noticing and could lead in attractive investments. Furthermore, as a rough estimate of the availability factor, the equiprobable contours of the event [4 m/s ≤ wind speed ≤ 25 m/s] are also estimated and presented. The selected lower and upper bounds of wind speed correspond to typical cut-in and cut-out wind speed thresholds, respectively, for commercial offshore wind turbines. Finally, for seven offshore wind farms that are at the concept/planning phase the main wind climate and wind power density characteristics are also provided.


    The Lotka-Volterra equations that represent the system of continuous-time non-linear differential equations whose solution is deterministic, evolved to determine the dynamics of biological systems. This system is a more generalized version of the Kolmogorov model because it only focuses on the predator-prey interactions but ignores competition, disease and mutualism, which the Kolmogorov model includes. As compared with other mathematical model, following assumptions are made in order to formulate the Lotka-Volterra equations:

    (i) For the prey population there is no shortage of food.

    (ii) The amount of food supplied to the prey is directly related to the size of the prey population.

    (iii) The rate of change of population is directly proportional to its size.

    (iv) The environment is constant and genetic adaptation is not assumed to be negligible.

    (v) Predators will never stop eating.

    Based on said assumptions a continuous-time, two-species Lotka-Volterra system takes the following form [1,2,3]:

    dxdt=axbxy, dydt=cxydy, (1.1)

    where x and y respectively denotes the number of prey and predator. The quantities: dxdt and dydt respectively denotes the instantaneous rates of prey and predator, and a, b, c and d are positive real parameters. In 1st equation of system (1.1), ax represent that prey are reproduced exponentially whereas bxy shows that the rate at which the predator kills prey. So, if there is no population of prey or no population of predator, no decrease in the population of prey can occur. The equation for prey can be summed up as the rate at which new prey is born, minus the rate at which prey is killed off. Now in 2nd equation of system (1.1), cxy represent the rate at which predators kill prey except that a different constant is used to describe this relationship since the rate at which predators kill and the rate at which they reproduce are not identical. Also, dy is assumed to be an exponential decay, which mean that prey cannot kill the predators and that the decrease in predator population is due to natural death or emigration.

    It is anticipated that as compared to continuous-time system, discrete-time system designated by difference equations are widely explored because these models are more reasonable, if populations have non-overlapping generations, as compared to differential type models, and also one can obtain more accurate numerical simulations [4,5,6]. During the last few years, discrete-time predator-prey models have been extensively investigated by studying fixed points, local and global dynamical properties, existence of possible bifurcation and chaos [7,8,9] and references cited therein. So, in this work, we explore the local dynamics, bifurcations and hybrid control in a 3×3 discrete-time predator-prey model, i.e., two prey and one predator, representing by the following system of difference equations [10]:

    xn+1=(1+a)xnbxnzn, yn+1=ryn(1yn)cynzn, zn+1=(1d)zn+exnzn+fynzn, (1.2)

    where xn, yn and zn respectively denotes population densities of prey and predator respectively, and a,b,c,d,e,f and r are real positive parameters. More precisely, r is the net reproductive rate of the prey population, a denotes the intrinsic rate of growth the prey, b is the rate of conversion of a consumed prey to a predator, c is the per capita rate of predation of the predator, d is death rate of the predator, e is the rate of conversion of a consumed prey to a predator and f is the per-capita rate of predation of the predator.

    It is important here to mention that in 2013 Sagayaraj et al. [10] studied the local dynamics about three boundary fixed points only: P000(0,0,0), P0y0(0,r1r,0) and P0yz(0,df,rffdrcf), out of five fixed points including the unique positive fixed point. So for the interested readers, these results regarding the local dynamics about fixed points: P000(0,0,0), P0y0(0,r1r,0) and P0yz(0,df,rffdrcf) are summarized in Table 1.

    Table 1.  Summarized results of Sagayaraj et al. [10] about fixed points P000(0,0,0), P0y0(0,r1r,0) and P0yz(0,df,rffdrcf) of the model (1.2).
    Fixed pointsCorresponding behavior for the model (1.2)
    P000sink if r<1,2<a<0 and 0<d<2; source if a>0,r>1 and d>2;
    saddle if a>0,r<1 and 0<d<2; non-hyperbolic if r=1 or d=2.
    P0y0sink if 1<r<3, 2<a<0 and r<ffd; source if r>3,a>0 and r>ffd;
    saddle if 1<r<3,a>0 and r<ffd.
    P0yz sink if r>f(b+ac)b(fd) and r>ffd; source if r<f(b+ac)b(fd) and r<ffd;
    saddle if ffd<r<f(2c+ac+b)b(fd).

     | Show Table
    DownLoad: CSV

    From the Table 1, it is clear that for the said fixed points: P000(0,0,0), P0y0(0,r1r,0) and P0yz(0,df,rffdrcf), Sagayaraj et al. [10] cannot give complete local dynamical classifications to interested readers, which is still further consideration and improvements towards local dynamical properties along with topological classifications and bifurcation analysis for the model under consideration. So, in this paper we will study the local dynamics about fixed points, bifurcations like fold bifurcation, Neimark-Sacker or Hopf bifurcation, supercritical or subcritical Neimark-Sacker bifurcation, period-doubling or flip bifurcation, supercritical or subcritical period-doubling bifurcation, and hybrid control in the model (1.2). More precisely, our key contributions in this paper are as follows:

    ● To explore the local dynamical properties along with different topological classifications about fixed points of the under consideration model.

    ● To explore the existence of prime period and periodic points of the model (1.2).

    ● To explore the existence of possible bifurcations about fixed points of 3×3 model.

    ● Presentation of numerical simulation to verify theoretical results.

    ● To control the N-S and P-D bifurcations about the positive fixed point by Hybrid control strategy.

    This paper is structured as follows: In Section 2, we explored the existence of fixed points of the model (1.2) algebraically. The linearized form of the model (1.2) is presented in Section 3. In Section 4, we explored the local dynamics about boundary fixed points and unique positive fixed point of the model. Existence of prime period and periodic points of the model (1.2) are explored in Section 5. In Section 6, we explored the existence of bifurcations about boundary fixed points of the model (1.2) whereas detailed bifurcation analysis about boundary fixed points and unique positive fixed point of the model (1.2) are given in Section 7. Theoretical results are verified numerically in Section 8. In Section 9, the hybrid control strategy is presented to control the N-S and P-D bifurcations about the unique fixed point of the model while a brief summary of the paper is presented in Section 10.

    The existence of fixed points in the interior of R3+ is explored in this Section. The result regarding the existence of fixed points for the 3×3 discrete-time model (1.2) can be summarized as a following Lemma:

    Lemma 2.1. In the interior of R3+, model (1.2) has four boundary fixed points and a unique positive fixed point. Precisely

    (i)  a,b,c,d,e,f,r, P000(0,0,0) is a boundary fixed point of model (1.2);

    (ii) P0y0(0,r1r,0) is a boundary fixed point of the model (1.2) if r>1;

    (iii) P0yz(0,df,rffdrcf) is a boundary fixed point of the model (1.2) if f>drr1 where r>1;

    (iv) Px0z(de,0,ab) is also boundary fixed point of the model (1.2),  a,b,c,d>0;

    (v) P+xyz(br(df)+f(b+ac)ber,brbacbr,ab) is the unique positive fixed point of (1.2) if d>bf(r1)afcbr and b>acr1 where r>1.

    Proof. If Pxyz(x,y,z) is a fixed point of model (1.2) then

    x=(1+a)xbxz, y=ry(1y)cyz, z=(1d)z+exz+fyz. (2.1)

    Since for P000(0,0,0), P0y0(0,r1r,0), P0yz(0,df,rffdrcf) and Px0z(de,0,ab), Eq (2.1) satisfied identically, and hence P000(0,0,0), P0y0(0,r1r,0), P0yz(0,df,rffdrcf) and Px0z(de,0,ab) are boundary fixed points of the model (1.2). In order for the existence of the unique positive fixed point one need to solve the following equation simultaneously:

    abz=0, r(1y)cz=1, d+ex+fy=0. (2.2)

    From first equation of system (2.2) one gets:

    z=ab. (2.3)

    From second equation of system (2.2) and Eq (2.3), one gets:

    y=brbacbr. (2.4)

    Finally, from third equation of system (2.2) and Eq (2.4), one gets:

    x=br(df)+f(b+ac)ber. (2.5)

    In view of Eqs (2.3)–(2.5) one can conclude that Pxyz(br(df)+f(b+ac)ber,brbacbr,ab) is the interior fixed point of model (1.2). Additionally, P+xyz(br(df)+f(b+ac)ber,brbacbr,ab) is the unique positive fixed point of the model (1.2) if d>bf(r1)afcbr and b>acr1 where r>1.

    The linearized form of 3×3 model (1.2) about fixed point Pxyz(x,y,z) is explored in this Section. So we have the following map to create the corresponding linearized form of model (1.2):

    (f,g,h)(xn+1,yn+1,zn+1), (3.1)

    where

    f=(1+a)xbxz, g=ry(1y)cyz, h=(1d)z+exz+fyz. (3.2)

    Finally, J|Pxyz(x,y,z) about Pxyz(x,y,z) under the map, which is depicted in Eq (3.1), is

    J|Pxyz(x,y,z)=(1+abz0bx0r2ryczcyezfz1d+ex+fy). (3.3)

    Local dynamical properties along with topological classifications about fixed points: P000(0,0,0), P0y0(0,r1r,0), P0yz(0,df,rffdrcf), Px0z(de,0,ab) and P+xyz(br(df)+f(b+ac)ber,brbacbr,ab) of 3×3 model (1.2) are explored in this section.

    The J|P000(0,0,0) about P000(0,0,0) is

    J|P000(0,0,0)=(1+a000r0001d). (4.1)

    The eigenvalues of J|P000(0,0,0) about P000(0,0,0) are

    λ1=1+a, λ2=r, λ3=1d. (4.2)

    So by stability theory [11], the local dynamics of 3×3 model (1.2) about P000(0,0,0) can be conclude as follows.

    Lemma 4.1. For P000(0,0,0) following conditions hold:

    (i) P000(0,0,0) is never sink;

    (ii) P000(0,0,0) is a source if

    a>0, r>1 and d>2; (4.3)

    (iii) P000(0,0,0) is a saddle if

    a>0, 0<r<1 and 0<d<2; (4.4)

    (iv) P000(0,0,0) is non-hyperbolic if

    r=1, (4.5)

    or

    d=2. (4.6)

    The J|P0y0(0,r1r,0) about P0y0(0,r1r,0) is

    J|P0y0(0,r1r,0)=(1+a0002rc(1r)r001d+f(r1)r). (4.7)

    The eigenvalues of J|P0y0(0,r1r,0) about P0y0(0,r1r,0) are

    λ1=1+a, λ2=2r, λ3=1d+f(r1)r. (4.8)

    We also summarized the local dynamics of 3×3 model (1.2) about P0y0(0,r1r,0) as follows:

    Lemma 4.2. For P0y0(0,r1r,0), following statements hold:

    (i) P0y0(0,r1r,0) is never sink;

    (ii) P0y0(0,r1r,0) is a source if

    a>0, and 3<r<ffd+2; (4.9)

    (iii) P0y0(0,r1r,0) is a saddle if

    a>0, and ffd+2<r<3; (4.10)

    (iv) P0y0(0,r1r,0) is non-hyperbolic if

    r=3, (4.11)

    or

    r=ffd+2. (4.12)

    The J|P0yz(0,df,rffdrcf) about P0yz(0,df,rffdrcf) is

    J|P0yz(0,df,rffdrcf)=(1+a+bc(1r+drf)0001drfdcfec(r1drf)fc(r1drf)1). (4.13)

    The eigenvalues of J|P0yz(0,df,rffdrcf) about P0yz(0,df,rffdrcf) are

    λ1=1+a+bc(1r+drf),λ2,3=2fdrf±Δ2, (4.14)

    where

    Δ=(dr2ff)24(fdr+drffdd2rf),=4df2+d2r24drf2+4d2frf2. (4.15)

    In the following Lemmas, we give the topological classifications about P0yz(0,df,rffdrcf) if Δ0 (resp. Δ<0) by existing literature [9,12].

    Lemma 4.3. If Δ=4d2fr+d2r24drf2+4d2frf20 then for P0yz(0,df,rffdrcf) of model (1.2), following statements hold:

    (i) P0yz(0,df,rffdrcf) is a stable node if

    f(ac+b)b(fd)<r<f(2c+ac+b)b(fd)  and  ffd<r<f(4d)d(d+2f); (4.16)

    (ii) P0yz(0,df,rffdrcf) is an unstable node if

    r>max{f(4d)d(d+2f),f(2c+ac+b)b(fd)}; (4.17)

    (iii) P0yz(0,df,rffdrcf) is saddle node if

    f(ac+b)b(fd)<r<f(2c+ac+b)b(fd)  and  r>f(4d)d(d+2f), (4.18)

    or

    r>f(2c+ac+b)b(fd)  and  ffd<r<f(4d)d(d+2f), (4.19)

    (iv) P0yz(0,df,rffdrcf) is non-hyperbolic if

    r=f(4d)d(d+2f), (4.20)

    or

    r=f(2c+ac+b)b(fd). (4.21)

    Lemma 4.4. If Δ=4d2fr+d2r24drf2+4d2frf2<0 then for P0yz(0,df,rffdrcf) of model (1.2), following statements hold:

    (i) P0yz(0,df,rffdrcf) is stable focus-node if

    f(ac+b)b(fd)<r<f(2c+ac+b)b(fd)  and  r<ff1d; (4.22)

    (ii) P0yz(0,df,rffdrcf) is an unstable focus-node if

    r>max{f(2c+ac+b)b(fd),ff1d}; (4.23)

    (iii) P0yz(0,df,rffdrcf) is saddle focus-node if

    f(ac+b)b(fd)<r<f(2c+ac+b)b(fd)  and  r>ff1d, (4.24)

    or

    r>f(2c+ac+b)b(fd)  and  r<ff1d; (4.25)

    (iv) P0yz(0,df,rffdrcf) is non-hyperbolic if

    r=ff1d, (4.26)

    or

    r=f(2c+ac+b)b(fd). (4.27)

    The J|Px0z(de,0,ab) about Px0z(de,0,ab) is

    J|Px0z(de,0,ab)=(10bde0racb0acbafb1). (4.28)

    The eigenvalues of J|Px0z(de,0,ab) about Px0z(de,0,ab) are

    λ1=racb,λ2,3=1±ιab. (4.29)

    Lemma 4.5. For Px0z(de,0,ab), following statements hold:

    (i) Px0z(de,0,ab) is never stable focus-node;

    (ii) Px0z(de,0,ab) is an unstable focus-node if

    r<acb1; (4.30)

    (iii) Px0z(de,0,ab) is a saddle focus-node if

    r>acb1; (4.31)

    (iv) Px0z(de,0,ab) is non-hyperbolic if

    r=acb1. (4.32)

    The J|P+xyz(br(df)+f(b+ac)ber,brbacbr,ab) about P+xyz(br(df)+f(b+ac)ber,brbacbr,ab) is

    J|P+xyz(br(df)+f(b+ac)ber,brbacbr,ab)=(10br(df)+f(b+ac)re02r+acbc(b+ac)brcaebafb1). (4.33)

    The characteristic equation of J|P+xyz(br(df)+f(b+ac)ber,brbacbr,ab) about P+xyz(br(df)+f(b+ac)ber,brbacbr,ab) is

    P(λ)=λ3+δ1λ2+δ2λ+δ3, (4.34)

    where

    δ1=r4acb,δ2=52r+adaf+2acbacfbra2c2fb2r+acfb+afr+a2cfbr,δ3=r2+adrafr2ad+3af+acb+acfbr+a2c2fb2racfb3a2cfbr+2a2cfba2cdb2afra3c2fb2r. (4.35)

    Now in the following lemma, we explore the necessary and sufficient condition under which the roots of (4.34) lies on open unit disk by utilizing Theorem 1.4 of [11].

    Lemma 4.6. P+xyz(br(df)+f(b+ac)ber,brbacbr,ab) of 3×3 model (1.2) is a sink iff the following conditions hold:

    |δ1+δ3|<1+δ2,  |δ13δ3|<3δ2,  δ23+δ2δ3δ1<1, (4.36)

    where δ1,δ2 and δ3 are depicted in (4.35).

    We explore that the fixed points: P000(0,0,0), P0y0(0,r1r,0), P0yz(0,df,rffdrcf), Px0z(de,0,ab) and P+xyz(br(df)+f(b+ac)ber,brbacbr,ab) of the model (1.2) are periodic points of prime period- 1. Additionally, we also explore that P000(0,0,0), P0y0(0,r1r,0), P0yz(0,df,rffdrcf), Px0z(de,0,ab) and P+xyz(br(df)+f(b+ac)ber,brbacbr,ab) of the model (1.2) are periodic points having period- 2,3,,n, by existing Literature [13,14,15].

    Theorem 5.1. Fixed points: P000(0,0,0), P0y0(0,r1r,0), P0yz(0,df,rffdrcf), Px0z(de,0,ab) and P+xyz(br(df)+f(b+ac)ber,brbacbr,ab) of model (1.2) are periodic points of prime period- 1.

    Proof. From model (1.2), one has

    F(x,y,z):=(f(x,y,z),g(x,y,z),h(x,y,z)), (5.1)

    where f,g and h are represented in system (3.2). After straightforward computation, from Eq (5.1) one gets:

    F|P000(0,0,0)=P000(0,0,0),F|P0y0(0,r1r,0)=P0y0(0,r1r,0),F|P0yz(0,df,rffdrcf)=P0yz(0,df,rffdrcf),F|Px0z(de,0,ab)=Px0z(de,0,ab),F|P+xyz(br(df)+f(b+ac)ber,brbacbr,ab)=P+xyz(br(df)+f(b+ac)ber,brbacbr,ab). (5.2)

    From Eq (5.2), one can conclude that P000(0,0,0), P0y0(0,r1r,0), P0yz(0,df,rffdrcf), Px0z(de,0,ab) and P+xyz(br(df)+f(b+ac)ber,brbacbr,ab) of model (1.2) is a periodic point of prime period- 1.

    Hereafter we proved that fixed points: P000(0,0,0), P0y0(0,r1r,0), P0yz(0,df,rffdrcf), Px0z(de,0,ab) and P+xyz(br(df)+f(b+ac)ber,brbacbr,ab) of 3×3 model (1.2) are periodic points having period- 2,3,,n.

    Theorem 5.2. The fixed point: P000(0,0,0) of 3×3 model (1.2) is a periodic point of period- 2,3,,n.

    Proof. From Eq (5.1), one has

    F2(x,y,z)=((1+a)f(x,y,z)bf(x,y,z)h(x,y,z), rg(x,y,z)(1g(z,y,z))cg(x,y,z)h(x,y,z), (1d)h(x,y,z)+ef(x,y,z)h(x,y,z)+ fg(x,y,z)h(x,y,z))F2|P000(0,0,0)=P000(0,0,0),F3(x,y,z)=((1+a)f2(x,y,z)bf2(x,y,z)h2(x,y,z), rg2(x,y,z)(1g2(z,y,z))cg2(x,y,z)h2(x,y,z),(1d)h2(x,y,z)+ef2(x,y,z)h2(x,y,z)+ fg2(x,y,z)h2(x,y,z))F2|P000(0,0,0)=P000(0,0,0), Fi(x,y,z)=((1+a)fi1(x,y,z)bfi1(x,y,z)hi1(x,y,z), rgi1(x,y,z)(1gi1(z,y,z))cgi1(x,y,z)hi1(x,y,z),(1d)hi1(x,y,z)+efi1(x,y,z)hi1(x,y,z)+ fgi1(x,y,z)hi1(x,y,z))Fi|P000(0,0,0)=P000(0,0,0). (5.3)

    Thus Eq (5.3) indicates that P000(0,0,0) of model (1.2) is a periodic point of period-2,3,,n.

    Theorem 5.3. The fixed points: P0y0(0,r1r,0), P0yz(0,df,rffdrcf), Px0z(de,0,ab) and P+xyz(br(df)+f(b+ac)ber,brbacbr,ab) of 3×3 model (1.2) are periodic points of period- 2,3,,n.

    Proof. Same as the proof of Theorem 5.2.

    Existence of bifurcations about boundary fixed points: P000(0,0,0), P0y0(0,r1r,0), P0yz(0,df,rffdrcf) and Px0z(de,0,ab) are explored in this Section. Moreover occurrence of bifurcations about boundary fixed points can be summarized as follows:

    (i) About P000(0,0,0), one of the eigenvalues λ2|(4.5)=1 but λ1,3=1+a,1d1 or 1. So 3×3 model (1.2) may undergoes a fold bifurcation if parameters Λ=(a,b,c,d,e,f,r) goes through the following curve:

    FB|P000(0,0,0)={Λ:r=1}. (6.1)

    For the interested readers, it is important here to mention that we have introduced new notation like λ2|(4.5) which stands that we have evaluated the eigenvalue λ2 at non-hyperbolic condition, which is depicted in Eq (4.5). Similarly we use such notation for same purpose in the rest of the paper.

    Moreover, about P000(0,0,0), one of the eigenvalues λ3|(4.6)=1 but λ1,2=1+a,r1 or 1. So 3×3 model (1.2) may undergoes a P-D bifurcation if parameters Λ goes through the following curve:

    PDB|P000(0,0,0)={Λ:d=2}. (6.2)

    It is noted that for d=2 the predator leaves a negative number of offspring, i.e., in the absence of prey and therefore its no make sense biology. Thus we will exclude this case from further study.

    (ii) About P0y0(0,r1r,0), one of the eigenvalues λ2|(4.11)=1 but λ1,3=1+a,1d+f(r1)r1 or 1. So 3×3 model (1.2) undergoes a P-D bifurcation if parameters Λ goes through the following curve:

    PDB|P0y0(0,r1r,0)={Λ:r=3}. (6.3)

    Also P0y0(0,r1r,0), one of the eigenvalues λ3|(4.12)=1 but λ1,2=1+a,2r1 or 1. So 3×3 model (1.2) undergoes a P-D bifurcation if Λ goes through the following curve:

    PDB|P0y0(0,r1r,0)={Λ:r=ffd+2}. (6.4)

    (iii) About P0yz(0,df,rffdrcf), pair of complex conjugate eigenvalues |λ2,3|(4.26)=1 but λ1=1+a+bc(1r+drf)1 or 1. So 3×3 model (1.2) undergoes a N-S bifurcation if parameters Λ goes through the following curve:

    NSB|P0yz(0,df,rffdrcf)={4d2fr+d2r24drf2+4d2frf2<0, r=ff1d}. (6.5)

    (iv) About Px0z(de,0,ab), λ1|(4.32)=1 but |λ2,3|=1+ad>1. So 3×3 model (1.2) may undergoes P-D bifurcation if Λ goes through the following curve:

    PDB|Px0z(de,0,ab)={Λ:r=acb1}. (6.6)

    For the completeness of this Section, the detailed bifurcation analysis about fixed points of the model (1.2) is presented as follows:

    About P000(0,0,0), 3×3 model (1.2) undergoes a fold bifurcation when parameters Λ goes through the curve, which is depicted in Eq (6.1). But by computation, fold bifurcation do not occur and therefore P000(0,0,0) is degenerate with high co-dimension for the original parameters Λ of 3×3 model (1.2).

    About P0y0(0,r1r,0), 3×3 model (1.2) undergoes P-D bifurcation when Λ goes through the curve, which is depicted in Eq (6.3). Since model (1.2) is invariant w.r.t x=z=0 then we can restrict the model (1.2) on x=z=0, to determine the bifurcation, where model (1.2) takes the form:

    yn+1=ryn(1yn). (7.1)

    Also the map corresponds to Eq (7.1) takes the form:

    f(y):=ryry2. (7.2)

    Denote r=r=3, y=y=r1r. By computation one gets: fy|r=r=3, y=y=r1r=1 which is non-hyperbolic condition for P-D bifurcation. Moreover, the non-degenerate conditions: fyy|r=r=3, y=y=r1r=60 and fr|r=r=3, y=y=r1r=290 hold, which implies that 3×3 model (1.2) undergoes the P-D bifurcation when Λ goes through the curve, which is depicted in Eq (6.3). Additionally, (fyyfr+2fyr)(2fyyy+3f2yy)|r=r=3, y=y=r1r=81<0, and hence by existing literature [15,16,17] model (1.2) undergoes the supercritical P-D bifurcation.

    In this Subsection, N-S bifurcation of the 3×3 model (1.2) is explored about P0yz(0,df,rffdrcf), if Λ goes through the curve, which is depicted in Eq (6.5) by assuming that rf(2c+ac+b)b(fd). Let

    un=xn, vn=yndf, wn=znrffdrcf, (7.3)

    then 3×3 model (1.2) takes the form:

    (un+1vn+1vn+1)=(1+a+bc(1r+drf)0001drfdcfec(r1drf)fc(r1drf)1)(unvnvn)+(bununrv2ncvnwneunwn+fvnwn), (7.4)

    whose linear part is same as J|P0yz(0,df,rffdrcf), which is depicted in Eq (4.13). Moreover at Eq (4.26), following non-degenerate conditions hold:

    |λ2,3|(4.26)=1 and d|λ2,3|dr|(4.26)=fddd22f0, (7.5)

    and additionally it is require that λm2,31 where m=1,2,3,4 which also true by computation. Now using the transformation:

    (unvnwn)=T(xtytzt), (7.6)

    Equation (7.4) takes the following form:

    (xt+1yt+1zt+1)=(2fdr2f4drf24df24d2frr2d22f04drf24df24d2frr2d22f2fdr2f0001+a+bc(1r+drf))×(xtytzt)+(Γ(xt,yt,zt)Φ(xt,yt,zt)Ψ(xt,yt,zt)), (7.7)

    where

    T=(00c2dfe(cdr+acfbfr+bf+bdr)bf+acf+bdrbfre(f+drfr)cdr2f(frfdr)c4drf24df24d2frd2r22f(frfdr)c2dbf+acf+bdr+cdrbfr101), (7.8)

    and

    Γ(xt,yt,zt)=bxtzt+bz2t+(xtzt+z2t)(c2dfbf+acf+bdr+cdrbfrbf+acf+bdrbfrf+drfr)cdr2(frfdr)x2t+c4drf24df24d2frr2d22(frfdr)xtytc2fdbf+acf+bdr+cdrbfrztxtcdr2(frfdr)xtzt+c4drf24df24d2frr2d22(frfdr)ytztc2fdbf+acf+bdr+cdrbfrz2t,Φ(xt,yt,zt)=bcdr2bc2df(frfdr)bf+acf+bdr+cdrbfrc4drf24df24d2frr2d2(xtzt+z2t)+2f(frfdr)c4drf24df24d2frr2d2×[r(cdr2f(frfdr)xt+c4drf24df24d2frr2d22f(frfdr)ytc2dbf+acf+bdr+cdrbfrzt)2+(xt+zt)(c2dr2f(frfdr)xtc24drf24df24d2frr2d22f(frfdr)yt+c3dbf+acf+bdr+cdrbfrzt)]+dr4drf24df24d2frr2d2×[zt(xt+zt)(c2dfbf+acf+bdr+cdrbfrbf+acf+bdrbfr(f+drfr))+(xt+zt)(cdr2(frfdr)xt+c4drf24df24d2frr2d22(frfdr)ytc2dfbf+acf+bdr+cdrbfrzt)],Ψ(xt,yt,zt)=bxtztbz2t. (7.9)

    Hereafter, we have the system (7.7) on the center manifold, i.e.,

    Wc(0)={(xt,yt,zt)zt=h(xt,yt),h(0,0)=0,Dh(0,0)=0}, (7.10)

    where

    h(xt,yt)=α1x2t+α2xtyt+α3y2t+O((|xt|+|yt|)3). (7.11)

    In view of Eqs (7.10) and (7.7) becomes:

    h(2fdr2fxt+4drf24df24d2frr2d22fyt+Γ(xt,yt,h),4drf24df24d2frr2d22fxt+2fdr2fyt+Φ(xt,yt,h))=(1+abc(r1drf))h(xt,yt)+Ψ(xt,yt,h). (7.12)

    After some manipulation, from Eq (7.12) one gets: α1=α2=α3=0. So, the system (7.7) restricted to Wc(0) becomes

    (xt+1yt+1)=(2fdr2f4drf24df24d2frr2d22f4drf24df24d2frr2d22f2fdr2f)(xtyt)+(ζ1(xt,yt)ζ2(xt,yt)), (7.13)

    where

    ζ1(xt,yt)=cdr2(frfdr)x2t+c4drf24df24d2frr2d22(frfdr)xtyt,ζ2(xt,yt)=2f(frfdr)c4drf24df24d2frr2d2[r(cdr2f(frfdr)xt+c4drf24df24d2frr2d22f(frfdr)yt)2+(c2dr2f(frfdr))x2tc24drf24df24d2frr2d22f(frfdr)xtyt]+dr4drf24df24d2frr2d2(cdr2(frfdr)x2t+c4drf24df24d2frr2d22(frfdr)xtyt). (7.14)

    In order for Eq (7.13) to undergoes a N-S bifurcation, one require that the following discriminatory quantity σ should be non-zero [12,18,19]:

    σ=Re((12λ)ˉλ21λL11L12)12|L11|2|L21|2+Re(ˉλL22), (7.15)

    where

    L11=14[ˉζ1xtxt+ˉζ1ytyt+ι(ˉζ2xtxt+ˉζ2ytyt)],L12=18[ˉζ1xtxtˉζ1ytyt+2ˉζ2xtyt+ι(ˉζ2xtxtˉζ2ytyt2ˉζ1xtyt)],L21=18[ˉζ1xtxtˉζ1ytyt2ˉζ2xtyt+ι(ˉζ2xtxtˉζ2ytyt+2ˉζ1xtyt)],L22=116[ˉζ1xtxtxt+ˉζ1xtytyt+ˉζ2xtxtyt+ˉζ2ytytyt+ι(ˉζ2xtxtxt+ˉζ2xtytytˉζ1xtxtytˉζ1ytytyt)]. (7.16)

    After straightforward calculation, one gets:

    L11=14[cdrfrfdr+ι(14drf24df24d2frr2d2(2c2drc2d2r3f(frfdr))cd2r2(frfdr)4drf24df24d2frr2d2rc4drf24df24d2frr2d2f(frfdr))],L12=18[cdrf+2cdr22cf2r+2cf2+cdrf(frfdr)+ι(1c4drf24df24d2frr2d2×(2c2drc2d2r3f(frfdr))cd2r2(frfdr)4drf24df24d2frr2d2+c4drf24df24d2frr2d2f(frfdr)(rf))],L21=18[2cf2r3cdrf2cdr22cf2cdrf(frfdr)+ι(1c4drf24df24d2frr2d2(2c2drc2d2r3f(frfdr))cd2r2(frfdr)4drf24df24d2frr2d2+c4drf24df24d2frr2d2f(frfdr)(r+f))],L22=0. (7.17)

    Based on N-S bifurcation analysis about P0yz(0,df,rffdrcf), one has

    Theorem 7.1. 3×3 model (1.2) undergoes a N-S bifurcation at P0yz(0,df,rffdrcf) if σ0. Moreover, if σ<0 (resp. σ>0), then an attracting (resp. repelling) closed curve bifurcates from P0yz(0,df,rffdrcf). Finally Supercritical (resp. subcritical) N-S bifurcation take place if σ<0 (resp. σ>0).

    Remark 1: Recall that about Px0z(de,0,ab), 3×3 model (1.2) may undergoes P-D bifurcation if Λ goes through the curve, which is depicted in Eq (6.6). But by computation, P-D bifurcation do not occur and therefore Px0z(de,0,ab) is degenerate with high co-dimension.

    Hereafter, we explore the N-S as well as P-D bifurcations of the 3×3 model (1.2) without calculating the eigenvalues of J|P+xyz(br(df)+f(b+ac)ber,brbacbr,ab) about P+xyz(br(df)+f(b+ac)ber,brbacbr,ab). Precisely, in the following two theorems, we explore N-S and P-D bifurcations by utilizing explicit criteria, i.e., the criteria which is formulated by using a set of inequalities or equalities that contain the coefficients of the characteristic equation derived from J|P+xyz(br(df)+f(b+ac)ber,brbacbr,ab) about the unique positive fixed point [20,21,22].

    Theorem 7.2. If

    1δ2+δ3(δ1δ3)=0,1+δ2δ3(δ1+δ3)>0,1+δ1+δ2+δ3>0,1δ1+δ2δ3>0,ddr(1δ2+δ3(δ1δ3))|r=r00,cos2πl11+δ1+δ2+δ32(1+δ2), l=3,4,, (7.18)

    then 3×3 model (1.2) undergoes a N-S bifurcation about P+xyz(br(df)+f(b+ac)ber,brbacbr,ab), where δ1,δ2 and δ3 are depicted in Eq (4.35).

    Proof. By utilizing explicit criterion, which is depicted in [20,21], for n=3 one gets:

    Δ2(r)=1δ2+δ3(δ1δ3)=0,Δ+2(r)=1+δ2δ3(δ1+δ3)>0,Pr(1)=1+δ1+δ2+δ3>0,(1)3Pr(1)=1δ1+δ2δ3>0,ddr(Δ2(r))|r=r0=ddr(1δ2+δ3(δ1δ3))|r=r00. (7.19)

    Finally

    10.5Pr(1)Δ0(r)Δ+1(r)=11+δ1+δ2+δ32(1+δ2).

    Theorem 7.3. If

    1δ2+δ3(δ1δ3)>0,1+δ2δ3(δ1+δ3)>0,1±δ2>0,1+δ1+δ2+δ3>0,1+δ1δ2+δ3=0,δ1δ2+δ332δ1+δ20, (7.20)

    then 3×3 model (1.2) undergoes a P-D bifurcation about P+xyz(br(df)+f(b+ac)ber,brbacbr,ab).

    Proof. By utilizing explicit criterion, which is depicted in [21,22], for n=3 one gets:

    Δ2(r)=1δ2+δ3(δ1δ3)>0,Δ+2(r)=1+δ2δ3(δ1+δ3)>0,Δ±1(r)=1±δ2>0,Pr(1)=1+δ1+δ2+δ3>0,Pr(1)=1+δ1δ2+δ3=0. (7.21)

    Finally

    ni=1(1)niδini=1(1)ni(ni+1)δi1=δ1δ2+δ332δ1+δ20.

    Theoretical results are illustrated numerically in this Section. In this regard, following cases are presented to discuss the correctness of obtained theoretical results about fixed points for the 3×3 model (1.2):

    Case I: If a=0.007,b=0.001,c=0.1,d=0.39,e=0.5,f=2.7 and varying r[1.2,10.4] then Figure 1a clearly indicates that about P0y0(0,0.6884735202492211,0), model (1.2) undergoes P-D bifurcation, which correctness the theoretical discussion in Subsection 7.2. More 3D bifurcation diagrams are plotted and depicted in Figure 1b, 1c.

    Figure 1.  P-D bifurcation diagram about P0y0(0,r1r,0) of model (1.2) by varying r[1.2,10.4].

    Case II: If a=0.45, b=0.9,c=0.99,d=0.29,e=0.3,f=2.3 then from Eq (4.26) one gets: r=2.2772277227722775. By Lemma 4.4, P0yz(0,df,rffdrcf) is locally asymptotically focus if r<2.2772277227722775, exchange of stability if r=2.2772277227722775, and unstable focus if r>2.2772277227722775. For instance, if r=1.4<2.2772277227722775 then it is clear from Figure 2a that P0yz(0,0.12608695652173912,0.2257356170399649) of model (1.2) is locally asymptotically focus. Similarly others values of r, if r<2.2772277227722775 then P0yz(0,df,rffdrcf) of model (1.2) is also locally asymptotically focus (Figure 2b, 2c). But if r>2.2772277227722775 then discrete-time model (1.2) loss its stability and meanwhile stable closed curve appears (Figure 3a, 3c). The appearance of these closed curve implies that model (1.2) undergoes supercritical N-S bifurcation when parameters goes through the curve, which is depicted in Eq (6.5). This correctness of theoretical results obtain in Subsection 7.3. Finally, the bifurcation diagrams along with Maximum Lyapunov Exponent in this case, are plotted and drawn in Figure 4.

    Figure 2.  Topological classification of model (1.2) about P0yz(0,df,rffdrcf) if r<2.2772277227722775.
    Figure 3.  Supercritical N-S bifurcation of model (1.2) about P0yz(0,df,rffdrcf) if r>2.2772277227722775.
    Figure 4.  Bifurcation diagram and their corresponding Maximum Lyapunov Exponent of the model (1.2) about P0yz(0,df,rffdrcf). (4a, 4b) Bifurcation diagram of the model if r[1.2,10.4] and initial condition (0,0.17,0.54). (4c) Maximum Lyapunov Exponent corresponding to Figure 4a, 4b.

    Case III: If a=0.45, b=0.93,c=0.8,d=0.95,e=0.15,f=0.45 and parameter r vary from 1.23 to 10.4, i.e., r[1.23,10.4] then at r=2.86317 the P+xyz(br(df)+f(b+ac)ber,brbacbr,ab) loss its stability. That means that P+xyz(br(df)+f(b+ac)ber,brbacbr,ab) is locally stable if r<2.86317. For instance, if r=1.75<2.86317 then it is clear from Figure 5a that unique positive fixed point P+xyz(19.044546850998465,0.2737327188940094,0.48387096774193544) of model (1.2) is locally stable. Similarly others values of r, if r<2.86317 then P+xyz(br(df)+f(b+ac)ber,brbacbr,ab) of model (1.2) is also locally stable (Figure 5b, 5c). Also if a=0.45,b=0.93,c=0.8,d=0.95,e=0.15,f=0.45,r=2.86317 then Eq (4.34) along with Eq (4.35) takes the following form:

    P(λ)=λ31.5239267741935483λ2+0.15468797005336343λ+0.4321468661912955=0. (8.1)
    Figure 5.  Topological classification of model (1.2) about P+xyz(br(df)+f(b+ac)ber,brbacbr,ab) if r<2.86317.

    By computer computation, from Eq (8.1) one gets:

    λ1=0.43214706529074226, λ2,3=0.9780369197421452±0.208430618911403ι

    with

    |λ2,3|=|0.9780369197421452±0.208430618911403ι|=1.

    This grantee that for above parametric values, discrete-time model (1.2) undergoes N-S bifurcation when parameter r goes through the bifurcation values, i.e., r=2.86317. Moreover the parametric conditions, which are depicted in Eq (7.18), under which the model (1.2) undergoes a N-S bifurcation hold, i.e.,

    1δ2+δ3(δ1δ3)=0,1+δ2δ3(δ1+δ3)=1.6264972358671579>0,1+δ1+δ2+δ3=0.06290806205111066>0,1δ1+δ2δ3=2.246467878055616>0,ddr(1δ2+δ3(δ1δ3))|r=2.86317=1.94299184399591770,cos2πl11+δ1+δ2+δ32(1+δ2)=0.9727597135838336,  l=3,4,5,,

    which again grantee that model (1.2) undergoes a N-S bifurcation if parameters goes through the bifurcation value. Hence if r>2.86317 then discrete-time model (1.2) loss its stability and meanwhile stable closed curve appears (Figure 6a6c). The appearance of these closed curve implies that model (1.2) undergoes supercritical N-S bifurcation when parameters goes through r>2.86317. This correctness of theoretical results obtain in Subsection 7.4. Finally, the bifurcation diagrams along with Maximum Lyapunov Exponent in this case, are plotted and drawn in Figure 7.

    Figure 6.  Supercritical N-S bifurcation of model (1.2) about P+xyz(br(df)+f(b+ac)ber,brbacbr,ab) if r>2.86317.
    Figure 7.  Bifurcation diagram and their corresponding Maximum Lyapunov Exponent of the model (1.2) about P+xyz(br(df)+f(b+ac)ber,brbacbr,ab). (7a, 7b) Bifurcation diagram of the model if r[1.1,10.4] and initial condition (0.2,0.2,0.2). (7c) Maximum Lyapunov Exponent corresponding to Figure 7a, 7b.

    Case IV: If a=0.3, b=1.99,c=3.75,d=0.75,e=0.235,f=1.0305 and parameter r vary from 1.6 to 5.2, i.e., r[1.6,5.2] then at r=3.9603484492658967 the P+xyz(br(df)+f(b+ac)ber,brbacbr,ab) loss its stability. That means that P+xyz(br(df)+f(b+ac)ber,brbacbr,ab) is stable if r<3.9603484492658967, unstable if r=3.9603484492658967 and P-D bifurcation appears if r>3.9603484492658967. For instance, if a=0.3,b=1.99,c=3.75,d=0.75,e=0.235,f=1.0305,r=3.97 then Eq (4.34) along with Eq (4.35) takes the following form:

    P(λ)=λ30.6049781838999324λ20.8791823027251384λ+0.7257958811747927=0. (8.2)

    From Eq (8.2) one gets:

    λ1=1, λ2,3=0.802489091949966±0.28601947219745705ι,

    with

    |λ2,3|=|0.802489091949966±0.28601947219745705ι|=0.85193654762241121.

    This grantee that for above parametric values, discrete-time model (1.2) undergoes P-D bifurcation when parameter r goes through the bifurcation values, i.e., r=3.9603484492658967. Moreover the parametric conditions, which are depicted in Eq (7.20), under which the model (1.2) undergoes a P-D bifurcation hold, i.e.,

    1δ2+δ3(δ1δ3)=0.9133119675196675>0,1+δ2δ3(δ1+δ3)=0.03312871021974505>0,1+δ2=0.12081769727486158>0,1δ2=1.8791823027251384>0,1+δ1+δ2+δ3=0.24163539454972183>0,1+δ1δ2+δ3=0,δ1δ2+δ332δ1+δ2=0.078240040995743830,

    which again grantee that model (1.2) undergoes a P-D bifurcation. The bifurcation diagrams along with Maximum Lyapunov Exponent in this case, are plotted and drawn in Figure 8. More figure are plotted and drawn in Figure 9 that indicates that model (1.2) complex dynamics of period-2, 10, 11, 12 respectively. For the qualitative behavior of the continuous-time model, we refer the reader [23,24,25] and references cited therein.

    Figure 8.  Bifurcation diagram and their corresponding Maximum Lyapunov Exponent of the model (1.2) about P+xyz(br(df)+f(b+ac)ber,brbacbr,ab). (8a, 8b) Bifurcation diagram of the model if r[1.5,5.2] and initial condition (0.5,0.15,0.15). (8c) Maximum Lyapunov Exponent corresponding to Figure 8a, 8b.
    Figure 9.  Complex dynamics of model (1.2) in the Case Ⅳ.

    In discrete-time models, chaos control can be obtained using various methods, some of them are the state feedback method, pole-placement technique and hybrid control method. Among all these methods, hybrid control technique is the simplest to implement. Moreover this method can be implemented for controlling chaos under the influence of both N-S and P-D bifurcations. To control the chaos in model (1.2), we will utilize Hybrid control strategy to control the N-S and P-D bifurcations about the positive fixed point: P+xyz(br(df)+f(b+ac)ber,brbacbr,ab). For more detailed study about the Hybrid control strategy, we refer the reader to [26,27,28] and reference cited therein.

    The Hybrid control strategy is applied for control of N-S bifurcation in this Subsection. Recall that if the uncontrolled discrete-time model (1.2) undergoes a N-S bifurcation about P+xyz(br(df)+f(b+ac)ber,brbacbr,ab) then the controlled discrete-time model correspond to model (1.2) becomes:

    xn+1=α((1+a)xnbxnzn)+(1α)xn,yn+1=α(ryn(1yn)cynzn)+(1α)yn,zn+1=α((1d)zn+exnzn+fynzn)+(1α)zn, (9.1)

    where 0<α<1. Here it is noted that control model (9.1) reduces to model (1.2) if α=1, and chaos for model (9.1) can be delayed, advance or even completely eliminate by choosing the suitable choice of the parametric value of α. Now if one reconsider the Case Ⅲ of Section 8 then J|P+xyz(19.044546850998465,0.20737327188940094,0.48387096774193544) about P+xyz(19.044546850998465,0.20737327188940094,0.48387096774193544) for the controlled system (9.1) becomes

    J|P+xyz(19.044546850998465,0.20737327188940094,0.48387096774193544)=(106.63919α011.47607α0.41243α0.0725806α1.42742α1). (9.2)

    The characteristic polynomial of Eq (9.2) becomes

    F(λ)=λ3+(1.47607α3)λ2+(0.106834α22.95215α+3)λ0.711285α30.106834α2+1.47607α1=0, (9.3)

    whose roots lies in an open disk |λ|<1 iff α(0,0.66675). For instance, if α=0.56 the phase portraits for controlled system (9.1) are plotted and drawn in Figure 10.

    Figure 10.  Phase portraits for controlled system (9.1).

    Recall that if Case Ⅳ of Section 8 hold then J|P+xyz(13.30555251494148,0.604750275583457,0.1507537688442211) about P+xyz(13.30555251494148,0.604750275583457,0.1507537688442211) for the controlled system, which is depiction in system (9.1), becomes

    J|P+xyz(13.30555251494148,0.604750275583457,0.1507537688442211)=(1010.4777α012.39502α2.26781α0.0354271α0.565327α1). (9.4)

    The characteristic polynomial of Eq (9.4) becomes

    F(λ)=λ3+(2.39502α3)λ2+(0.910861α24.79004α+3)λ0.889018α30.9108661α2+2.39502α1=0, (9.5)

    whose roots lies in an open disk |λ|<1 iff α(0,0.543543). For instance, if α=0.50 the phase portraits for controlled system (9.1) are plotted and drawn in Figure 11.

    Figure 11.  Phase portraits for controlled system (9.1).

    This work deals with the study of local dynamics along with different topological classifications, bifurcations and Hybrid control of the model (1.2) in R3+. In the interior of R3+, it is explored that the model has boundary fixed points: P000(0,0,0), P0y0(0,r1r,0), P0yz(0,df,rffdrcf), Px0z(de,0,ab) and a unique positive fixed point: P+xyz(br(df)+f(b+ac)ber,brbacbr,ab) under certain parametric conditions. We have explored the local dynamics along with different topological classifications about fixed points: P000(0,0,0), P0y0(0,r1r,0), P0yz(0,df,rffdrcf) and Px0z(de,0,ab) by applying the method of Linearization and stability theory, and conclusions are presented in Table 2. We have also investigated the local dynamics about positive fixed point: P+xyz(br(df)+f(b+ac)ber,brbacbr,ab), and proved that it is a sink if and only if Eq (4.36) along with Eq (4.35) hold. We have also investigated the existence of prime period and periodic points for the model (1.2). Further, about each fixed point of the under consideration model, we have discussed the possible bifurcation scenarios. In general, bifurcation theory is the mathematical study of changes in the qualitative or topological structure of a given family, such as the integral curves of a family of vector fields and the solutions of family of differential equations. It is most commonly applied to the mathematical study of dynamical systems, a bifurcation occurs when a small smooth change made to the parameter values namely the bifurcation parameters of the system causes a sudden qualitative or topological change in its behavior. So, in the current study, we proved that about P000(0,0,0) there exist a fold bifurcation as well as P-D bifurcation if Λ respectively goes through the curves: FB|P000(0,0,0)={Λ:r=1} and PDB|P000(0,0,0)={Λ:d=2}. But by computation it goes to degenerate to high dimension and hence there does not such bifurcations about P000(0,0,0). About P0y0(0,r1r,0), it is proved that there exist a P-D bifurcation when Λ goes through the curve: PDB|P0y0(0,r1r,0)={Λ:r=3}. We have also proved that about P0yz(0,df,rffdrcf) the model (1.2) undergoes a N-S bifurcation when Λ goes through the curve: NSB|P0yz(0,df,rffdrcf)={4d2fr+d2r24drf2+4d2frf2<0, r=ff1d}. By explicit criterion(without finding eigenvalues) it is proved that about P+xyz(br(df)+f(b+ac)ber,brbacbr,ab) the model (1.2) undergoes N-S as well as P-D bifurcations. Numerical simulations are presented to verify theoretical results. Further, we have computed maximum Lyapunov exponents numerically. At the end, hybrid control strategy is applied to stabilize chaos existing in the model (1.2). Global dynamics about fixed points: P000(0,0,0), P0y0(0,r1r,0), P0yz(0,df,rffdrcf), Px0z(de,0,ab) and P+xyz(br(df)+f(b+ac)ber,brbacbr,ab), existence of global bifurcations and construction of forbidden set for the model (1.2) are our next aim to study.

    Table 2.  Fixed points with corresponding behavior of the 3×3 model (1.2).
    Fixed pointsCorresponding behavior
    P000never sink; source if r>1 and  d>2; saddle if 0<r<1 and 0<d<2;
    non-hyperbolic if r=1  or  d=2.
    P0y0never sink; source if 3<r<ffd+2; saddle if ffd+2<r<3;
    non-hyperbolic if r=3 or r=ffd+2.
    P0yzstable node if
    f(ac+b)b(fd)<r<f(2c+ac+b)b(fd)  and  ffd<r<f(4d)d(d+2f);
    unstable node if
    r>max{f(4d)d(d+2f),f(2c+ac+b)b(fd)};
    saddle node if
    f(ac+b)b(fd)<r<f(2c+ac+b)b(fd)  and  r>f(4d)d(d+2f);
    or
    r>f(2c+ac+b)b(fd)  and  ffd<r<f(4d)d(d+2f);
    non-hyperbolic if
    r=f(4d)d(d+2f) or  r=f(2c+ac+b)b(fd);
    stable focus-node if
    f(ac+b)b(fd)<r<f(2c+ac+b)b(fd)  and  r<ff1d;
    unstable focus-node if
    r>max{f(2c+ac+b)b(fd),ff1d};
    saddle focus-node if
    f(ac+b)b(fd)<r<f(2c+ac+b)b(fd)  and  r>ff1d;
    or
    r>f(2c+ac+b)b(fd)  and  r<ff1d;
    non-hyperbolic if
    r=ff1d or  r=f(2c+ac+b)b(fd).
    Px0znever stable focus-node; unstable focus-node if r<acb1;
    saddle focus-node if r>acb1; non-hyperbolic if r=acb1.
    Pxyzsink if Eq (4.36) along with Eq (4.35) hold.

     | Show Table
    DownLoad: CSV

    This research is partially supported by the Higher Education Commission of Pakistan.

    The authors declare that they have no conflicts of interest regarding the publication of this paper.

    [1] European Wind Energy Association, The European offshore wind industry-key trends and statistics 2015. European Wind Energy Association, 2016. Available from: https://www.ewea.org/fileadmin/files/library/publications/statistics/EWEA-European-Offshore-Statistics- 2015.pdf.
    [2] Bilgili M, Yasar A, Simsek E (2011) Offshore wind power development in Europe and its comparison with onshore counterpart. Renew Sust Energ Rev 15: 905–915. doi: 10.1016/j.rser.2010.11.006
    [3] Perveen R, Kishor N, Mohanty SR (2014) Off-shore wind farm development: present status and challenges. Renew Sust Energ Rev 29: 780–792.
    [4] Soukissian TH, Papadopoulos A (2015) Effects of different wind data sources in offshore wind power assessment. Renew Energ 77: 101–114. doi: 10.1016/j.renene.2014.12.009
    [5] Colmenar SA, Perera PJ, Borge DD, et al. (2016) Offshore wind energy: a review of the current status, challenges and future development in Spain. Renew Sust Energ Rev 64: 1–18. doi: 10.1016/j.rser.2016.05.087
    [6] European Wind Energy Association, Wind in power: 2015 European statistics. European Wind Energy Association, 2016. Available from: http://www.ewea.org/fileadmin/files/library/ publications/ statistics/ EWEA-Annual-Statistics-2015.pdf.
    [7] Soukissian T, Reizopoulou S, Drakopoulou P, et al. (2016) Greening offshore wind with the smart wind chart evaluation tool. Web Ecol 16: 73–80. doi: 10.5194/we-16-73-2016
    [8] Kaldellis JK, Apostolou D, Kapsali M, et al. (2016) Environmental and social footprint of offshore wind energy. Comparison with onshore counterpart. Renew Energ 92: 543–556.
    [9] Brownlee MTJ, Hallo JC, Jodice LW, et al. (2015) Attitudes toward offshore wind energy development. J Leisure Res 47: 263–284.
    [10] Westerberg V, Jacobsen JB, Lifran R (2013) The case for offshore wind farms, artificial reefs and sustainable tourism in the French Mediterranean. Tourism Manag 34: 172–183. doi: 10.1016/j.tourman.2012.04.008
    [11] DIRM Méditerranée, Document de planification: Le développement de l'éolien en mer Méditerranée. France: Préfecture maritime de la Méditerranée, Préfecture de région Provence Alpes Côte d'Azur, 2015. Available from: http://www.dirm.mediterranee.developpement-durable.gouv.fr/ IMG/pdf/ Document_de_planification_pour_transmission.pdf.
    [12] 4C Offshore, Two more French Floaters get approved! 4C Offshore, 2016. Available from: http://www.4coffshore.com/windfarms/two-more-french-floaters-get-approved!-nid4813.html.
    [13] Rodrigues S, Restrepo C, Kontos E, et al. (2015) Trends of offshore wind projects. Renew Sust Energ Rev 49: 1114–1135. doi: 10.1016/j.rser.2015.04.092
    [14] Westerberg V, Jacobsen JB, Lifran R (2015) Offshore wind farms in Southern Europe-determining tourist preference and social acceptance. Energ Res Soc Sci 10: 165–179. doi: 10.1016/j.erss.2015.07.005
    [15] Zountouridou EI, Kiokes GC, Chakalis S, et al. (2015) Offshore floating wind parks in the deep waters of Mediterranean Sea. Renew Sust Energ Rev 51: 433–448. doi: 10.1016/j.rser.2015.06.027
    [16] Onea F, Deleanu L, Rusu L, et al. (2016) Evaluation of the wind energy potential along the Mediterranean Sea coasts. Energ Explor Exploit 34: 766–792. doi: 10.1177/0144598716659592
    [17] Balog I, Ruti PM, Tobin I, et al. (2016) A numerical approach for planning offshore wind farms from regional to local scales over the Mediterranean. Renew Energ 85: 395–405. doi: 10.1016/j.renene.2015.06.038
    [18] Kotroni V, Lagouvardos K, Lykoudis S (2014) High-resolution model-based wind atlas for Greece. Renew Sust Energ Rev 30: 479–489. doi: 10.1016/j.rser.2013.10.016
    [19] Emmanouil G, Galanis G, Kalogeri C, et al. (2016) 10-year high resolution study of wind, sea waves and wave energy assessment in the Greek offshore areas. Renew Energ 90: 399–419. doi: 10.1016/j.renene.2016.01.031
    [20] Soukissian T, Karathanasi F, Axaopoulos P (2017) Satellite-based offshore wind resource assessment in the Mediterranean Sea. IEEE J Oceanic Eng 42: 73–86.
    [21] Soukissian TH (2014) Probabilistic modeling of directional and linear characteristics of wind and sea states. Ocean Eng 91: 91–110. doi: 10.1016/j.oceaneng.2014.08.018
    [22] Song M, Chen K, Zhang X, et al. (2016) Optimization of wind turbine micro-siting for reducing the sensitivity of power generation to wind direction. Renew Energ 85: 57–65. doi: 10.1016/j.renene.2015.06.033
    [23] Watson SJ (2014) Quantifying the variability of wind energy. Wires Energ Environ 3: 330–342. doi: 10.1002/wene.95
    [24] EMODnet, EMODnet Bathymetry portal. EMODnet, 2016. Available from: http://www.emodnet-hydrography.eu/.
    [25] Caralis G, Chaviaropoulos P, Ruiz Albacete V, et al. (2016) Lessons learnt from the evaluation of the feed-in tariff scheme for offshore wind farms in Greece using a Monte Carlo approach. J Wind Eng Ind Aerod 157: 63–75. doi: 10.1016/j.jweia.2016.08.008
    [26] Greek Parliament, Governmental Gazette, A' No. 149/9-8-2016, L. 4414/2016. Official Government Gazette of the Hellenic Republic.
    [27] European Commission, Official Journal of the European Union, Guidelines on State aid for environmental protection and energy 2014-2020 (2014/C 200/01). European Commission, 2014. Available from: http://eur-lex.europa.eu/legal-content/EN/TXT/PDF/?uri=CELEX:52014XC 0628(01)&from=EN.
    [28] Karathanasi FE, Soukissian TH, Axaopoulos PG (2016) Calibration of wind directions in the Mediterranean Sea. In: Proceedings of the 26th International Ocean and Polar Engineering Conference; 2016; Rhodes, Greece, 491–497.
    [29] Fisher N (1995) Statistical analysis of circular data. 1st ed. Cambridge: Cambridge University Press, 294.
    [30] Jammalamadaka R, SenGupta A (2001) Topics in circular statistics. Singapore: World Scientific Publishing Co. Pte. Ltd., 334.
    [31] Hansen FV (1993) Surface roughness lengths. White Sands Missile Range, New Mexico: U.S. Army Research Laboratory, 1–40.
    [32] Shu ZR, Li QS, He YC, et al. (2016) Observations of offshore wind characteristics by Doppler-LiDAR for wind energy applications. Appl Energ 169: 150–163. doi: 10.1016/j.apenergy.2016.01.135
    [33] Papadopoulos A, Katsafados P (2009) Verification of operational weather forecasts from the POSEIDON system across the Eastern Mediterranean. Nat Hazards Earth Syst Sci 9: 1299–1306. doi: 10.5194/nhess-9-1299-2009
    [34] Papadopoulos A, Korres G, Katsafados P, et al. (2011) Dynamic downscaling of the ERA-40 data using a mesoscale meteorological model. Mediterranean Mar Sci 12: 183–198.
    [35] Ferrier BS, Jin Y, Lin Y, et al. (2002) Implementation of a new grid-scale cloud and precipitation scheme in the NCEP Eta Model. 19th Conference on weather analysis and forecasting/15th Conference on numerical weather prediction. San Antonio: Am Meteorol Soc, 280–283.
    [36] Janjic ZI, Gerrity JP, Nickovic S (2001) An alternative approach to nonhydrostatic modeling. Mon Weather Rev 129: 1164–1178.
    [37] Janjić ZI (1994) The step-mountain Eta coordinate model: further developments of the convection, viscous sublayer, and turbulence closure schemes. Mon Weather Rev 122.
    [38] Chen F, Janjić Z, Mitchell K (1997) Impact of atmospheric surface-layer parameterizations in the new land-surface scheme of the NCEP mesoscale eta model. Bound Lay Meteorol 85: 391–421. doi: 10.1023/A:1000531001463
    [39] Lacis AA, Hansen J (1974) A parameterization for the absorption of solar radiation in the earth's atmosphere. J Atmos Sci 31: 118–133.
    [40] Schwarzkopf MD, Fels SB (1991) The simplified exchange method revisited: an accurate, rapid method for computation of infrared cooling rates and fluxes. J Geophys Res 96: 9075–9096. doi: 10.1029/89JD01598
    [41] Soukissian T, Chronis G (2000) Poseidon: a marine environmental monitoring, forecasting and information system for the Greek Seas. Mediterranean Mar Sci 1: 71–78.
  • This article has been cited by:

    1. Sabrina H. Streipert, Gail S. K. Wolkowicz, 2023, Chapter 22, 978-3-031-25224-2, 473, 10.1007/978-3-031-25225-9_22
    2. Subhadipa Das, Amar Nath Chatterjee, Sourav Rana, Complex Dynamics of a Discrete Time One Prey Two Predator System with Prey Refuge, 2025, 11, 2349-5103, 10.1007/s40819-025-01896-9
  • Reader Comments
  • © 2017 the Author(s), licensee AIMS Press. This is an open access article distributed under the terms of the Creative Commons Attribution License (http://creativecommons.org/licenses/by/4.0)
通讯作者: 陈斌, bchen63@163.com
  • 1. 

    沈阳化工大学材料科学与工程学院 沈阳 110142

  1. 本站搜索
  2. 百度学术搜索
  3. 万方数据库搜索
  4. CNKI搜索

Metrics

Article views(9793) PDF downloads(1857) Cited by(23)

Article outline

Figures and Tables

Figures(8)  /  Tables(6)

/

DownLoad:  Full-Size Img  PowerPoint
Return
Return

Catalog