Bifurcation Analysis and Its Applications

Continuous dynamical systems that involve differential equations mostly contain parameters. It can happen that a slight variation in a parameter can have significant impact on the solution. The main questions of interest in this chapter are: How to continue equilibria and periodic orbits of dynamical systems with respect to a parameter? How to compute stability boundaries of equilibria and limit cycles in the parameter space? How to predict qualitative changes in system’s behavior (bifurcations) occurring at these equilibrium points? This chapter will also cover the classification of bifurcations in terms of equilibria and periodic orbits. Especially it will present the specific bifurcation called ”Hopf bifurcation” which refers to the development of periodic orbits from stable equilibrium point, as a bifurcation parameter crosses a critical value. Since the theory of bifurcation from equilibria based on center manifold reduction and Poincaré-Normal forms, the direction of bifurcations for the mathematical models will also be explained using this theory. Finally, by introducing several software packages and numerical methods this chapter will also cover the techniques to determine and continue in some control parameters all local bifurcations of periodic orbits of dynamical systems and relevant normal form computations combined with the center manifold theorem, including periodic normal forms for periodic orbits.


Introduction
Continuous dynamical systems that involve differential equations mostly contain parameters. It can happen that a slight variation in a parameter can have significant impact on the solution. The main questions of interest in this chapter are: How to continue equilibria and periodic orbits of dynamical systems with respect to a parameter? How to compute stability boundaries of equilibria and limit cycles in the parameter space? How to predict qualitative changes in system's behavior (bifurcations) occurring at these equilibrium points? This chapter will also cover the classification of bifurcations in terms of equilibria and periodic orbits. Especially it will present the specific bifurcation called "Hopf bifurcation" which refers to the development of periodic orbits from stable equilibrium point, as a bifurcation parameter crosses a critical value. Since the theory of bifurcation from equilibria based on center manifold reduction and Poincaré-Normal forms, the direction of bifurcations for the mathematical models will also be explained using this theory. Finally, by introducing several software packages and numerical methods this chapter will also cover the techniques to determine and continue in some control parameters all local bifurcations of periodic orbits of dynamical systems and relevant normal form computations combined with the center manifold theorem, including periodic normal forms for periodic orbits.
In general, in a dynamical system, a parameter is allowed to vary, then the differential system may change. An equilibrium can become unstable and a periodic solution may appear or a new stable equilibrium may appear making the previous equilibrium unstable. The value of parameter at which these changes occur is known as "bifurcation value" and the parameter that is varied is known as the "bifurcation parameter".
In this chapter, we also discuss several types of bifurcations, saddle node, transcritical, pitchfork and Hopf bifurcation. Among these types, we especially focus on Hopf bifurcation. The first three types of bifurcation occur in scalar and in systems of differential equations. The fourth type called Hopf bifurcation does not occur in scalar differential equations because this type of bifurcation involves a change to a periodic solution. Scalar autonomous differential equations can not have periodic solutions. Hopf bifurcation occurs in systems of differential equations consisting of two or more equations. This type is also referred to as a "Poincare-Andronov-Hopf bifurcation". For a given system of differential equations first we shall consider the stability and the local Hopf bifurcation. By using the Hopf bifurcation theorem we prove the occurrence of the Hopf bifurcation. And then, based on the normal form method and the center manifold reduction introduced by Hassard et al., [10], we derive the formulae determining the direction, stability and the period of the bifurcating periodic solution at the critical value of the bifurcation parameter. To verify the theoretical analysis, numerical simulations for bifurcation analysis are given in this chapter. For references see [1]- [22].
We also introduce the Hopf bifurcation for continuous dynamical systems and state the Hopf bifurcation theorem for these models. As it is well known, Hopf bifurcations occur when a conjugated complex pair of eigenvalues crosses the boundary of stability. In the time-continuous case, a limit cycle bifurcates. It has an angular frequency which is given by the imaginary part of the crossing pair. In the discrete case, the bifurcating orbit is generally quasi-periodic, except that the argument of the crossing pair times an integer gives just 2π. If we consider an ordinary differential equation (ODE) that depends on one or more parameters α where, for simplicity, we assume α to be the only parameter. There is the possibility that under variation of α nothing interesting happens to Equation (1). There is only a quantitatively different behavior. Let us define Equation (1) to be structurally stable in the case there are no qualitative changes occurring. However, the ODE (Ordinary Differential Equation) might change qualitatively. At that point, bifurcations will have occurred.
Many of the basic principles for one dimensional systems apply also for two-dimensional systems. Let us define a two-dimensional system where biologically we mostly interpret x as prey or resource and y as predator or consumer. Equilibria can be found by taking the equations equal to zero, i.e., We have three possibilities for the stability of an equilibrium. Next to the stable and unstable equilibrium, there is the saddle equilibrium. A two-dimensional stable equilibrium is attracting in two directions, while a two-dimensional unstable equilibrium is repelling in two directions. A saddle point is attracting in one direction and repelling in the other direction.
In the less formal literature saddles are often considered just unstable equilibria. A second remark is that also the dynamics of the system around the equilibria can differ.
The attracting or repelling can occur via straight orbits (a node) or via spiralling orbits (a spiral or focus). Note that, it is not possible to have a saddle focus in two dimensions. It is possible though in three of higher dimensional systems.
To prove the existence of Hopf bifurcation, we first obtain the Hopf bifurcation theorem hypothesis, i.e., the existence of purely imaginary eigenvalues of the corresponding characteristic equation with respect to the parameter α and also we prove the transversality condition at the critical value α 0 where Hopf bifurcation occurs. Then based on the normal form approach and the center manifold theory introduced by Hassard et. al, [10], we derive the formula for determining the properties of Hopf bifurcation of the model.
Finally in this chapter, to support these theoretical results, we illustrate them by numerical simulations. In numerical analysis, generally MATLAB solver packages are used to analyze the dynamics of nonlinear models. In these solvers, differential equation systems are simulated by difference equations. In this chapter, we also give some examples from biology (such as well known predator-prey models with time delay) with numerical simulations and by graphing the solutions in two or three dimensions, we illustrate the occurrence of periodic solutions.

