## Abstract

The dynamics of a single breathing localized structure in a three-component reaction–diffusion system subjected to time-delayed feedback is investigated. It is shown that variation of the delay time and the feedback strength can lead either to stabilization of the breathing or to delay-induced periodic or quasi-periodic oscillations of the localized structure. A bifurcation analysis of the system in question is provided and an order parameter equation is derived that describes the dynamics of the localized structure in the vicinity of the Andronov–Hopf bifurcation. With the aid of this equation, the boundaries of the stabilization domains as well as the dependence of the oscillation radius on delay parameters can be explicitly derived, providing a robust mechanism to control the behaviour of the breathing localized structure in a straightforward manner.

## 1. Introduction

Starting with the work of Ott *et al.* [1], a variety of different techniques for controlling unstable or chaotic states in complex systems have been developed within the past decade. Among other control methods, the time-delayed feedback scheme [2] (also referred to as Pyragas control or time-delay autosynchronization) has proved itself to be an efficient tool, allowing a non-invasive stabilization of unstable periodic orbits of dynamical systems (see [3] and references therein). Time-delayed feedback control has been successfully applied to a broad variety of spatially extended systems, including plasma physics [4], nonlinear optics [5–8], electrochemical [9] and neural systems [10,11].

In particular, control of the dynamics of spatio-temporal patterns in reaction–diffusion systems has been of increasing interest in recent years. I mention only delay-induced turbulent structures in a diffusive Hutchinson equation [12], delay-modulated oscillatory hexagonal superlattices and stripes in a Brusselator model [13], control of spatio-temporal patterns in a Gray–Scott model [14], spatio-temporal patterns in a prey–predator plankton system [15], dynamics of Turing patterns in the Lengyel–Epstein system under time-delayed feedback [16] or moving localized and periodic structures in the FitzHugh–Nagumo model [17]. Quite recently, the influence of delayed feedback on the stability properties of a single stationary localized structure in a three-component reaction–diffusion system was studied in detail [18]. It was shown that the presence of the feedback force can induce complex dynamical behaviour of the localized solution, leading among other things to the formation of moving or breathing structures. However, the opposite problem of stabilization of a certain unstable localized state is still not understood to a large extent.

In this paper, I am interested in the influence of time-delayed control on the dynamics of breathing localized structures in a three-component reaction–diffusion system with one activator and two inhibitors. I shall show that a variation of the delay time and the feedback strength can indeed lead to stabilization of the breathing. In addition, more complex delay-induced periodic or quasi-periodic oscillations of the localized structure can also be found. In order to understand the impact of the delayed feedback term on the dynamics of the breathing localized structure, I derive an order parameter equation in the vicinity of the bifurcation point where oscillatory dynamics sets in. The desired equation is subject to a nonlinear delay–differential equation, explicitly describing the temporal evolution of the localized structure.

## 2. Linear stability analysis

As mentioned above, here I am interested in the dynamics of breathing localized structures in a three-component reaction–diffusion system with one activator and two inhibitors:
2.1
Here *u*=*u*(** r**,

*t*) is the activating component, whereas

*v*=

*v*(

**,**

*r**t*) and

*w*=

*w*(

**,**

*r**t*) denote the inhibiting components, . The coefficient λ in the polynomial nonlinear function

*f*(

*u*)=

*λu*−

*u*

^{3}is positive as well as the diffusion coefficients

*D*

_{u},

*D*

_{v}and

*D*

_{w}of the corresponding components and the dimensionless constants

*η*and

*θ*, representing the ratios of the characteristic time scales of both inhibitors

*v*and

*w*with respect to that of the activator. The constants

*κ*

_{3}and

*κ*

_{4}are also positive, whereas

*κ*

_{1}violates the inversion symmetry and has arbitrary sign. Finally,

*τ*denotes the delay time, whereas the parameter

