Control of traveling localized spots

TL;DR

This paper develops an analytical open-loop control method to guide stable traveling localized spots in reaction-diffusion systems along desired paths without altering their shape, confirmed by numerical optimal control solutions.

Contribution

It provides an explicit analytical control strategy based on Goldstone modes for guiding spots, bypassing complex reaction kinetics details and enhancing computational efficiency.

Findings

Control signal expressed analytically in terms of Goldstone modes.

Numerical confirmation of the control's optimality.

Analytical expressions serve as effective initial guesses for optimal control.

Abstract

Traveling localized spots represent an important class of self-organized two-dimensional patterns in reaction-diffusion systems. We study open-loop control intended to guide a stable spot along a desired trajectory with desired velocity. Simultaneously, the spot's concentration profile does not change under control. For a given protocol of motion, we first express the control signal analytically in terms of the Goldstone modes and the propagation velocity of the uncontrolled spot. Thus, detailed information about the underlying nonlinear reaction kinetics is unnecessary. Then, we confirm the optimality of this solution by demonstrating numerically its equivalence to the solution of a regularized, optimal control problem. To solve the latter, the analytical expressions for the control are excellent initial guesses speeding-up substantially the otherwise time-consuming calculations.

Tables

Table 1. Table 1: Parameter values used in the numerical simulations. The parameters θ = 1 , κ 2 = 2 , formulae-sequence 𝜃 1 subscript 𝜅 2 2 \theta=1,\,\kappa_{2}=2, κ 3 = 1 , subscript 𝜅 3 1 \kappa_{3}=1, and κ 4 = 8.5 subscript 𝜅 4 8.5 \kappa_{4}=8.5 are the same for set 1 and set 2.

	$D_u$	$D_v$	$D_w$	$\tau$	${\kappa }_1$	$v_0^x$
set $1$	$1.0\cdot {10}^{-4}$	$1.86\cdot {10}^{-4}$	$9.6\cdot {10}^{-3}$	$48.0$	-$6.92$	$2.599\cdot {10}^{-3}$
set $2$ [44]	$0.9\cdot {10}^{-4}$	$1.00\cdot {10}^{-3}$	$1.0\cdot {10}^{-2}$	$40.0$	-$7.30$	$1.776\cdot {10}^{-3}$

Equations

$${\partial }_{t}u(\mathbf{r},t)=$$

$$\tau {\partial }_{t}v(\mathbf{r},t)=$$

$$\theta {\partial }_{t}w(\mathbf{r},t)=$$

$${\partial }_t\mathbf{U}(\mathbf{r},t)-\mathbb{D}\Delta \mathbf{U}(\mathbf{r},t)-\mathbf{R}(\mathbf{U}(\mathbf{r},t))=\mathbb{B}\mathbf{f}(\mathbf{r},t).$$

$$\mathbf{U}(\mathbf{r},t_0)={\mathbf{U}}_0(\mathbf{r}),$$

$$\mathbf{U}(x_a,y,t)=$$

$$\mathbf{U}(x,y_a,t)=$$

$$\mathbf{f}(\mathbf{r},t)={\mathbb{B}}^{-1}\{{\partial }_t{\mathbf{U}}_d(\mathbf{r},t)-\mathbb{D}\Delta {\mathbf{U}}_d(\mathbf{r},t)-\mathbf{R}({\mathbf{U}}_d(\mathbf{r},t))\}.$$

$$0=$$

$$\mathbf{U}_{d}(\mathbf{r},t)=\mathbf{U}_{c}(\mathbb{A}(-\Phi_{\varphi}(t))\left(\mathbf{r}-\mbox{\boldmath$\Phi$}(t)\right)).$$

$$\mathbf{f}(\mathbf{r},t)=$$

$$\displaystyle\mathbf{f}_{\mathrm{Gold}}(\mathbf{r},t)=\,\left[\left(\begin{pmatrix}v_{0}^{x}\\ v_{0}^{y}\\ 0\end{pmatrix}-\mathbb{A}_{z}\left(-\Phi_{\varphi}(t)\right)\begin{pmatrix}\dot{\Phi}_{x}(t)\\ \dot{\Phi}_{y}(t)\\ \dot{\Phi}_{\varphi}(t)\end{pmatrix}\right)\cdot\tilde{\nabla}_{\mbox{\boldmath$\xi$}}\right]\mathbf{U}_{c}(\mbox{\boldmath$\xi$})\Bigg{|}_{\mbox{\boldmath$\xi$}=\bar{\mbox{\boldmath$\xi$}}(t)}$$

$$\displaystyle\mathcal{L}=\,\mathbb{D}\Delta_{\mbox{\boldmath$\xi$}}+\mathbf{v}_{0}\cdot\nabla_{\mbox{\boldmath$\xi$}}+\mathcal{D}\mathbf{R}(\mathbf{U}_{c}\left(\mbox{\boldmath$\xi$}\right)),$$

$$J(\mathbf{U},\mathbf{f})=$$

$$\left[\underset{0}{\overset{T}{\int }}\underset{\Omega }{\int }\left(\left(\mathbf{S}({\mathbf{f}}_{\mathrm{opt}})\right)-{\mathbf{U}}_d\right)\cdot \left({\mathbf{S}}^{\prime }({\mathbf{f}}_{\mathrm{opt}})\mathbf{h}\right)\text{rmd}\mathbf{r}\text{rmd}t+\nu \underset{0}{\overset{T}{\int }}\underset{\Omega }{\int }{\mathbf{f}}_{\mathrm{opt}}\cdot \mathbf{h}\text{rmd}\mathbf{r}\text{rmd}t\right]=0\forall \mathbf{h}.$$

$$\underset{0}{\overset{T}{\int }}\underset{\Omega }{\int }\left(\mathbf{P}+\nu {\mathbf{f}}_{\mathrm{opt}}\right)\cdot \mathbf{h}\text{rmd}\mathbf{r}\text{rmd}t=0\forall \mathbf{h},$$

$$-{\partial }_t\mathbf{P}(\mathbf{r},t)-\mathbb{D}\Delta \mathbf{P}(\mathbf{r},t)-\mathcal{D}{\mathbf{R}}^T({\mathbf{U}}_{\mathrm{opt}}(\mathbf{r},t))\mathbf{P}(\mathbf{r},t)={\mathbf{U}}_{\mathrm{opt}}-{\mathbf{U}}_d\text{mbox}inQ,$$

$$\mathbf{P}+\nu {\mathbf{f}}_{\mathrm{opt}}=0.$$

$${\Phi }_x(t)=r\sin (4\pi t/T)),{\Phi }_y(t)=r\sin (6\pi t/T)),\text{mbox}and{\Phi }_{\phi }(t)=0,$$

$$\Vert h(t){\Vert }_{L^1(\Omega )}=\underset{\Omega }{\int }|h\left(\mathbf{r},t\right)|\text{rmd}\mathbf{r}.$$