Basic concepts of bifurcation analysis
As it is stated above, in dynamical systems, a bifurcation occurs when a small smooth change made to the parameter values (the bifurcation parameters) of a system causes a sudden "qualitative" or topological change in its behavior. Generally, at a bifurcation, the local stability properties of equilibria, periodic orbits or other invariant sets changes. It has two types; Local bifurcations, which can be analyzed entirely through changes in the local stability properties of equilibria, periodic orbits or other invariant sets as parameters cross through critical thresholds; and Global bifurcations, which often occur when larger invariant sets of the system "collide" with each other, or with equilibria of the system. They cannot be detected purely by a stability analysis of the equilibria (fixed or equilibrium points, see the next section).

Equilibrium points
In dynamical systems, only the solutions of linear systems may be found explicitly. The problem is that in general real life problems may only be modeled by nonlinear systems. The main idea is to approximate a nonlinear system by a linear one (around the equilibrium point).
Of course, we do hope that the behavior of the solutions of the linear system will be the same as the nonlinear one. But this is not always true. Before the linear stability analysis, we give some basic definitions below.
Definition (Equilibrium Point): Consider a nonlinear differential equation where f is a function mapping R n xR 3 → R n . A pointx is called an equilibrium point if there is a specificū R m such that f (x(t), u(t)) = 0 n .
Supposex is an equilibrium point (with the inputū). Consider the initial condition x(0) =x, and applying the input u(t) =ū for all t ≥ t 0 , then resulting solution x(t) satisfies for all t ≥ t 0 . That is why it is called an equilibrium point or solution.
Example: As an example, consider the logistic growth equation (the rate of population density) where x(t) denotes the population density at time t, r and K are positive constants, K is the carrying capacity. Then by setting right hand side function equal to zero, we obtain two equilibrium points x = 0 and x = K.

Linear stability analysis
Linear stability of dynamical equations can be analyzed in two parts: one for scalar equations and the other for two dimensional systems;