*α*is the delay strength. Note that the time-delayed feedback term is introduced in such a way that the corresponding coupling matrix is a unit one [18]. In the absence of the delayed feedback, system (2.1) was first introduced [19,20] as an extension of the phenomenological model for a planar DC gas-discharge system with high ohmic semiconductor electrode. On the other hand, system (2.1) can be considered in a more general context as a three-component extension of the FitzHugh–Nagumo system for nerve pulse transmission [11,21–23] or a model system of a pattern-forming chemical reaction in a microemulsion [24,25]. However, in this study I have no specific application in mind and provide a general analysis of the system in question, which can be later addressed to the specific applications mentioned above. Note that, from a practical perspective, the usage of an identity control scheme may be quite restrictive. Notwithstanding that, a close analytical treatment of the simplest case of the control scheme enables one to gain deeper insights into the underlying stabilization mechanisms and to get some ideas about the impact of the time-delayed feedback on the dynamical properties of the localized structure.

From now on, I consider the general form of the reaction–diffusion system (2.1),
2.2
where ** q**=

**(**

*q***,**

*r**t*)=(

*u*(

**,**

*r**t*),

*v*(

**,**

*r**t*),

*w*(

**,**

*r**t*))

^{T}is a vector function, , is a nonlinear reaction–diffusion operator and

*E*denotes an identity coupling matrix. I am interested in the dynamics of a two-dimensional stationary single localized solution

*q*_{0}(

**) of system (2.2), which exists and which is stable in an appropriate range of parameters [26]. In the absence of the delay term, i.e. for**

*r**α*=0, this stationary localized structure can lose its stability with the change of one or more control parameters, e.g.

*η*or

*θ*. Typical destabilization scenarios include drift bifurcation, leading to a motion of localized structures [27,28], the Hopf bifurcation (also called Andronov–Hopf bifurcation), where the so-called breathing localized structures are formed [29], or a non-trivial combination of both instabilities in the vicinity of a codimension-two bifurcation point [30].

For *α*=0, linear stability of the stationary solution *q*_{0}(** r**) can be analysed by means of the ansatz

**(**

*q***,**

*r**t*)=

*q*_{0}(

**)+**

*r***(**

*φ***) exp(**

*r**μt*), leading to the linear eigenvalue problem 2.3 where the linear operator stands for a linearization of the operator around the stationary solution

*q*_{0},

*μ*is the set of eigenvalues and

**(**

*φ***) are the corresponding eigenfunctions. Note that generally the operator is not self-adjoint, so its eigenvalues and eigenfunctions are typically complex. The destabilization scenario I am interested in here is that a pair of complex-conjugated eigenvalues pass through the imaginary axis as a control parameter exceeds some critical value. In what follows, I use the time constant**

*r**θ*as a control parameter. If

*θ*is below a critical value

*θ*

_{c}, both inhibitors almost instantaneously adapt to the current distribution of the activator and the stationary solution remains stable. Therefore, a mechanism for imposing instabilities is delayed inhibition, where one or several inhibitors are too slow to adapt to activator changes [23]. Indeed, if

*θ*>

*θ*

_{c}, the addition of a symmetric perturbation in the form of the unstable breathing mode to the localized structure leads to slightly increased concentration in the centre, that is, both concentrations become narrower than their counterparts in the stationary solution. Therefore, they spread to the sides due to diffusion. While the inhibitor is slow, it strongly increases its concentration to overcompensate the losses caused by diffusion, whereas the activator decreases only slightly. That is, after some time one gets an excess amount of the inhibitor, distributed over the whole localized structure. It causes a decay of the whole solution, especially on the tails of the activator. Hence both peaks become narrower. On the other hand, the losses due to diffusion cause the decay of the inhibitor in the centre. As a consequence, the activator can again slightly increase its concentration, and new breathing cycles begin. Figure 1 shows the typical behaviour of the localized structure beyond the Andronov–Hopf bifurcation point

*θ*=

*θ*

_{c}calculated by a numerical integration of system (2.1). In figure 1

*a*, a stationary localized solution bifurcates to a breathing localized structure [29], which oscillates with a constant amplitude. Figure 1