$${\dot{\Phi }}_i(t)=$$

$${\Phi }_x(t)=r\sin \left({\Phi }_{\phi }(t)\right),{\Phi }_y(t)=-r\cos \left({\Phi }_{\phi }(t)\right),{\Phi }_{\phi }(t)=2\pi t/T.$$

$${\dot{\Phi }}_x(t)=v_0^x\cos \left({\Phi }_{\phi }(t)\right),{\dot{\Phi }}_y(t)=v_0^x\sin \left({\Phi }_{\phi }(t)\right).$$

$$\mathbf{f}_{\mathrm{Gold}}(\mathbf{r},t)=\,-\dot{\Phi}_{\varphi}(t)\partial_{\varphi}\mathbf{U}_{c}(\mathbb{A}(-\Phi_{\varphi}(t))\left(\mathbf{r}-\mbox{\boldmath$\Phi$}(t)\right)).$$

$${\Phi }_{\phi }(t)=\frac{\pi }{4}\sin \left(\frac{2\pi t}{T}\right),$$

$${\partial }_{t}u(\mathbf{r},t)=$$

$$\tau {\partial }_{t}v(\mathbf{r},t)=$$

$$\theta {\partial }_{t}w(\mathbf{r},t)=$$

$$\tilde{v}(\mathbf{r},t)={\mathcal{K}}_v^0\circ v_0+\frac{1}{\tau }{\mathcal{K}}_v\circ \left(u+f_v\right),\tilde{w}(\mathbf{r},t)={\mathcal{K}}_w^0\circ w_0+\frac{1}{\theta }{\mathcal{K}}_w\circ \left(u+f_w\right),$$

$${\partial }_{t}u(\mathbf{r},t)=$$

$$+f_u-\frac{{\kappa }_3}{\tau }{\mathcal{K}}_v\circ f_v-\frac{{\kappa }_4}{\theta }{\mathcal{K}}_w\circ f_w.$$

$$f_{u,\mathrm{Gold}}^{\mathrm{sing}}\left(\mathbf{r},t\right)=f_{u,\mathrm{Gold}}\left(\mathbf{r},\mathbf{t}\right)-\frac{{\kappa }_3}{\tau }{\mathcal{K}}_v\circ f_{v,\mathrm{Gold}}-\frac{{\kappa }_4}{\theta }{\mathcal{K}}_w\circ f_{w,\mathrm{Gold}},$$

Full text

††: New J. Phys.

Keywords: open-loop control, optimal control, reaction-diffusion system, spots, coherent structures

1 Introduction

Localized spots, sometimes referred to as auto-solitons [1], dissipative solitons [2], or bumps [3], are a subclass of traveling waves that spontaneously evolve in two-dimensional (2D) dissipative nonlinear systems far from thermodynamic equilibrium. In a co-moving reference frame, spots are stationary solutions to coupled nonlinear partial differential equations (PDE), such as reaction-diffusion (RD) or neural field equations, for example. The characteristic length and time scales of the spots, i.e., their wave profile, propagation velocity, etc., are selected by the experimental conditions or the parameters of the model.

Experimentally, localized spots have been observed as current filaments in gas-discharge [4], as bright intensity spots in nonlinear optics and laser physics [5], as moving localized regions of increased concentration in chemical reactions [6], or coverage of adsorbed species in heterogeneous catalysis [7]. Further examples include temperature spots in fixed-bed catalytic reactors [8], actin conformation in dictyostelium discoideum [9], neural activity in head-direction cells [10], vegetation patterns [11], and many others.

Although control of self-organized patterns attracted considerable attention over the last decades, compare [12, 13] and references therein, it is still a challenging problem in applied nonlinear science. Often, one distinguishes between open-loop, closed-loop, and optimal control.

Open-loop control is independent of the instantaneous state of the system. As a consequence, it is inherently susceptible to perturbations in the initial conditions as well as to parameter uncertainty. Thus, detailed knowledge of the system’s dynamics and in-depth stability analysis are pre-conditions for reliable open-loop schemes. Typical examples of open-loop control are space-time dependent external forcing [14, 15] or control by imposed geometric constraints [16, 17].

On the other hand, in closed-loop or feedback control, the controlled state is permanently monitored to adjust the control signal accordingly [18, 19, 20]. Particularly, time-delayed feedback can induce pattern forming instabilities in addition to the pattern to be controlled [21].

Optimal control reformulates control problems in terms of the minimization of a cost functional [22, 23]. The cost functional measures the distance in function space between a desired target state and the actual controlled state of the system. If a control signal is the unique solution to an optimal control problem, then no other control, be it open- or closed-loop, will be able to enforce a controlled state closer to the target state. Conditions for sufficiency and uniqueness of optimal control are discussed extensively in the mathematical literature [24, 25, 26]. Optimal control of self-organized patterns requires complete knowledge of the PDEs governing the system’s evolution in time and space. Numerical solutions to optimal control of PDEs often base on computationally expensive iterative algorithms restricted to relatively small spatial domains and short time intervals. Clearly, the convergence to the target state sensitively depends on an appropriate initial guess for the control signal.

For traveling wave patterns, a general control task is position control aimed at guiding the pattern according to a given protocol of motion (POM), i.e., moving it with desired velocity along a desired trajectory through a spatial domain. In some technical applications like catalytic reactors, it is necessary to avoid the collision of high-temperature spots with the reactor walls or their pinning at heterogeneities of the catalyst’s support to maintain operational safety [8]. Another example of open-loop position control is the enhancement of the $\mathrm{CO}_{2}$ production rate during the low-pressure catalytic oxidation of $\mathrm{CO}$ on $\mathrm{Pt}(110)$ single crystal surfaces by dragging of reaction pulses and fronts using a focused laser beam with speeds differing from their natural propagation velocity in the absence of control [27, 28]. In a photosensitive Belousov-Zhabotinsky (BZ) medium, periodic variations of the applied light intensity forces a spiral wave tip to describe a wide range of hypocycloidal and epicycloidal trajectories [29, 30]. In optical bistable media like dye-doped liquid crystals and Kerr cavities, interface dynamics can be controlled by spatially inhomogeneous forcing [31]. Position control of traveling wave patterns can be tackled by feedback control as well. For example, the spiral wave core in a photosensitive BZ medium was steered around obstacles using feedback signals obtained from wave activity measured at a point detector, from tangential crossing of wavefronts with detector lines, or a spatially extended control domain [14, 32, 33]. Two feedback loops were used to stabilize and guide unstable traveling wave segments along pre-given trajectories [34]. Furthermore, feedback-mediated control loops were employed to stabilize plane waves undergoing transversal instabilities [35].