Linear stability analysis for scalar equations
To analyze the ODE locally about the equilibrium point x =x, we expand the function f (x) in a Taylor series about the equilibrium pointx. To emphasize that we are doing a local analysis, it is customary to make a change of variables from the dependent variable x to a local variable. Now let where it is assumed that ε(t) 1, so that we can justify dropping all terms of order two and higher in the expansion. Substituting x(t) =x + ε(t) into the RHS of the ODE yields; and dropping higher order terms, we obtain Note that dropping these higher order terms is valid since ε(t) The goal is to determine if we have growing or decaying solutions. If the solutions grows, then the equilibrium point is unstable. If the solution decays, then the fixed point is stable.
To determine whether or not the solution is stable or unstable we simply solve the ODE and get the solution as where ε 0 is a constant. Hence, the solution is growing if f (x) > 0 and decaying if f (x) < 0. As a result, the equilibrium point is stable if f (x) < 0, unstable if f (x) > 0 as it is stated in the following theorem.
Theorem: Suppose for scalar differential equation the derivative function f is continuous on an open interval I where the equilibrium pointx I.
If the equilibrium point is stable and in addition lim t→∞ x(t) =x, then it is called asymptotically stable equilibrium point.

Example:
Determine the stability of the fixed points to the Logistic growth equation We first find the equilibrium points by setting f (N) = 0, which yields two pointsN = 0, K.
Next we compute f (N) and evaluate it at the equilibrium points.

Linear stability analysis for systems
Consider the two dimensional nonlinear system and suppose that (x,ȳ) is a steady state (equilibrium point), i.e., f (x,ȳ) = 0 and g(x,ȳ) = 0. Now let's consider a small perturbation from the steady state (x,ȳ) where u and v are understood to be small as u 1 and v 1. Is is natural to ask whether u and v are growing or decaying so that x and y will move away form the steady state or move towards the steady states. If it moves away, it is called unstable equilibrium point, if it moves towards the equilibrium point, then it is called stable equilibrium point.As in scalar equations, by expanding the Taylor's series for f (x, y) and g(x, y); Since u and v are assumed to be small, the higher order terms are extremely small, we can neglect the higher order terms and obtain the following linear system of equations governing the evolution of the perturbations u and v, where the matrix f x f y g x g y is called Jacobian matrix J of the nonlinear system.The above linear system for u and v has the trivial steady state (u, v) = (0, 0), and the stability of this trivial steady state is determined by the eigenvalues of the Jacobian matrix, as follows: Theorem : An equilibrium point (x,ȳ) of the differential equation is stable if all the eigenvalues of J, the Jacobian evaluated at (x,ȳ) have negative real parts. The equilibrium point is unstable if at least one of the eigenvalues has a positive real part.
As a summary, Asymptotically stable : A critical point is asymptotically stable if all eigenvalues of the jacobian matrix J are negative, or have negative real parts.

Unstable:
A critical point is unstable if at least one eigenvalue of the jacobian matrix J is positive, or has positive real part.

Stable (or neutrally stable) :
Each trajectory move about the critical point within a finite range of distance.

Definition(Hyperbolic point):
The equilibrium is said to be hyperbolic if all eigenvalues of the jacobian matrix have non-zero real parts.
Hyperbolic equilibria are robust(i.e., the system is structurally stable): Small perturbations of order do not change qualitatively the phase portrait near the equilibria. Moreover, local phase portrait of a hyperbolic equilibrium of a nonlinear system is equivalent to that of its linearization. This statement has a mathematically precise form known as the Hartman-Grobman. This theorem guarantees that the stability of the steady state (x,ȳ) of the nonlinear system is the same as the stability of the trivial steady state (0, 0) of the linearized system.

Definition(Non-Hyperbolic point):
If at least one eigenvalue of the Jacobian matrix is zero or has a zero real part, then the equilibrium is said to be non-hyperbolic.
Non-hyperbolic equilibria are not robust (i.e., the system is not structurally stable). Small perturbations can result in a local bifurcation of a non-hyperbolic equilibrium, i.e., it can change stability, disappear, or split into many equilibria. Some refer to such an equilibrium by the name of the bifurcation (See the section below).
Example: Consider the following nonlinear autonomous system The equilibria are the points (x,ȳ) = (0, 0) and (x,ȳ) = (1, 2) and the Jacobian matrix is We compute the Jacobian at the equilibrium point (0,0) where J(0, 0) = 0 −1 2 0 which implies that the eigenvalues are purely imaginary by solving the characteristic equation Since the point is non-hyperbolic, the linearized system can not tell about the stability. Later on we will show that this is a center.
For the equilibrium point (1, 2), the Jacobian J(1, 2) = 2 0 0 −1 thus the point is locally where one of the eigenvalues is strictly positive). Since it is a hyperbolic equilibrium point, the stability of fixed point is the same as in the linearized system. So it is also unstable.