*b*shows another possible scenario, where the amplitude of oscillations of the localized structure increases with time. Here, slow inhibitors cannot suppress a strong increase in the activator concentration in the course of time, leading to a collapse of the solution. The former case is referred to as supercritical Hopf bifurcation, whereas the latter scenario is called subcritical. I focus first on the supercritical case, that is, an increase of the control parameter beyond the critical value

*θ*

_{c}leads to the formation of the breathing localized structure. The goal now is to investigate the behaviour of the breathing localized structure in the presence of the time-delayed feedback term and explore whether this type of oscillating solution can be stabilized with this kind of control. The stationary localized solution

*q*_{0}exists for all values of the delay strength

*α*and is not affected by the delay term. However, its stability may change [18]. For

*α*≠0, the linear stability of

*q*_{0}is given by the eigenvalue problem 2.4 with the same set of eigenfunctions

**as in (2.3), as the linearization operator commutes with the identity coupling matrix. The complex eigenvalues λ can be found in terms of the Lambert function**

*φ**W*

_{m}, [6,31–33], as

A stabilization of the unstable localized solution can be achieved for such values of *α* and *τ* where Re(λ) is negative. Separating the real and imaginary parts of the last equation and solving the obtained system for Re(λ)=0, the following solvability condition for the instability threshold can be derived:
2.5
As Re(*μ*)>0 beyond the bifurcation point *θ*=*θ*_{c}, equation (2.5) possesses non-trivial solutions only for negative values of the delay strength *α*, whereas a maximal possible value of *α* to stabilize the unstable periodic localized solution is given by *α*=−Re(*μ*)/2. In order to find the shape of the domains of stabilization, the solvability condition (2.5) is solved numerically for different values of *τ* and *α* as shown in figure 2. The obtained bifurcation diagram clearly indicates the influence of the time-delayed feedback term on the dynamics of the breathing localized solution: the coloured domains appertain to *Re*(λ)≤0, where the stabilization is successful, whereas outside of these domains no stabilization takes place and the induced instabilities are caused by eigenvalues with Re(λ)>0 and (in general) with non-vanishing imaginary parts.

## 3. Direct numerical simulations

In order to verify results obtained from linear stability analysis, direct numerical simulations of the time evolution of the localized breathing solution of system (2.1) have been performed. The system parameters were chosen to be in the supercritical regime (as shown in figure 1*a*). The calculations were performed on the rectangular domain *Ω*=[−*L*,*L*]×[−*L*,*L*] with periodic boundary conditions using a pseudospectral method, whereas a fourth-order Runge–Kutta scheme is employed for the time stepping. At the start of the simulation, the feedback term was switched off, that is, a breathing localized structure, located in the centre of the domain *Ω* and widely separated from its boundaries, emerges. As soon as the amplitude of oscillations reaches a constant value, the time-delayed term was switched on. Note that, as the localized solution is situated far from domain boundaries, boundary interaction effects are infinitesimal and play no role in the dynamical behaviour. In the case of Neumann boundary conditions, boundary effects are negligible, too, for the same geometrical settings. Three simulation examples are presented in figure 3, where space–time plots for different values of *α* and *τ* are shown (figure 3*a*–*c*). The white line indicates the time moment at which the time-delayed term was switched on. Figure 3*d*–*f* represent zoomed views of the regions where delayed feedback is applied. One can see that, for parameters within stabilization domains, the amplitude of oscillations decreases after the activation of the time-delayed term and stationary localized structure emerges (figure 3*a*,*d*). On the other hand, outside of these domains, no stabilization takes place and solutions that oscillate with slightly different amplitude can be found (figure 3*b*,*e*). However, the induced oscillatory instability scenario can be more complex, as shown in figure 3*c*,*f*, where a quasi-periodic oscillating localized structure arises as the time-delayed term is switched on. In order to characterize the behaviour of the obtained breathing and quasi-periodic breathing solutions, the power spectrum in terms of Fourier transform *F* was calculated for the activator component *u*(** r**,

*t*) (figure 3

*g*,

*h*). Here, one can clearly see the difference between the power spectra of a single breathing solution, presented in figure 3