Recently, we proposed an open-loop control that acts solely via the Goldstone modes of wave patterns [36] and provides analytical expressions for the amplitude of the control signal to be applied; it is coined Goldstone mode control. We demonstrated that this control is able to accelerate or decelerate 1D traveling front and pulse solutions to RD equations [36, 37, 38] without changing their spatial profile. The stability of the control loop with respect to small changes in the initial conditions was discussed in [39]. Goldstone mode control also applies to move the core of a spiral wave at desired velocity along a pre-given trajectory through a 2D spatial domain, or to shape iso-concentration lines of 2D traveling pulses [40]. Interestingly enough, the control turned out to be equivalent to the solution of an appropriately formulated optimal control problem [36, 41].

In this paper, we extend Goldstone mode control to spatially localized moving spots. We introduce a three-component RD model supporting stable traveling spot solutions in section 2 and derive analytical expressions for position and orientation control of traveling spots in the fully-actuated case in section 3.1. The corresponding optimal control problem with an objective functional involving a Tikhonov regularization term is formulated explicitly in section 3.2. Here, we discuss the relation between Goldstone mode control derived in section 3.1 and the solutions to the optimal control problem. In section 4, after a brief description of the numerical methods being used, we discuss examples for fully-actuated position and orientation control of spots in subsection 4.1 and 4.2, respectively, as well as for under-actuated position control by a single control signal 4.3. Finally, we conclude the results in section 5.

2 Three-component spot model

Throughout this work, we consider the following three-component RD system exhibiting immobile and traveling stable spot solutions in $2$D [2, 42, 43]

[TABLE]

Here, $\Delta=\partial_{x}^{2}+\partial_{y}^{2}$ represents the Laplacian in Cartesian coordinates, $\mathbf{r}$ is the position vector in the spatial domain $\Omega$, $\mathbf{r}=(x,y)^{T}\in\Omega\subset\mathbb{R}^{2}$, and $t$ indicates time. $D_{u},D_{v}$, and $D_{w}$ denote the diffusion coefficients of components $u,v$, and $w$ while $\tau$ and $\theta$ set the time scales for the $v$ and $w$ kinetics, respectively. Beside spots, the model (1) is capable to support peanut patterns [44], breathing solitons [45], and jumping oscillons [46], for example.

The three-component RD system (1) was first introduced by Purwins and co-workers to model the dynamics of current filaments in planar gas-discharge [2]. In this context, activator $u$ and inhibitor $v$ represent the current density and the voltage drop over the high-ohmic electrode, respectively. The second inhibitor $w$ is linked to the surface charge, and the additive bifurcation parameter $\kappa_{1}$ is related to the supply voltage. Replacing the constant additive parameter in (1) by a space-dependent quantity $\kappa_{1}(\mathbf{r})$ breaks the translation and the rotation Euclidean symmetries of the equations. The interaction of traveling spots with different types of parameter heterogeneities in $1$D and $2$D has been studied in detail by many authors, see [44] and references therein. Penetration, rebound, annihilation, oscillation, as well as stationary or oscillatory pinning of spots were observed. Figure 2 illustrates different outcomes of the interaction between a traveling spot and a localized circular defect formed by a finite jump $\delta k_{1}=\kappa_{1}^{\mathrm{het}}-\kappa_{1}^{\mathrm{back}}$ from a background value $\kappa_{1}^{\mathrm{back}}$ and a higher value $\kappa_{1}^{\mathrm{het}}$ inside the circular heterogeneity. Additionally, other scenarios of spot-defect interaction have been found like repeated creation of spots inside the heterogeneity as well as spots orbiting both inside and outside of the defect’s boundary. If the activator describes the temperature in a catalytic packed-bed reactor, resting hot spots [47] or those pinned to local heterogeneities can damage the catalyst support. In particular, collision of hot spots with the reactor walls must be prevented for safety reasons. Consequently, guidance of a traveling spot with given velocity along a desired trajectory through a bounded spatial domain might be particular challenge in chemical engineering applications.

3 Controlling position and orientation of traveling spots

3.1 Analytical expression for control amplitudes in position control

Let us consider a controlled RD system according to

[TABLE]

The space-time dependent control signals $\mathbf{f}(\mathbf{r},t)=(f_{1}(\mathbf{r},t),\dots,f_{m}(\mathbf{r},t))^{T}$, $m\in\mathbb{N}$, on the right hand side of (2a) are assumed to act for all times $t$ everywhere within $\Omega$.

The constant $n\times m$ matrix $\mathbb{B}$ determines which components are directly affected by the control signals. A system with strictly less independent control signals than components, $m<n$, is underactuated. For $m=n$ and $\mathbb{B}$ invertible, the system is fully actuated. In what follows, we focus on fully actuated systems and set $\mathbb{B}$ equal to the identity matrix $\mathds{1}$. The limiting case of single component control, i.e., $\mathbb{B}\mathbf{f}(\mathbf{r},t)\propto(f_{1}(\mathbf{r},t),0,\dots,0)^{T}$, we consider in subsection 4.3.

The partial differential equations (2a) describe the evolution of the components $\mathbf{U}(\mathbf{r},t)$ in the presence of spatio-temporal perturbations $\mathbf{f}(\mathbf{r},t)$ that break the translation and rotation invariance of the unperturbed equations. In this interpretation, the response of the unperturbed solution to a given small input $\mathbf{f}$ can be calculated perturbatively, see Ref. [48, 49, 50], and the SI.

In this paper, following [36], for given desired spot dynamics, we perceive (2a) as a conditional equation for the perturbations which now are considered as control inputs. The goal of the control $\mathbf{f}$ is to enforce a state $\mathbf{U}$ to follow a given desired distribution $\mathbf{U}_{d}(\mathbf{r},t)=\left(u_{1,d}(\mathbf{r},t),\ldots,u_{n,d}(\mathbf{r},t)\right)^{T}$ as closely as possible everywhere in the spatial domain $\Omega$ and for all times $0\leq t\leq T$. We call a desired distribution $\mathbf{U}_{d}$ exactly realizable if there exists a control $\mathbf{f}$ such that the controlled state $\mathbf{U}$ equals $\mathbf{U}_{d}$ everywhere in the space-time cylinder $Q=\Omega\times[0,T]$.

Inserting $\mathbf{U}_{d}$ for $\mathbf{U}$ in (2a) yields for the control

[TABLE]

For $\mathbf{U}_{d}$ to be exactly realizable, three more conditions must be satisfied: First, the initial condition for the controlled state, (2b), must coincide with the initial state of the desired distribution, $\mathbf{U}(\mathbf{r},0)=\mathbf{U}_{d}(\mathbf{r},0)$. Second, all boundary conditions for the desired distribution $\mathbf{U}_{d}$ have to comply with the boundary conditions for $\mathbf{U}$, (2c). Third, $\mathbf{U}_{d}$ must be sufficiently smooth in the space-time cylinder $Q=\Omega\times[0,T]$ such that the derivatives $\partial_{t}\mathbf{U}_{d}$ and $\Delta\mathbf{U}_{d}$ are continuous.