Bifurcation in one dimension
We may assume that n = 1 so that f : R 2 → R 1 , and x 0 (λ) is a real valued analytic function of λ provided are solutions of The most common bifurcation types are illustrated by the following examples.

Saddle-node bifurcation
A saddle-node bifurcation or tangent bifurcation is a collision and disappearance of two equilibria in dynamical systems. In autonomous systems, this occurs when the critical equilibrium has one zero eigenvalue. This phenomenon is also called fold or limit point bifurcation.
Now consider the dynamical system defined by An equilibrium solution (where if a > 0, then we have two real solutions.
We now consider each of the two solutions for a > 0, and examine their linear stability in the usual way. First, we add a small perturbation: Substituting this into the equation yields and since the term in brackets on the RHS is trivially zero, therefore which has the solution ε(t) = A exp(−2xt).
As sketched in the "bifurcation diagram" below, therefore, the saddle node bifurcation at a = 0 corresponds to the creation of two new solution branches. One of these is linearly stable, the other is linearly unstable.

Transcritical bifurcation
In a transcritical bifurcation, two families of fixed points collide and exchange their stability properties. The family that was stable before the bifurcation is unstable after it. The other fixed point goes from being unstable to being stable.  Again, a and b are control parameters. We can find two steady states (x = 0) to this system We now examine the linear stability of each of these states in turn, following the usual procedure.
For the statex 1 , we add a small perturbation Therefore, perturbations grow for a > 0 and decay for a < 0. So It can be easily seen that the bifurcation point a = 0 corresponds to an exchange of stabilities between the two solution branches.

The pitchfork bifurcation
In pitchfork bifurcation one family of fixed points transfers its stability properties to two families after or before the bifurcation point. If this occurs after the bifurcation point, then pitchfork bifurcation is called supercritical. Similarly, a pitchfork bifurcation is called subcritical if the nontrivial fixed points occur for values of the parameter lower than the bifurcation value. In other words, the cases in which the emerging nontrivial equilibria are stable are called supercritical whereas the cases in which these equilibria are called subcritical.
Consider the dynamical system As usual, a and b are external control parameters. Steady states, for which x = 0 are as follows: Note that the equilibrium pointsx 2 andx 3 only exist when a > 0 if b > 0 and for a < 0 if b < 0.
As usual, we now examine the linear stability of each of these steady states in turn. (This can be done for a general b). First we write the perturbation forx 1 = 0, with the solution e(t) = A exp(at).
So we see that the statex 1 = 0 is linearly unstable if a > 0, the statex 2 = 0 is linearly stable if a < 0.
For the states x =x 2 and x =x 3 , settinḡ Thus it is obvious that the statesx 2 andx 3 are linearly stable if a > 0, the statesx 2 andx 3 are linearly unstable if a < 0,

Hopf bifurcation
Definition : A Hopf or Poincare-Andronov-Hopf bifurcation is a local bifurcation in which a fixed point of a dynamical system loses stability as a pair of complex conjugate eigenvalues of linearization around the fixed point cross the imaginary axis of the complex plane.

Hopf bifurcation theorem
Consider the two dimensional system where τ is the parameter and suppose that (x (τ) , y (τ)) is the equilibrium point and α (τ) ± iβ (τ) are the eigenvalues of the Jacobian matrix which is evaluated at the equilibrium point.
In addition let's assume that the change in the stability of the equilibrium point occurs at τ = τ * where α (τ * ) = 0. First the system is transformed so that the equilibrium is at the origin and the parameter τ at τ * = 0 gives purely imaginary eigenvalues. System (7) is rewritten as follows; The linearization of the system (7) about the origin is given by dX dt = J(τ)X, where X = x y and is the Jacobian matrix evaluated at origin.
Then in any open set U containing the origin in R 2 and for any τ 0 > 0, there exists a valuē τ, |τ| < τ 0 such that the system of differential equations (8) has a periodic solution for τ =τ in U. (Allen, L.J.S).
Note: The Hopf bifurcation requires at least a two dimensional differential equation system to appear.

Definition:
The bifurcation stated in the Hopf bifurcation theorem is called "supercritical" if the equilibrium point (0, 0) is asymptotically stable when τ = 0 (at the bifurcation point) and it is called "subcritical" if the equilibrium point (0, 0) is negatively asymptotically stable (as t → −∞) when τ = 0.
In a supercritical Hopf bifurcation, the limit cycle grows out of the equilibrium point. In other words, right at the parameters of the Hopf bifurcation, the limit cycle has zero amplitude, and this amplitude grows as the parameters move further into the limit-cycle. (See the figure below) However in a subcritical Hopf bifurcation, there is an unstable limit cycle surrounding the equilibrium point, and a stable limit cycle surrounding that. The unstable limit cycle shrinks down to the equilibrium point, which becomes unstable in the process. For systems started near the equilibrium point, the result is a sudden change in behavior from approach to a stable focus, to large-amplitude oscillations.(See the figure below).
Example: Consider the two dimensional system where α is a parameter. When we compute the equilibrium points depending on parameter α; If α < 0, then there is no x-nullclines, hence the system has no equilibrium points.
If α > 0, then the system has two equilibrium points (− √ α, 0) and ( √ α, 0). Then the Jacobian matrix is where at the equilibrium points For α = 0, there will be a line equilibrium (since one of the eigenvalues is zero) and for α > 0, the point (− √ α, 0) is a sink and ( √ α, 0) is a saddle point so that α = 0 is the bifurcation point for this differential equation system.
where α is the bifurcation parameter. We can easily show that the conditions of the Hopf Bifurcation theorem hold.
In this system f 1 and f 2 are zero. Then the Jacobian matrix is for which the eigenvalues are λ 1,2 = α ± i where Reλ(α) = α and the imaginary part Imλ(α) = ±1.It follows that Reλ(0) = 0 and Imλ(0) = 0. and also Hence, we conclude that there exists a periodic solution for α = 0 in every neighborhood of origin.

Center Manifold teorem
Let f (0) = 0, for the dynamical system and let the eigenvalues of the Jacobian matrix be λ 1 , λ 2 , ..., λ n . Suppose that, the real parts of the eigenvalues are zero and if not, suppose there are n + numbers of eigenvalues with Re λ > 0, n 0 number of eigenvalues with Re λ = 0 and n − number of eigenvalues with Re λ < 0. Let T c be the eigenspace on imaginary axis corresponding to n 0 eigenvalues. The eigenvalues on the imaginary axis (Re λ = 0) are called the critical eigenvalues as on the eigenspace T c . And suppose the function ϕ t denote the flow corresponding to the equation (11).
With these assumption, we state the Center Manifold theorem as follows;

Theorem: (Center Manifold theorem)
There exists a locally invariant C ∞ center manifold W c loc (0) such that such that the dynamics of the system (where A c and A s are are the blocks in the canonical form whose diagonals contain the eigenvalues with Re λ = 0 and Re λ < 0; respectively) restricted to the center manifold are given by And the manifold W c loc is called center manifold. Remark: Center manifolds are not unique.

Center Manifold reduction for two dimensional systems
Example: Consider the two dimensional system of differential equations The only equilibrium point is (0, 0), we linearize around that and obtain and on the other hand, Comparing the two expressions, we deduce that a = −1, b = 0, c = −2, d = 0 and the center manifold reduction takes the form Hence for the last equation x = 0 is asymptotically stable and therefore (0, 0) is asymptotically stable for the original system.
Example: Consider the two dimensional system Again (0, 0) is an equilibrium point and the jacobian matrix for the linearized system is from which we deduce that a = 1 and b = 0 and For this reduced equation, x = 0 is unstable and hence, (0, 0) is also unstable for the original system.

Center Manifold reduction for Hopf bifurcation
The aim of this section is to give a formal framework for the analytical bifurcation analysis of Hopf bifurcations in delay differential equations with a single fixed time delay τ to be chosen as a bifurcation parameter. Characteristic equations of the delay differential equation form (14) are often studied in order to understand changes in the local stability of equilibria of certain delay differential equations. It is therefore important to determine the values of the delay at which there are roots with zero part. We give a general formalization of these calculations and determine closed form algebraic equations where the stability and amplitude of periodic solutions close to bifurcation can be calculated.
We shall determine the direction of Hopf bifurcation and the stability of the bifurcating periodic solutions by applying the normal form theory and the center manifold theorem by Hassard et al., [10], and throughout this section, we assume that the three dimensional system of delay differential equations (14) undergoes Hopf bifurcations at the positive equilibrium (N * 0 , P * 0 , S * 0 ) at τ = τ k , and iω 1 is the corresponding purely imaginary root of the characteristic equation at the positive equilibrium (N * 0 , P * 0 , S * 0 ). For the sake of simplicity, we use the notation iω for iω 1 .
We first consider the system (14) by the transformation where By Riesz representation theorem, there exists a function η(θ, μ) of bounded variation for θ [−1, 0], such that Indeed we may take where δ is the Dirac delta function.
To determine g 21 , we need to compute W 20 (θ) and W 11 (θ). By (15) and (18), we havė Note that on the center manifold C 0 near to the origin, Thus we obtain, By using (19), for θ [−1, 0), Comparing the coefficients with (20), we obtain the following From (22) and (24) and the definition of A, we geṫ 1 , E 1 ) R 3 is a constant vector. From the definition of A and (22), we obtain where dη(θ) = η(θ, 0). Next we compute W 20 (θ) and W 11 (θ) from (25) and (26) and determine the following values to investigate the qualities of bifurcating periodic solution in the center manifold at the critical value τ k . For this purpose, we express the direction of Hopf bifurcation in terms of g ij s and eigenvalues λτ k . And then we can evaluate the following values; and we state this as in the following theorem.
Theorem : μ 2 determines the direction of Hopf bifurcation; if μ 2 > 0, then the Hopf bifurcation is supercritical and the bifurcating periodic solutions exist for τ > τ 0 , if μ 2 < 0, then the Hopf bifurcation is subcritical and the bifurcating periodic solutions exist for τ < τ 0 . β 2 determines the stability of the bifurcating periodic solutions; bifurcating periodic solutions are stable if β 2 < 0, unstable if β 2 > 0. T 2 determines the period of the bifurcating solution; the period increases if T 2 > 0, it decreases if T 2 < 0.
In the following section, we shall give a numerical example to verify the theoretical results.
By the theory of Hassard et al., [10], as it is discussed in previous section, we also determine the direction of Hopf bifurcation and the other properties of bifurcating periodic solutions. From the formulae in Section 5.2 we evaluate the values of μ 2 , β 2 and T 2 as from which we conclude that Hopf bifurcation of system (29) occurring at τ 0 = 1.5663 is supercritical and the bifurcating periodic solution exists when τ crosses τ 0 to the right, and also the bifurcating periodic solution is stable.
In computer simulations, the initial conditions are taken as (N 0 , P 0 , S 0 ) = (50, 25, 25) and MATLAB DDE (Delay Differential Equations) solver is used to simulate the system (29). We first take τ = 1.5 < τ 0 and plot the density functions N(t), P(t) and S(t) in Figs.9,10,11 respectively which shows the positive equilibrium is asymptotically stable for τ < τ 0. . Moreover in Fig.12, we illustrate the asymptotic stability in three dimension.
However in Figs.13,14,15 and 16 below, we take τ = 2 > τ 0 sufficiently close to τ 0 which illustrates the existence of bifurcating periodic solutions from the equilibrium point E * .