*b*,

*e*, and a quasi-periodic breathing (figure 3

*c*,

*f*). In the first case, the periodic signal gives peaks at the fundamental frequency

*ω*and its harmonics, as shown in figure 3

*g*. In the case of quasi-periodic breathing, several peaks at linear combinations of two or more irrationally related frequencies can be observed (figure 3

*h*).

## 4. Order parameter equation

The next goal is to try to understand the results obtained from the linear stability analysis and supported by numerical simulations of the system in question from the point of view of bifurcation theory. As mentioned above, the instability scenario of interest corresponds to the situation where a pair of complex-conjugated eigenvalues pass through the imaginary axis as one gradually changes the control parameter *θ*. That is, for some critical value *θ*=*θ*_{c}, the corresponding eigenvalues are purely imaginary, i.e. λ=±i*ω*. Now, if we increase the control parameter *θ*=*θ*_{c}+*ε*, *ε*≪1, the real-valued vector function ** q**(

**,**

*r**t*) can be represented as [29,30] 4.1 Here,

**(**

*ξ**t*) is a slowly varying complex amplitude of the critical eigenfunction

*φ*(

**), corresponding to the eigenvalue λ=±i**

*r**ω*at the bifurcation point

*θ*=

*θ*

_{c}. In addition,

*ξ*_{2}(

**,**

*r**t*) and

*ξ*_{0}(

**,**

*r**t*) are for the contributions of the second and zero harmonics, respectively. The specific choice of the perturbation in the ansatz (4.1) becomes apparent if one takes a look at the power spectrum calculated by means of Fourier transform

*F*(see figure 1

*c*, where the power spectrum for the activator perturbation

*u*−

*u*

_{0}is shown). Here, one can clearly see that, apart from the fundamental frequency

*ω*, only zero and second harmonics impact on the spectrum of the oscillating solution. The goal is to write down an ordinary differential equation for the amplitude

**(**

*ξ**t*) of the unstable mode

**(**

*φ***) that describes the behaviour of the single localized structure in the vicinity of the Andronov–Hopf bifurcation point. For this purpose, substitute equation (4.1) into equation (2.2), equate the terms with frequency**

*x**ω*and obtain 4.2 Here, and , so that whereas the overline stands for complex conjugate. In addition, equating the terms with frequencies 2

*ω*and 0, respectively, and neglecting contributions from the delayed terms, one gets two equations describing the time evolution of the amplitudes

*ξ*_{2}and

*ξ*_{0}of the stable modes: 4.3 and 4.4 As all critical modes except for the breathing mode

**are assumed to be stable, one can apply the adiabatic approximation to equations (4.3) and (4.4), yielding two solvability conditions with respect to unknown functions**

*φ*

*X*_{2}and

*X*_{0}: 4.5 and 4.6 Note that equations (4.3) and (4.4) were obtained neglecting contributions of the delayed terms. However, this simple approximation may be reasonable, while adiabatic elimination excludes from consideration the time evolution of stable modes. The functions

*X*_{2}and

*X*_{0}live in the same space as

**and are connected to the amplitudes**

*φ*

*ξ*_{2}and

*ξ*_{0}as [34,35] Now, in order to write the desired equation for

**, one needs to project onto the corresponding mode**

*ξ*

*φ*^{†}of the adjoint operator . The projection yields 4.7 where complex coefficients

*a*

_{1}and

*a*

_{2}can be expressed as Here, 〈⋅|⋅〉 denotes the scalar product defined in terms of full spatial integration over the considered domain.