Next, we formulate the control goal for spot solutions to the uncontrolled RD equations(2). These solutions propagate with constant velocity $\mathbf{v}_{0}=(v_{0}^{x},v_{0}^{y})^{T}$ and wave profile $\mathbf{U}_{c}$ through the spatial domain. In a co-moving frame of reference, \mbox{\boldmath\xi}=(\xi_{x},\xi_{y})^{T}\equiv\mathbf{r}-\mathbf{v}_{0}t, $\mathbf{U}_{c}$ satisfies the equation

[TABLE]

where, \nabla_{\mbox{\boldmath\xi}}=(\partial_{\xi_{x}},\partial_{\xi_{y}})^{T} and \Delta_{\mbox{\boldmath\xi}}=\partial_{\xi_{x}}^{2}+\partial_{\xi_{y}}^{2} denote the component-wise gradient and Laplacian, respectively. We emphasize that resting localized spots, $\mathbf{v}_{0}=\mathbf{0}$, are rotationally-symmetric solutions while traveling localized spots are axis-symmetric with the symmetry axis directed tangentially to the trajectory of motion, cf. figure 1(a) and (b), respectively. We characterize the current position of a spot by the $x$- and $y$-coordinates of the maximum value of the activator concentration along its symmetry axis at a given time, \mbox{\boldmath\Phi}(t)=(\Phi_{x}(t),\Phi_{y}(t))^{T}, and its orientation by the angle $\Phi_{\varphi}(t)$ between the spot’s symmetry axis and the $x$-axis, compare figure 3.

A distribution following a prescribed POM \mbox{\boldmath\Xi}(t)=(\mbox{\boldmath\Phi}(t),\Phi_{\varphi}(t))^{T}, while simultaneously preserving the profile of the uncontrolled spot $\mathbf{U}_{c}$, reads

[TABLE]

Here, $\mathbb{A}(\alpha)=\left[\,\cos(\alpha),-\sin(\alpha);\,\sin(\alpha),\cos(\alpha)\,\right]$ is the clockwise rotation matrix in $2$D. For the desired distribution (5) to be exactly realizable, the initial condition must be a spot solution of the form $\mathbf{U}(\mathbf{r},t_{0})=\mathbf{U}_{c}(\mathbb{A}(-\phi_{0})\left(\mathbf{r}-\mathbf{r}_{0}\right))$, which yields for the initial values of the POM \mbox{\boldmath\Phi}(t_{0})=\mathbf{r}_{0} and $\Phi_{\varphi}(t_{0})=\phi_{0}$, respectively. Inserting the desired distribution (5) into the general control solution (3) leads to

[TABLE]

with \bar{\mbox{\boldmath\xi}}(t)=\mathbb{A}(-\Phi_{\varphi}(t))\left(\mathbf{r}-\mbox{\boldmath\Phi}(t)\right). For the sake of a compact notation, we introduced the differential operator \tilde{\nabla}_{\mbox{\boldmath\xi}}=(\partial_{\xi_{x}},\partial_{\xi_{y}},\partial_{\varphi})^{T} with the angular derivative $\partial_{\varphi}=-\xi_{y}\partial_{\xi_{x}}+\xi_{x}\partial_{\xi_{y}}$. The dot denotes the derivative with respect to time $t$, and $\mathbb{A}_{z}(\alpha)$ is the clockwise rotation matrix around the $z$-axes in $3$D, $\mathbb{A}_{z}(\alpha)=\mathrm{diag}(\mathbb{A}(\alpha),1)$. Using equation (4) for the uncontrolled spot profile, we end up with the expression

[TABLE]

for our Goldstone mode control.

Remarkably, any reference to the nonlinear functions $\mathbf{R}$ drops out from the result (7). This is of great advantage in all applications where the details of the underlying reaction kinetics $\mathbf{R}$ are largely unknown or difficult to identify. Once propagation velocity $\mathbf{v}_{0}$ and wave profile $\mathbf{U}_{c}$ of the uncontrolled spot are measured with an accuracy sufficient to calculate the Goldstone modes $\partial_{\xi_{x}}\mathbf{U}_{c}$, $\partial_{\xi_{y}}\mathbf{U}_{c}$, and $\partial_{\varphi}\mathbf{U}_{c}$, the control signals can be computed in advance for the complete time interval $[0,T]$. Consequently, in contrast to feedback control, a continuous recording of the system is not required.

One notices that (7) equals the sum of Goldstone modes with time-dependent prefactors, \mathbf{f}_{\mathrm{Gold}}(\mathbf{r},t)=\,P_{1}(t)\,\partial_{\xi_{x}}\mathbf{U}_{c}(\mbox{\boldmath\xi})+P_{2}(t)\,\partial_{\xi_{y}}\mathbf{U}_{c}(\mbox{\boldmath\xi})+P_{3}(t)\,\partial_{\varphi}\mathbf{U}_{c}(\mbox{\boldmath\xi}). The Goldstone modes are the right eigenvectors to the linear stability operator $\mathcal{L}$ of (4)

[TABLE]

to the eigenvalue zero. They are associated with the translational and rotational invariance of equation (2a) in $\mathbb{R}^{2}$ for $\mathbf{f}(\mathbf{r},t)=0$. Clearly, the prefactors’ magnitudes are proportional to the difference between the intrinsic velocity, $\mathbf{v}_{0}$, and the current prescribed spot velocity projected onto the $x$- and $y$-axes. If the prescribed POM \mbox{\boldmath\Xi}(t) coincides with the spot’s natural motion, then all prefactors vanish identically and $\mathbf{f}_{\mathrm{Gold}}$ disappears everywhere in $Q$. Importantly, the control signal is localized around the spot position and vanishes far from it because the spatial derivatives of its profile decay sufficiently fast, \lim_{\lVert\mbox{\boldmath\xi}\rVert\to\infty}\nabla_{\mbox{\boldmath\xi}}\mathbf{U}_{c}=\mathbf{0}.

Alongside with these advantages, limitations in the applicability of Goldstone mode control (7) exist as well. For instance, the magnitude of the applied control may locally attain values that are unfeasible to realize physically because $\mathbf{f}_{\mathrm{Gold}}$ is proportional to the slope of the controlled wave profile $\mathbf{U}_{c}$. The stability of the control scheme depends sensitively on how precise the Goldstone modes can be calculated. Further, the complete spatial domain $\Omega$ accessible by the spot has to be available for the control as well. Additionally, as already mentioned above, $\mathbf{f}_{\mathrm{Gold}}$ cannot be applied to desired trajectories $\mathbf{U}_{d}$ which do not comply with initial as well as boundary conditions or which are non-smooth. While all these cases cannot be treated within the analytical approach proposed here, optimal control can deal with many of these complications.