Consider the Van der Pol equation
and convert this into a system and check the stability of fixed points.
6. Find the equilibrium points and identify the bifurcation in the following system, and sketch the appropriate bifurcation diagram and phase portraits: Then compute the extended center manifold near the bifurcation point by choosing τ as bifurcation parameter.

7.
Show that the following system is structurally unstable, dx dt = x + y − x 2 + y 2 , dy dt = −2x − y + xy, and the following system is structural stable. Explain our reason.

9.
For r ∈ R, consider the differential equation on the real line.
(a) Show that x * = 0 is a fixed point for any value of the parameter r, and determine its stability. Hence identify a bifurcation point r 1.
(b) Show that for certain values of the parameter r there are additional fixed points. For which values of r do these fixed points exist? Determine their stability and identify a further bifurcation points r 2 .
(c)Using a Taylor expansion of the differential equation above, determine the normal form of the bifurcation at r 1 . What type of bifurcation takes place.  (c) Sketch the phase portrait in the (x, y) plane, including trajectories through (1, 0) and (2, 0). Which fixed point does the trajectory through (2, 0) approach? 11. Consider the Lorenz system (the model of heat convection by Rayleigh-Benard occurring in the earth's atmosphere) in three dimension as follows; where σ, r, b are constants. Perfom the stability analysis of this nonlinear system.

y (t) = −cy(t) + βx(t)y(t),
with time delay τ and positive constants r, b, α, c and β. By choosing τ as bifurcation parameter, check whether bifurcating periodic solutions occur around the equilibrium points or not.