Note that, in general, the analytical calculation of the eigenfunctions *φ*^{†} of the adjoint operator is difficult, but in the case of the reaction–diffusion system (2.1) it is possible using the relation [28,36]
4.8
where *M*(*η*,*θ*) is a diagonal matrix, defined as
For example, the coefficient standing in the numerator of *a*_{1} can be calculated as [29]
where *M*_{c}=*M*(*η*,*θ*_{c}) and *M*_{θ}=∂*M*(*η*,*θ*)/∂*θ*|_{θ=θc}. That is, using relation (4.8) one can calculate all scalar products in equation (4.7) in terms of the critical eigenfunction ** φ**=(

*φ*

_{u},

*φ*

_{v},

*φ*

_{w})

^{T}of the linearization operator .

A nonlinear delay–differential equation (4.7) is the desired order parameter equation, describing the behaviour of a single localized structure in the vicinity of the bifurcation point *θ*=*θ*_{c}. Note that a similar equation was obtained [37] for a model of Doppler's autodyne, described by the van der Pol–Duffing generator with additional delayed feedback.

For the vanishing feedback force (*α*=0), equation (4.7) is reduced to a normal form of a Hopf bifurcation [29,38]. That is, for *α*=0 a trivial solution ** ξ**=0, corresponding to the case of a stable stationary localized structure, is stable for Re(

*a*

_{1})<0, whereas for Re(

*a*

_{1})>0 a non-trivial periodic solution

**(**

*ξ**t*)=

*R*

_{0}e

^{iν0t}can be found, which is stable for Re(

*a*

_{2})<0 and unstable for Re(

*a*

_{2})>0 [29,38]. The former case corresponds to supercritical Hopf bifurcation, where the growth of the unstable mode is stabilized by a limit cycle as shown in figure 1

*a*, whereas in the latter case the bifurcation is subcritical (figure 1

*b*).

As indicated above, interest lies in the supercritical case, that is, in what follows Re(*a*_{1})>0 and *Re*(*a*_{2})<0. Now the aim is to find out the influence of the time-delayed feedback term on the dynamics of the periodic solution given by the latter relation. For non-vanishing feedback strength, *α*≠0, the exponential ansatz for the complex amplitude, ** ξ**(

*t*)=

*R*e

^{iϕ}, leads to the following system of coupled nonlinear real-valued delay–differential equations: 4.9 and 4.10 This system possesses a stationary solution, corresponding to a periodic motion on a limit cycle of the radius

*R*=

*R*

_{c}. Here,

*ϑ*

_{0}is an arbitrary constant and 4.11 and 4.12 One can easily show that for any fixed value of the delay strength

*α*, the maximal value of the radius

*R*

_{c}as a function of the delay time

*τ*is achieved at

*τ*=

*πk*/(

*ϑ*+

*ω*), , which corresponds to the radius

*R*

_{0}for vanishing delay term. That is, for these values of

*τ*, the time-delayed feedback control procedure fails. However, for other values of

*τ*and

*α*, equations (4.11) and (4.12) allow the control of the value of

*R*

_{c}in a straightforward manner: indeed, the value of

*R*

_{c}is smaller than

*R*

_{0}and can be directly calculated from system (4.11) and (4.12) (figure 4

*a*), where the dependence of the radius

*R*=

*R*

_{c}on the delay time

*τ*is shown for

*α*=−0.2.

Here, non-trivial values of *R*_{c} correspond to a localized structure, breathing with a constant amplitude (cf. figure 3*b*), whereas the vanishing *R*_{c} corresponds to stabilization of the breathing due to time-delayed feedback (figure 3*a*). Solving system (4.11) and (4.12) for *R*_{c}=0, the critical delay time *τ*, which is necessary to achieve stabilization, can be explicitly found as
4.13
One can see that, for a given value of *α*, the critical delay time depends on the parameters *a*_{1} and *ω*, defined at the bifurcation point *θ*_{c}, as well as on the distance to the bifurcation point *ε*. That is, the coefficients of the order parameter equation (4.7), derived at the bifurcation point, provide the full information concerning stabilization threshold (4.13) in the vicinity of the bifurcation point without any need to calculate the spectrum of the linearization operator as was done for the solvability condition (2.5).