3.2 Optimal control

An optimal control minimizes a so-called objective functional $J$ defined as a non-negative tracking-type functional

[TABLE]

$\mathbf{U}$ satisfies the controlled state equation associated to $\mathbf{f}$ with respect to given initial and boundary condition, cf. equations (2). The first term appearing in $J$ measures the distance between the actual and the desired solution $\mathbf{U}$ and $\mathbf{U}_{d}$ up to the terminal time $T$ in an $L^{2}(Q)$-sense. In the second, so-called Tikhonov regularization term, a small but finite, positive value $\nu$ guarantees the existence of an optimal control $\mathbf{f}_{\mathrm{opt}}$ that minimizes the objective functional $J$ (9) for $\Omega\subset\mathbb{R}^{q},q=1,2,3$, see Ref. [51].

For exactly realizable desired states, $\mathbf{U}=\mathbf{U}_{d}$, the solution to (3) equals the solution to the unregularized optimal control problem for $\nu=0$. If $\mathbf{U}_{d}$ is not exactly realizable, the controlled state $\mathbf{U}$ must be obtained as part of the solution to the optimal control problem. The minimization of $J$ must be performed with respect to state $\mathbf{U}$ and control $\mathbf{f}$. Expressing $\mathbf{U}$ in terms of $\mathbf{S}(\mathbf{f})$, where $\mathbf{S}:\mathbf{f}\mapsto\mathbf{U}$ is the solution operator to (2) in $Q$, justifies the definition of a reduced objective functional $J(\mathbf{f}):=J(\mathbf{S}(\mathbf{f}),\mathbf{f})$. In order to minimize $J(\mathbf{f})$, its first directional derivative with respect to $\mathbf{f}$ has to equal zero in all directions $\mathbf{h}$; yielding the necessary optimality conditions

[TABLE]

The state $\mathbf{U}$ is constrained to satisfy the controlled state equation together with given initial and boundary conditions, cf. equations (2). Similar as in ordinary minimization problems, a constrained minimization can be transformed to an unconstrained one by introducing Lagrange multipliers $\mathbf{P}(\mathbf{r},t)=(p_{1}(\mathbf{r},t),\ldots,p_{n}(\mathbf{r},t))^{T}$, also called the adjoint state. By means of the latter, (10) can be reformulated

[TABLE]

whereby the adjoint state is the solution of the adjoint equation

[TABLE]

subject to terminal condition $\mathbf{P}(\cdot,T)=\mathbf{0}$ in $\Omega$ and periodic boundary conditions in $\partial\Omega$. Here, $\mathcal{D}\mathbf{R}^{T}$ denotes the transposed Jacobian matrix of $\mathbf{R}$ with respect to $\mathbf{U}$. It is rather obvious that the condition (11) is equivalent to the condition

[TABLE]

This is nothing more than the well-known condition that, in a minimum, the gradient of the function to be minimized is zero.

Due to the mixed initial and terminal conditions for $\mathbf{U}$ and $\mathbf{P}$ it is rarely possible to find numerical solutions to optimal control by a direct integration method. To reduce numerical costs, we employ Model Predictive Control and divide our optimal control problem in subproblems with a $4$ time-step small time-horizon [41]. Thereby, each subproblem is solved with a gradient-type method. Details on the iteration scheme are discussed in the supplementary information (SI), paragraph S1.

4 Examples

In the following, we discuss three examples for position control of traveling spot solutions to the three-component RD model (1). Mainly, we compare Goldstone mode control $\mathbf{f}_{\mathrm{Gold}}$ with optimal control $\mathbf{f}_{\mathrm{opt}}$. If not stated otherwise, the state equation (2a) and the adjoint equation (12) are solved on a squared domain $\Omega=(-0.5,0.5]\times(-0.5,0.5]$ with periodic boundary conditions (2c) both in $x$ and $y$. The domain size is sufficiently large to avoid self-interaction of the spots in the periodic simulation domain. Without loss of generality, we fix the spots’ direction of motion to coincide with the $x$-axis, i.e., $v_{0}^{x}\neq 0$ and $v_{0}^{y}=0$. Any numerical simulation of equation (2a) is initialized with the profile $\mathbf{U}_{c}$ of the uncontrolled spot. This profile and the corresponding natural velocity $\mathbf{v}_{0}$ are obtained by solving the nonlinear eigenvalue problem (4) with adequate accuracy. Further details on the used numerical methods, the spatial and temporal resolution, and the initial conditions are presented in the SI, S1.

4.1 Translational position control of spots

In our first example, we aim to shift the spot’s position along a Lissajous curve without controlling its orientation, i.e., the spot’s symmetry axis is kept frozen to the $x$-axis. Thus, the POM \mbox{\boldmath\Xi}(t)=\left(\Phi_{x}(t),\Phi_{y}(t),\Phi_{\varphi}(t)\right)^{T} is given by

[TABLE]

with radius $r=0.2$ and protocol duration $T=200$. The video [SI_video1] shows the complete dynamics of all three state components $\mathbf{U}$ as well as $\mathbf{f}_{\mathrm{Gold}}$ and $\mathbf{f}_{\mathrm{opt}}$. In figure 4(a), we depict the time evolution of the activator distribution $u(\mathbf{r},t)$ under the action of the control $\mathbf{f}_{\mathrm{Gold}}$ shown in figure 4(b). One observes that the spot follows the desired trajectory indicated by the dashed line. The orientation of its symmetry axis remains fixed while the control signal realigns at any instants of time. Indeed, in the absence of orientation control, $\Phi_{\varphi}(t)=\,0$, $\mathbf{f}_{\mathrm{Gold}}$ can be expressed by the projection of $\nabla\mathbf{U}_{c}$ onto the tangential vector to the Lissajous curve $\mathbf{T}(t)$, \mathbf{f}_{\mathrm{Gold}}(\mathbf{r},t)=\left[\left(v_{0}^{x}-\dot{\Phi}_{x}(t)\right)\partial_{\xi_{x}}\mathbf{U}_{c}\left(\mbox{\boldmath\xi}(t)\right)-\dot{\Phi}_{y}(t)\partial_{\xi_{y}}\mathbf{U}_{c}\left(\mbox{\boldmath\xi}(t)\right)\right]\propto\mathbf{T}(t)\cdot\nabla\,\mathbf{U}_{c}, with \mbox{\boldmath\xi}(t)=\mathbf{r}-\mbox{\boldmath\Phi}(t). Obviously, the control is localized at the current spot position \mbox{\boldmath\Phi}(t) and vanishes far away from it. Despite that the average speed $\bar{v}=L_{\mathrm{curve}}/T\approx 6v_{0}^{x}$ along the studied Lissajous curve (14) with arc length $L_{\mathrm{curve}}$ is almost five times larger than the propagation velocity of the uncontrolled spot, the magnitude of $f_{u,\mathrm{Gold}}$ is of the same order as the local reaction terms (1). The control signals applied to the inhibitors $v$ and $w$ are one and two magnitudes smaller [SI_video1] than the activator’s control, respectively.

