PD control at Neimark–Sacker bifurcations in a Mackey–Glass system

This article focuses on a proportional-derivative (PD) feedback controller to control a Neimark–Sacker bifurcation for a Mackey–Glass system by the Euler method. It has been shown that the onset of the Neimark–Sacker bifurcation can be postponed or advanced via a PD controller by choosing proper control parameters. Finally, numerical simulations are given to confirm our analysis results and the effectiveness of the control strategy. Especially, a derivative controller can significantly improve the speed of response of a control system.


Introduction
In Mackey and Glass [1], Mackey and Glass described a physiological system as the delay differential equation (DDE): Here, p(t) denotes the density of mature cells in blood circulation, τ is the time delay between the production of immature cells in the bone marrow and their maturation for release in circulating bloodstreams, β, θ , γ , and n ≥ 3 are all positive constants and β γ > 1. (1.2) In [1], Michael C. Mackey and Leon Glass associated the onset of a disease with bifurcations in Eq. (1.1). The fluctuations were caused in peripheral, while blood cell counts in chronic granulocytic leukemia (CGL). So studying the dynamics of Eq. (1.1) is significant to medical research. Various bifurcations exist in nonlinear dynamical systems, such as complex circuits, networks, and so on. Bifurcations can be important and beneficial if they are under appropriate control. In [2], a nonlinear feedback control method was given, in which a polynomial function was used to control the Hopf bifurcation. Bifurcation control refers to the task of designing a controller to suppress or reduce some existing bifurcation dynamics of a given nonlinear system, thereby achieving some desirable dynamical behaviors [3][4][5][6][7].
In [8], a hybrid control strategy was proposed, in which state feedback and parameter perturbation were used to control the bifurcations. In [9], a discrete-time delayed feedback controller was presented and studied. Proportional-integral-derivative (PID) controller is a widely used control method for dynamics control in a nonlinear system for its superior performance [10][11][12][13][14]. The results show one can delay or advance the onset of bifurcations by changing the control parameters, including the proportional control parameter and the derivative control parameter.
The discussion on the importance of discrete-time analogues in preserving the properties of stability and bifurcation of their continuous-time counterparts has been studied by some authors. For example in [15], the authors considered the numerical approximation of a class of DDEs undergoing a Hopf bifurcation by using the Euler forward method. Motivated by the works on bifurcation control, we adopt a proportional-derivative (PD) feedback control Euler scheme in which a proportional control parameter and a derivative control parameter are used to control the Neimark-Sacker bifurcations. By selecting appropriate control parameters, we obtain that the dynamic behavior of a controlled system can be changed. Especially, a derivative controller can significantly improve the speed of response of a control system. As far as we know, the numerical controlled dynamics for DDES by a PD feedback controller has been rarely studied.
The rest of the paper is organized as follows. In the next section, for a DDE of Mackey-Glass system with PD controller, we analyze the local stability of the equilibria and existence of the Hopf bifurcation. The main results are obtained in Sect. 3. Among them, by using a PD controlled Euler algorithm, the dynamics of the numerical discrete systems are derived according to the Neimark-Sacker bifurcation theorem. In Sect. 4, by applying the theories of discrete bifurcation systems, the direction and stability of bifurcating periodic solutions from the Neimark-Sacker bifurcation of controlled delay equation are confirmed. Section 5 gives numerical examples to illustrate the validity of our results.

Hopf bifurcation in DDE via a PD controller
Under the transformation p(t) = θ x(t), Eq. (1.1) becomeṡ x * is a positive fixed point to Eq. (2.1), and x * satisfies Apply a PD controller to system (2.1) as follows: where k p is the proportional control parameter, and k d is the derivative control parameter. (2.4) The linearization of Eq. (2.4) atx = 0 is whose characteristic equation is For τ = 0, the only root of Eq. (2.6) is here k d < 1 and Let u(t) = x(τ t), then Eq. (2.1) can be written aṡ u * is a positive fixed point to Eq. (2.9), and u * satisfies Apply a PD controller to system (2.9) as follows: The linearization of Eq. (2.12) at z = 0 is whose characteristic equation is with λ =λτ for τ = 0.
Let iω be a root of Eq. (2.14) if and only if Separating the real and imaginary parts, we have According to (2.6), (2.8), (2.14), and (2.16), this is impossible if and only if k d < 1 and We obtain the following result.
Proof By differentiating both sides of Eq. (2.14) with respect to τ , we obtain Therefore, This implies that The result is confirmed. (2.3) and (2.11), we give the following statements:

Neimark-Sacker bifurcation analysis of the PD control Euler method
This section concerns the stability and bifurcation of the numerical discrete PD control system. We implement the PD control strategy [10][11][12][13]. Set Employing the Euler method to Eq. (3.1) yields the difference equation here h = 1 m , m ∈ Z, y n is an approximate value to y(nh). Providing a new variable Y n = (y n , y n-1 , . . . , y n-m ) T , we can rewrite (3.2) as where F = (F 0 , F 1 , . . . , F m ) T , and The characteristic equation of A is The equation has an m-fold root λ = 0 and a simple root λ = 1. Consider the root λ(τ ) such that |λ(0)| = 1. This root is a C 1 function of τ . For Eq. (3.7), we have If k d < 1 and (3.8) are satisfied, then Consequently, all roots of Eq. (3.7) lie in |λ| < 1 for sufficiently small τ > 0.

Lemma 3 If the step-size h is sufficiently small, k d < 1 and
then Eq. (3.7) has no root with modulus one for all τ > 0.
Proof When two roots of characteristic equation (3.7) pass through the unit circle, a Neimark-Sacker bifurcation occurs. Assume that there exists τ * such that e iω * , ω * ∈ (-π, π] is the root of characteristic equation (3.7). Then Hence, separating the real and the imaginary parts gives We get . (3.13) By virtue of the step-size h being sufficiently small, k d < 1 and (3.10) being satisfied, we obtain cos ω * > 1, which yields a contradiction. So Eq. (3.7) has no root with modulus one for all τ > 0.
from Lemmas 2 and 3, we know that Eq. (3.7) has no root with modulus one for all τ > 0. Applying Corollary 2.4 in [16], all roots of Eq. (3.7) have modulus less than one for all τ > 0. The conclusion follows.
Remark 1 Through the above conclusions of Lemmas 2-4 and Theorem 2, for the stepsize h is sufficiently small, due to β γ > nγ k p +(n-2)γ , here k d < 1, (2n)γ < k p < 0 or 0 < k p < γ (when k p = 0, β γ > n n-2 ), we can delay (or advance) the onset of a Neimark-Sacker bifurcation by choosing appropriate control parameters k p and k d .

Direction and stability of the Neimark-Sacker bifurcation in a discrete control model
In the previous section, we have verified the conditions for the Neimark-Sacker bifurcation to occur when τ = τ * k for k = 0, 1, 2, . . . , [ m- 1 2 ]. In this section we continue to study the direction of the Neimark-Sacker bifurcation and the stability of the bifurcating periodic solutions when τ = τ * 0 using the techniques from normal form and center manifold theory [17,18].
So, we can rewrite system (3.2) as (4.1) Let q = q(τ * 0 ) ∈ C m+1 be an eigenvector of A corresponding to e iω * 0 , then We also introduce an adjoint eigenvector q * = q * (τ ) ∈ C m+1 having the properties and satisfying the normalization q * , q = 1, where q * , q = m i=0 q * i q i .
Let T center denote a real eigenspace corresponding to e ±iw * 0 , which is two-dimensional and is spanned by { (q), Im(q)}, and let T stable be a real eigenspace corresponding to all eigenvalues of A T , other than e ±iw * 0 , which is (m -1)-dimensional. All vectors x ∈ R m+1 can be decomposed as where v ∈ C, vq +vq ∈ T center , and y ∈ T stable . The complex variable v can be viewed as a new coordinate on T center , so we have v = q * , x , Let a(λ) be a characteristic polynomial of A and λ 0 = e iw * 0 . Following the algorithms in [17] and using a computation process similar to that in [15,19], we have g 20 = q * , B(q, q) , g 11 = q * , B(q, q) , g 02 = q * , B(q, q) , where So, we can get an expression for the critical coefficient c 1 (τ * 0 ) (4.6) By (4.1), (4.2), and Lemma 6, we get a(e 2w * 0 i ) We obtain the stability of the closed invariant curve by applying the Neimark-Sacker bifurcation theorem [20]. The results are as follows.

Numerical simulations
The purpose of this section is to validate the effectiveness of the PD control Euler method in Sects. 2-4 by numerical examples.   Table 1, we could argue that the PD control Euler method enlarges the stable region by choosing control parameters k p < 0, k d < 1. From Figs. 1 and 5, we could obtain that the PD control Euler method narrows the stable region by choosing control parameters k p > 0, k d < 1.
At the same time, for the purpose of comparison, we choose the same the proportional control parameter k p = 0.09, the different derivative control parameters k d = 0.2 ( Fig. 1)