Figure 4*b* showsQ1 the critical delay time (4.13) as a function of the delay strength *α* (blue solid lines), whereas red filled dots correspond to the stabilization threshold obtained from the solvability condition (2.5). Here, coloured regions indicate areas in the (*α*,*τ*) plane where stabilization of the breathing localized structure is possible, i.e. the amplitude of the oscillations *R* is zero inside of coloured domains and equals *R*_{c} outside. However, outside of the stabilization domains, the stationary solution (4.11) and (4.12) of system (4.9) and (4.10) can become unstable as one changes the delayed feedback parameters. This situation corresponds to the instability scenario shown in figure 3*c*, where the switching on of the control force leads to complex quasi-periodic oscillations of the localized structure.

So far, the control and stabilization of the supercritical regime (Re(*a*_{1})>0, Re(*a*_{2})<0 for *α*=0) has been discussed. However, the stationary solution (4.11) and (4.12) can also be found for a proper choice of *α* and *τ* even if Re(*a*_{2})>0 for *α*=0. This situation is presented in figure 5, where a direct numerical simulation of equation (2.1) is performed for parameters in the subcritial regime (see also figure 1*b*). Here, a single localized structure in the absence of nonlinear stabilization is prevented from collapse, which leads either to a stabilization of the oscillations (figure 5*a*) or to a breathing with a constant amplitude (figure 5*b*). As was mentioned above, the collapse of the localized structure in the absence of control is caused by the delayed inhibition, that is, both inhibitors are too slow to suppress a strong increase in the activator concentration in the course of time. In this framework, the time-delayed feedback control can be seen as an additional delayed inhibition, forcing slow inhibitors to be fast enough to follow the activator distribution. Using the concept of additional induced delayed inhibition, all results obtained for the supercritical case can be explained in the same way: the time-delayed feedback control provides a mechanism of additional delayed inhibition, which can lead either to the stabilization of the solution if controlled inhibitors are fast enough, or to oscillatory dynamics otherwise.

## 5. Conclusion

In this paper, I have discussed the dynamics of breathing localized structures in a three-component reaction–diffusion system under time-delayed feedback control. I have analytically studied the linear stability problem of the delayed system, yielding an explicit expression for the boundary of stabilization domains within which a breathing localized structure can be effectively stabilized. However, I have found that outside of these domains more complex delay-induced periodic or quasi-periodic oscillations can be obtained. In order to understand the influence of the delayed feedback term on the behaviour of the breathing localized structure, an order parameter equation for the amplitude of the breathing mode was derived. Information about the dynamics of the system is contained in the complex coefficients of this equation. Using this information, the dependence of the oscillation radius on delay parameters can be explicitly derived, providing a robust mechanism to control the behaviour of the breathing localized structure in a straightforward manner. Note that the order parameter equation obtained is different from the normal form of the (subcritical) Hopf bifurcation subjected to time-delayed feedback, which contains a complex phase ahead of the delayed term [39,40]. However, it falls into place if one recollects that the reaction–diffusion system (2.1) is a real-valued system with real-valued delay rate. That is, one may obtain the normal form, mentioned above, by considering complex-valued model systems (e.g. a complex Ginzburg–Landau equation) with complex delay rates. To conclude, note that, although all results are obtained for the identity coupling term, the form of the order parameter equation still stands for other coupling matrices. However, explicit calculation of the corresponding coefficients in this case is intricate and will be treated elsewhere. Moreover, all results are derived in general form and can be applied to a wide class of spatial extended non-variational systems admitting localized structures, such as control of breathing and spatio-temporal chaotic localized states discussed in the framework of the non-variational Swift–Hohenberg equation and also observed in liquid crystal light valve experiments [41,42].

## Footnotes

One contribution of 19 to a Theme Issue ‘Localized structures in dissipative media: from optics to plant ecology’.

- © 2014 The Author(s) Published by the Royal Society. All rights reserved.