On the scale of [SI_video1], there is no distinguishable difference between $\mathbf{f}_{\mathrm{Gold}}$ and $\mathbf{f}_{\mathrm{opt}}$. Both are always localized close to the current spot position, and their magnitudes change proportional to |\dot{\mbox{\boldmath\Phi}}(t)|. For a quantitative comparison, we compute the relative errors between $f_{u,\mathrm{Gold}}$ and $f_{u,\mathrm{opt}}$ measured by the $L^{1}(\Omega)$-norm

[TABLE]

Here, $|h\left(\mathbf{r},t\right)|$ indicates the absolute value of $h$ at position $\mathbf{r}$ and time $t$.

In figure 4(c), we depict solely the normalized error for the first half of the protocol because it starts to repeat after $T/2$, $\Phi_{y}(t)=-\Phi_{y}(t+T/2)$. The relative error between $f_{u,\mathrm{Gold}}$ and optimal control $f_{u,\mathrm{opt}}$ (solid line) is satisfactory and ranges between $2\%$ and $8\%$. As reported in S1 of the SI, the limiting error is dominated by the time step chosen in the implicit Euler-scheme. Albeit the scheme is A-stable, the error at a specific time $t$ is of the order of $\mathcal{O}(dt)$. Consequently, we observe that $\lVert f_{u,\mathrm{Gold}}(t)-f_{u,\mathrm{opt}}(t)\rVert_{L^{1}(\Omega)}$ is bounded from above by $dt$; $dt=0.1$ in the studied example.

The dashed line in figure 4(c) shows the relative error between the activator distribution obtained by Goldstone mode control and the one calculated under optimal control, $\lVert u_{\mathrm{Gold}}(t)-u_{\mathrm{opt}}(t)\rVert_{L^{1}(\Omega)}/\lVert u_{\mathrm{Gold}}(t)\rVert_{L^{1}(\Omega)}$. At any time, this error is less than $10^{-3}$, i.e., both controlled states agree remarkably well, despite that $\lVert f_{u,\mathrm{Gold}}(t)-f_{u,\mathrm{opt}}(t)\rVert_{L^{1}(\Omega)}$ is of the order $10^{-1}$. Additionally, the relative errors between the desired distribution $\mathbf{U}_{d}$ and the state solutions $\mathbf{U}_{\mathrm{Gold}}$ and $\mathbf{U}_{\mathrm{opt}}$ are less than $10^{-7}$ in both cases (not shown explicitly). This confirms that the Goldstone mode control (7), within numerical accuracy, is indeed the solution to the corresponding unregularized optimal control problem. Similar conclusions had been obtained in our previous study of position control of front solutions in one spatial dimension, see [41].

The gradient-type method, used to solve the optimal control problem, relies on an initial guess for the control signal. The closer the starting guess is to the final solution, the fewer iteration steps are necessary to converge for most established optimization methods. Starting every iteration with an initial zero control, it takes on average $\bar{n}_{\mathrm{iter}}\simeq 23$ iterations per time step for position control along the Lissajous curve (14). Using the control solution of the previously solved subproblem as initial guess reduces the average number of iterations to $\bar{n}_{\mathrm{iter}}\simeq 14$. Taking advantage of the similarity between $\mathbf{f}_{\mathrm{Gold}}$ and $\mathbf{f}_{\mathrm{opt}}$, see figure 4(c), the computational costs reduce even further. The most substantial computational speed-up is obtained by initiating every optimization subproblem with (3). Then, the iteration stops on average after the first step, $\bar{n}_{\mathrm{iter}}\simeq 1$.

4.1.1 Stability of position control

Any open-loop controls is sensitive against perturbations of the initial conditions, data uncertainty, or numerical roundoff errors. To test the stability of our Goldstone mode control for position control $\mathbf{f}_{\mathrm{Gold}}$, we accelerate or decelerate a single spot from its initial, intrinsic velocity $\mathbf{v}_{0}$ to a final velocity $\mathbf{v}_{1}$ using a translational POM \mbox{\boldmath\Xi}(t)=\left(\mbox{\boldmath\Phi},0\right)^{T} with velocity

[TABLE]

for $i\in\{x,y\}$. Note that both the protocol’s velocity \dot{\mbox{\boldmath\Phi}}(t) and acceleration \ddot{\mbox{\boldmath\Phi}}(t) are continuous functions within the interval $[0,T_{i}]$. $T_{i}$ denotes the duration of the protocol. The maximum acceleration $\pi\left(v_{i}^{1}-v_{i}^{0}\right)/(2T_{i})$ is proportional to the prescribed velocity difference $v_{i}^{1}-v_{i}^{0}$ and inversely proportional to $T_{i}$.

A sketch of the protocol is depicted in figure 5(b). Since the proposed control scheme is an open-loop control, deviations between the current spot position \mbox{\boldmath\Phi}_{\mathrm{curr}}(t) and the POM \mbox{\boldmath\Phi}(t) will grow unbounded in time if the difference between them exceeds a critical value [39]. A specific protocol is called stable and marked by green boxes in figure 5 if and only if the Euclidean distance is bounded as \lVert\mbox{\boldmath\Phi}_{\mathrm{curr}}(t)-\mbox{\boldmath\Phi}(t)\rVert<L/2 for all times $t\in[0,t_{\mathrm{end}}]$. Otherwise, it is called unstable (red boxes). Note that a protocol is also considered to be unstable if the control leads to the nucleation of additional spots. In order to make the results comparable for different protocol durations, we adjust the terminal simulation time $t_{\mathrm{end}}$ according to $t_{\mathrm{end}}=\mathrm{max}\left(10\,t_{\mathrm{drift}},T_{i}+10\,L/|v_{1}^{i}|\right)$ with drift time $t_{\mathrm{drift}}=L/v_{0}^{x}$. We stress that all simulation results presented in figure 5 have been computed for sufficiently long time intervals and do not alter upon an increase of the total simulation time.

Figure 5(a) depicts the numerically evaluated region of stable position control (green boxes) in $x$-direction as a function of the ratio of terminal spot velocity $v_{1}^{x}$ to the initial one $v_{0}^{x}$ and the ratio of the control duration $T_{x}$ to the drift time $t_{\mathrm{drift}}$. The translational POM in $y$-direction is set to zero, $\Phi_{y}(t)=0$. As expected, the numerical algorithm is stable in the absence of control, $v_{1}^{x}/v_{0}^{x}=1$. Further, it turns out that the control scheme is mostly stable for rapid, $T_{x}\ll t_{\mathrm{drift}}$, to moderately slow POMs, $T_{x}\lesssim 10\,t_{\mathrm{drift}}$, regardless of the velocity change, $|v_{1}^{x}-v_{0}^{x}|$. The stability regions exhibit an asymmetry with respect to the sign of the velocity change. Weakly accelerating protocols, $1<v_{1}^{x}/v_{0}^{x}\lesssim 2$, are unstable (red colored region) while decelerating ones, $v_{1}^{x}<v_{0}^{x}$, are always stable for $T_{x}\lesssim 10\,t_{\mathrm{drift}}$. This finding is in agreement with [39]. The instability for $v_{0}^{x}<v_{1}^{x}\lesssim 2v_{0}^{x}$ is caused by an undesired rotation of the spot induced by numerical truncation errors. These accumulate during the simulation and eventually result in an asymmetric perturbation (with respect to $y$) acting on the spot pattern. Once the spot starts to rotate and eventually drifts away from the centerline $y=0$, the proposed open-loop control $\mathbf{f}_{\mathrm{Gold}}$ can neither respond nor correct the undesired rotation. The impact of the numerical truncation error becomes more pronounced with growing protocol’s duration $T_{x}$ and results in a broad unstable region for long protocols, $T_{x}>10\,t_{\mathrm{drift}}$.

The situation changes if one aims to move the spot pattern perpendicular to its intrinsic direction of propagation, here in $y$-direction. In figure 5(c), we keep the motion in $x$ unchanged, $\Phi_{x}(t)=v_{0}^{x}\,t$, and accelerate the spot according to (16) along the $y$-direction. Because the controlled spot solution is symmetric with respect to the centerline $y=0$, position control in $y$ might be inherently unstable [39]. One notices immediately that regions with unstable position control are much larger compared to figure 5(a). Nevertheless, the control is stable for weak acceleration, $v_{1}^{y}\simeq 0.1v_{0}^{x}$, independent of the protocol’s duration. Increasing the terminal velocity further, Goldstone mode control starts to fail. Once a certain deviation between the current spots’ position and the proposed POM is attained, the pattern cannot follow the applied control anymore and starts to move freely. With further growing terminal velocity $v_{1}^{y}$, the control’s magnitude increases as well and thus successful position control can be re-stabilized. Longer protocols $T_{y}$ result in an accumulating of numerical truncation errors.

4.1.2 Orientation control with speed adjustment

In the previous paragraph, we’ve demonstrated that the stability of position control can be enhanced if in any current position of the spot its symmetry axis, given by $\Phi_{\varphi}(t)$, points tangentially to the direction of motion. Therefore, in our next example, we propose to shift the spot pattern along a circular trajectory by simultaneously controlling its orientation

[TABLE]

Here, $r$ denotes the radius of the circle and $T$ the protocol’s duration. For experimental realization compare [28], for example.

In figure 6, we present the temporal evolution of the activator distribution $u$ (a) controlled by $\mathbf{f}_{\mathrm{Gold}}$ (b). In line with the POM, the spot always keeps its symmetry axis at the tangent to the desired trajectory of motion. The control $\mathbf{f}_{\mathrm{Gold}}$ remains localized and is dominated by the translational Goldstone mode $\partial_{\xi_{x}}\mathbf{U}_{c}$ due to the acceleration along \mbox{\boldmath\Xi}(t) (17); the average speed is $\bar{v}\simeq 2.4v_{0}^{x}$. Notably, the maximum value of the control magnitude is half as strong compared to position control without adjusting the orientation, $\Phi_{\varphi}(t)=0$ (not explicitly shown). In panel (c), the temporal behavior of the relative error $\lVert f_{u,\mathrm{Gold}}(t)-f_{u,\mathrm{opt}}(t)\rVert_{L^{1}(\Omega)}/\lVert f_{u,\mathrm{Gold}}(t)\rVert_{L^{1}(\Omega)}$ measured by the $L^{1}\left(\Omega\right)$ norm is shown (solid line). They are large compared to pure position control along a Lissajous curve, cf. figure 4. Stronger deviations are caused by interpolation errors arising during numerical rotation of spot patterns by $\Phi_{\varphi}(t)$. The relative error attains a maximum at $\Phi_{\varphi}(t)=m\,45^{\circ},\,m\,\mbox{odd}$. At these angles, the distance between the nodes of the rotated grid and the underlying one is the largest, viz., $dx/\sqrt{2}$, and, hence, numerical interpolation errors become significant. Contrarily, the relative error minimizes at $\Phi_{\varphi}(t)=m\,90^{\circ},\,m\in\mathbb{Z}$, at which both grids coincide. Remarkable, the normalized error $\lVert u_{\mathrm{Gold}}(t)-u_{\mathrm{opt}}(t)\rVert_{L^{1}(\Omega)}/\lVert u_{\mathrm{Gold}}(t)\rVert_{L^{1}(\Omega)}$ (dashed line), is still less than $10^{-3}$ at any instants of time despite that the deviation of the associated controls rises up to $\sim 25\%$.

4.2 Orientation control

If the uncontrolled spot propagates at non-zero velocity \mathbf{v}_{0}\neq\mbox{\boldmath0}, the simplest way to navigate it through a spatial domain is to control exclusively its current orientation $\Phi_{\varphi}(t)$. If so, the translational components of the POM \mbox{\boldmath\Xi}(t) are determined by

[TABLE]

Clearly, one loses the possibility to control separately the $x$- and $y$-position of the pattern by limiting the speed to \lVert\dot{\mbox{\boldmath\Phi}}(t)\rVert=v_{0}^{x}. Inserting (18) into (7), the translational Goldstone modes drop out and we obtain

[TABLE]

Now, we pick up the problem formulated in section 2, namely, how to prevent pinning of a spot at a local heterogeneity in the domain. The heterogeneity is viewed as circular region where the parameter $\kappa_{1}$ jumps from a background value $\kappa_{1}^{\mathrm{back}},\,\forall\mathbf{r}\notin\Omega_{\circ}$ to a defect value $\kappa_{1}^{\mathrm{het}},\,\forall\mathbf{r}\in\Omega_{\circ}$ whereby $\Omega_{\circ}=\{\left(x,y\right)\in\mathbb{R}^{2}:(x+R)^{2}+y^{2}<R^{2}\}$ with radius $R=0.05$. The orientational POM for avoiding the heterogeneity is set to

[TABLE]

with duration $T=L_{x}/v_{0}^{x}$. Note that the corresponding prescribed positions $\left(\Phi_{x}(t),\Phi_{y}(t)\right)^{T}$ have to be calculated numerically.

In figure 7, we present the temporal evolution of the activator distribution $u$ in panel (a) and the corresponding control signal $f_{u,\mathrm{Gold}}$ in panel (b). The prescribed translational POM is indicated by the dashed lines. At first glance, the control signal possesses a more complicated shape and its magnitude is significantly reduced, $|f_{u,\mathrm{Gold}}|\lesssim 10^{-3}$, as compared to $|f_{u,\mathrm{Gold}}|\lesssim 10^{0}$ and $|f_{u,\mathrm{Gold}}|\lesssim 10^{-1}$ in the previous examples, cf. figure 4 and figure 6. Thus, orientation control is less invasive than position control. In return, we lose the ability for fast intervention into spot dynamics as well as for mayor increase in the speed of the spot. Additionally, orientation control is much more susceptible to fail. The small control magnitudes are too weak to suppress the impact of numerical round-off errors which may result in undesired spot rotation, cf. section 4.1.1. Caused by the small propagation velocity, the duration $T$ of the POM grows as compared to position control, see section 4.1, and therefore the probability of failure increases as well.

4.3 Position control by a single control signal

So far, we have discussed examples of fully actuated systems for which the number of state components equals the number of independent control signals. If the coupling matrix $\mathbb{B}$ is not invertible, expression (7) for $\mathbf{f}_{\mathrm{Gold}}$ cannot be used. The question arises how to extend our approach to underactuated systems [36, 52]. In the following example we assume a control acting on the activator $u$ only while inhibitors $v$ and $w$ remain uncontrolled, i.e., $f_{v,\mathrm{Gold}}(\mathbf{r},t)=f_{w,\mathrm{Gold}}(\mathbf{r},t)=0$. Control via an inhibitor has been discussed in detail for the Hodgkin-Huxley model and the three-component Oregonator model for photosensitive BZ reaction, compare supplemental information to [36].

To derive an expression for $f_{u,\mathrm{Gold}}(\mathbf{r},t)$, we start with the fully actuated system

[TABLE]

Equations (21b)-(21c) are linear, inhomogeneous PDEs with initial conditions $v(\mathbf{r},t_{0})=v_{0}(\mathbf{r})$ and $w(\mathbf{r},t_{0})=w_{0}(\mathbf{r})$, respectively. Their solutions can be written as

[TABLE]

where $\mathcal{K}_{i}^{0}$ and $\mathcal{K}_{i}$, $i\in\{v,w\}$, are integral operators involving Green’s functions to the homogeneous equations corresponding to (21b)-(21c) with associated initial conditions and to the inhomogeneous equations with zero initial conditions. Plugging (22) into (21a) gives

[TABLE]

From the last line of (4.3) we identify the expression for $f_{u,\mathrm{Gold}}(\mathbf{r},t)$ to be

[TABLE]

whereby the component of $\mathbf{f}_{\mathrm{Gold}}$ are determined by (7).

As an example for position control by a single control signal, we guide a spot along the Lissajous curve given by (14) with radius $r=0.2$ and protocol duration $T=200$. The spot’s orientation $\Phi_{\varphi}(t)=0$ remains uncontrolled. The relative errors between desired and controlled states are shown in figure 8(c). All states are obtained from numerical simulation of (1)-(2a) with control $\mathbf{f}_{\mathrm{Gold}}^{\mathrm{sing}}=(f_{u,\mathrm{Gold}}^{\mathrm{sing}},0,0)^{T}$ given by (24). The relative error for the activator $u$ (solid line) is less than $10^{-3}$ at any time $t$ and thus the controlled activator pattern agrees satisfactorily well with the desired distribution. This finding is corroborated by snapshots of $u$ at different instants of time in figure 8(a). In contrast to the activator, the profile of the inhibitor $v$ is not preserved under control but deformed considerably, see figure 8(b). In particular, an elongated region of activity becomes apparent along the Lissajous curve due to time scale separation in the RDS (1). The concentration of the slow inhibitor $v$, produced in the wake of the activator, decays exponentially to the rest state on a time scale $\tau=48\approx T/4$. Consequently, the relative error of $v$ (dashed line) attains relatively large values of the order $10^{-2}$. On the other hand, the fast inhibitor $w$ and the activator $u$ vary on the same characteristic time scale as $\theta=1$ was chosen in the considered example. Thus, we expect only small changes in both profiles in the presence of the control. In fact, the values of the relative error for $w$ turn out to be less than $10^{-4}$ which is even one magnitude smaller than the relative error of $u$, cf. the dash-dotted line in figure 8(c).

5 Conclusion

Localized traveling patterns are ubiquitous in spatially extended nonlinear systems driven far from thermodynamic equilibrium. These structures are often coined dissipative solitons or shortly spots and have been observed in various chemical, chemo-mechanical, electrical or neural systems. Hence, to control the position, orientation and velocity of a traveling spot is a key challenge.

Exploiting the translational and rotational symmetries of the governing equations, we’ve demonstrated that the control signal, which one has to apply to solve these tasks, is constituted by the Goldstone modes with time dependent prefactors. Intriguingly, for the latter analytic expressions have been derived. To deduce the control signal – coined Goldstone mode control – for a given protocol of motion, it is adequate to measure the spot’s profile and the corresponding propagation velocity with sufficient accuracy. In particular, Goldstone mode control is realized by external spatio-temporal forcing, i.e., it is an open-loop control. Contrary to closed-loop or feedback control, continuous monitoring of the system is not required. On the downside, as any open-loop control, the method is sensitive to perturbations. Therefore, the range of applicability has been checked by a stability analysis. Most importantly, our approach requires no detailed knowledge about the underlying reaction kinetics as opposed to standard open-loop control. Although the control signal is invasive, it is designed to preserve simultaneously the shape of the controlled pattern.

Remarkably, in all examples considered so far, Goldstone mode control is, within numerical accuracy, equal to solutions of an equivalent, non-regularized optimal control problem. Consequently, our control turns out to be optimal, i.e., no other control enforces the system closer to the desired target state according to the protocol of motion. Furthermore, these control signals have been proven to be excellent initial conditions for regularized optimal control problems; achieving a substantial computational speed-up. Generally, Goldstone mode control approach might serve as consistency check for numerical optimal control algorithms as well. We emphasize that optimal control is not only computationally demanding but requires full knowledge of the nonlinear kinetics. On the other hand, the scope of optimal control can be extended in many ways like sparse control or for inequality conditions for the control amplitudes’ upper and lower bounds [24, 41].

Due to the underlying symmetry considerations, Goldstone mode control is widely applicable. Already, the method have been successfully used to guide traveling interfaces and excitation pulses in $1$D [36, 41] and spiral waves [41] as well as to shape iso-concentration lines of traveling wave patterns [40] in $2$D. Recently, we successfully applied Goldstone mode control to spot solutions of neural field equations [53] that phenomenologically describe the dynamics of synaptically coupled neurons [54].

We thank Alexander Ziepke for helpful discussions as well as for critical reading of the manuscript and acknowledge financial support from the German Science Foundation DFG through the SFB 910 “Control of Self-Organizing Nonlinear Systems”.