Turing Patterns in Forced Open Two‐Side‐Fed‐Reactor

The mechanism suggested by Turing for reaction‐diffusion systems is widely used to explain pattern formation in biology and in many other areas. The persistence of patterns in altering environments is an important property in many natural cases. The experimental study of these phenomena can be done in chemical systems using appropriately designed reactors, e.g., in two‐side‐fed open gel reactors. This configuration allows for testing the effect of time‐periodic boundary conditions that generate periodic feeding of chemicals on the dynamics of Turing patterns. The numerical approach is based on a chemically realistic mechanism and a 2D description of the reactor that reproduces the feeding from the boundaries and the corresponding concentration gradients. Depending on the amplitude and the frequency of the forcing, two basic regimes are observed, spatiotemporal oscillations and pulsating spot pattern. In between them, a mixed‐mode pattern can also develop. Spot patterns can survive large amplitude forcing. The dynamics of the spot pulsation are analyzed in detail, considering the effect of the tanks and the chemical gradients that localize the patterns. These findings suggest that periodic feeding effectively controls pattern formation in chemical systems.


Introduction
Understanding the origin, the symmetries, and the structure of natural patterns is at the forefront of nonlinear science. The mutual interaction of reactions and transport processes, e.g., diffusion, is a simple mechanism to create such patterns.

DOI: 10.1002/adts.202300091
Due to their tractability and versatility, chemical systems with nonlinear feedback can be conveniently used in the experimental investigation of reaction-diffusion patterns. [1,2] The appearance of positive feedback, e.g., autocatalysis, is relatively common in inorganic redox reactions. Although the detailed mechanism of these reactions is often quite complex, the numerical investigation of the experimentally observed pattern formation can be successfully made by simplified models. Among the different spatiotemporal patterns observed in reaction-diffusion systems, the stationary ones deserve special attention according to their relevance in developmental biology, as it was revealed by Alan Turing. [3][4][5] These stationary patterns develop in the presence of short-range activation and long-range inhibition and have an intrinsic (chemical) wavelength. [6] A successful method to ensure this condition is based on using a low-mobility compound to bind the activator species, which results in an apparent decrease in its diffusivity. [4,5,7] For example, starch was used to bind triiodide ions in the first experimental observation of Turing patterns. Performing these reactions in conditions where pure reaction-diffusion patterns can develop also requires proper reactor design.
Open spatial reactors are convenient tools for the experimental study of sustained reaction-diffusion patterns. [8] The central part of these reactors is a porous material, e.g., a hydrogel, where the reaction-diffusion patterns can develop, avoiding any advection of the fluid. The continuous feed of the reactant necessary for sustained pattern formation is provided from the outer surface of the gel. The two most often used configurations are the oneside-feed (OSFR) and the two-side-feed (TSFR) open reactors. In the first case, the gel is connected through a single surface to the content of a continuous-stirred-tank reactor (CSTR). The CSTR contains all the chemicals but is kept in a state where the extent of the reaction is low, e.g., by operating at a short residence time. This way, the gel is fed by a mixture of fresh reactants. In a TSFR, two tanks are in contact with the gel part on two opposite sides. The reactants are distributed between the two tanks to avoid any chemical reactions. The reactions can occur only in the gel where the chemicals -diffused from the opposite sides -meet. This feeding mode creates counter gradients of the concentrations of the chemicals inside the gel, and the conditions of pattern formation meet in a localized region of the gel. The formation of localized patterns in the presence of cross-gradients has been studied both theoretically and experimentally in a configuration where two separate pools of reactants are connected. [9,10]  Turing patterns can be realized in different chemical reactions using OSFR and TSFR configurations. [7,[11][12][13][14][15] The typical experiments are performed at constant feeding compositions suitable for stationary pattern formation.
The effect of periodic forcing on Turing patterns is a challenging theoretical and experimental problem. Many experiments were carried out by using periodic illumination of the light-sensitive chlorine-dioxide-iodine-malonic acid (CDIMA) in OSFR, [16][17][18] , where suppression of patterns, [16] formation of entrained and oscillating patterns [19] and development of superlattice Turing patterns were observed. [14] This method creates forcing on the local kinetics all along the gel part of the OSFR. Recently, another method, periodic feeding from the boundary, was applied in a TSFR in similar chlorite-iodide-malonic acid (CIMA) reaction. [20] The experiments showed that patterns may appear and persist under periodic environmental conditions. Theoretical studies on the time-periodic forcing of Turing patterns near a Turing-Hopf bifurcation point presented that forcing may enhance or suppress Turing patterns, in agreement with the experimental observations. [21] Buceta and Lindberg showed that global forced alternation may lead to pattern formation in reactiondiffusion systems, and Turing instability develops when switching is sufficiently rapid. [22] In their corresponding numerical simulations, they used a reaction-diffusion model with periodic variations both in the reaction and diffusion terms and observed stationary and localized oscillatory patterns. The study of timedependent boundary conditions, by using periodic temperature variation in the environment, show spot splitting phenomena. [23] Here, we present 2D numerical simulations to investigate the effect of periodic feed in TSFR (Figure 1). Creating periodic boundary feed conditions is more straightforward in a TSFR, where the primary chemical reaction only takes place in the gel and not in the feeding tank, like in an OSFR. The applied description of the TSFR allows for studying the effect of gradients caused by the feeding method. As a chemical model, we use an extended version of the Rábai model of pH oscillators (reaction (R1)-(R4)). [24] pH oscillators produce many patterns, including Turing patterns, in experiments performed both in OSFRs and TSFRs. [12,13,15] The main aspects of these experimental observations have been successfully simulated by the Rabai model.
In this model B denotes an oxidant (bromate, iodate, or hydrogen peroxide), A − is the unprotonated form of a weak acid (e.g., sulfite ion), C and S denote a second substrate (e.g., ferrocyanide) and low mobility hydrogen ion binding agent, and P and Q are products. In this model, reaction (R2) represents positive feedback, whereas reaction (R3) provides negative feedback. The corresponding rate equations are written as follows: The role of reaction (R4) is to incorporate the presence of a low mobility agent (S − ) that reversibly binds the activator (H + ). This chemical reaction is responsible for the differential diffusion necessary for Turing pattern formation. The forcing on spatiotemporal oscillation has already been studied by 1D numerical simulations with the Rabai model in a TSFR and showed forced bursting and synchronization. [25] Here, we explore the forcing of Turing patterns in a TSFR focusing on the damping effect of the tank, the effect of the amplitude and frequency of the forcing on patterns dynamics, and the role of gradients caused by the TSFR configuration.

Turing Patterns in the Absence of Forcing
Dúzs and Szalai have experimentally and numerically demonstrated the formation of spatiotemporal patterns in pH oscillators performed in a TSFR. [26] Here, we recall some relevant aspects of the dynamics. To avoid reaction (R2) between reactants B and HA, they are separated as B is fed into tank A and HA is fed into tank B (Figure 1).
The necessary differential diffusion for Turing instability is regulated by reaction (R4). In the absence of S − (s tot = 0), the model shows spatiotemporal oscillations at the applied conditions. The oscillations are localized by the gradients caused by the feeding and appear along a line parallel to the feeding surfaces ( Figure 2). By increasing s tot the oscillation disappear at s tot = 0.14. The amplitude and the frequency of the oscillations are finite at this bifurcation point (Figures S1 and S2, Supporting Information). Above this critical s tot value localized spot pattern form in the gel ( Figure 2). There is a region (0.06 < s tot < 0.15) where the stability of spatiotemporal oscillations and Turing patterns overlap.
Here, we focus on the dynamics of Turing patterns at s tot = 0.5, far from the Turing-Hopf bifurcation point, where oscillations are entirely suppressed even in a spatially homogeneous reactor, e.g., in a continuous stirred-tank reactor ( Figure S3, Supporting Information). At these conditions, the spot pattern is stable in the range of the input flow concentration of HA, that is 0.43 ≤ a h0 ≤ 0.46. Below a h0 = 0.43, there is no autocatalytic outbreak of H + in the gel, and above a h0 = 0.46, a stationary stripe stays in the gel in parallel to the feeding surfaces. When the forcing was applied on the input concentration of HA according to Equation (3), a h0 = 0.5 was used.

Effect of Forcing on the Tank
In TSFR experiments, a straightforward control parameter is the composition of the input flow. At constant input flow, in the absence of any reaction in the tank, and if the feedback of the chemical composition of the gel content is negligible, the composition in the tank is equal to that of the input flow. Therefore, boundary conditions at the tank/gel surface are directly set by the compositions in the input flow. That is not the case when forcing is applied to the compositions in the input flow.
In the absence of any reaction, a sinusoidal forced ideal continuous stirred tank reactor (CSTR) is described by the following equation: where , c 0 , A f , and f are the residence time, the concentration in the input flow, the amplitude and frequency of the perturbation, respectively. The residence time of the CSTR is the ratio of the volume of the reactor and the volumetric flow rate. The term describing the chemical concentration in the input flow is c 0 + A f sin ( f t). The solution of Equation (1) is: where C is a constant. [25] In Equation (2), the amplitude of the damped concentration oscillations in the tank is described by the first term. In our simulations, a sinusoidal forcing was applied to the input flow concentration of HA according to Equation (3).
where h st B = h B at A f = 0. A rel 0 indicates the variation of the activator concentration during the forcing in tank B compared to its value at no forcing.
In Figure 3a the simulated values of A rel 0 are plotted as a function of the forcing frequency ( f ) at increasing values of A f . The general characteristics of the damping are clearer in Figure 3b where the ratio of A 0 ∕A max 0 is plotted as a function f . At low frequencies the 1∕( 2 f + 1) and at high frequencies the f ∕( 2 f + 1) term of Equation (2) determines the damping. As the forcing frequency on the input feed concentration increases, the amplitude of the oscillations at the tank/gel surface significantly decreases.
It is a significant effect that must be considered when analyzing the forced dynamics of patterns in TSFR.

Dynamics of the Forced Gel Content
In TSFR experiments, the patterns are localized along a line (a plane in 3D) parallel to the feeding surfaces inside the gel. As we have discussed, when the forcing is applied to the input of one of the tanks, the amplitude of the forced oscillations of the concentration at the tank/gel surface (the boundary conditions) strongly depends on f . The spread of this harmonic concentration oscillations at the boundary inside the gel content, in the absence of reactions, can be described as a 1D diffusion problem: with c = 0 ∀x at t = 0 and c = A f cos ( f t − ) at x = 0. The solution to this problem is the following: where k = √ f ∕(2D). [27] Equation (6) describes a wave whose amplitude decreases exponentially along x. The phase shift increases with , and x, and the velocity of the forced wave is √ 2D f . Accordingly, the forcing applied on the input feed of the TSFR generates a wave in the gel content that interacts with the Turing patterns. The effect of forcing is represented in the parameters of the sinusoidal input flow ( f and A f ) in Figure 4a. Three types of behavior are observed: pulsating spots (Figure 5c), spatiotemporal oscillations (Figure 5a), and mixed mode patterns (Figure 5b). Below A f = 0.08, only pulsating spots were found. Above this, forcing amplitude at low frequencies forced spatiotemporal oscillation to develop, while at high frequencies, pulsating spots formed. The forced oscillations and pulsation frequency match the forcing one, and the synchronization is perfect.
The amplitude of forcing (A f ) is an obvious operational control parameter, but what directly determines the dynamics of the gel content is the amplitude of the oscillations in the forced tank (A rel 0 ). Figure 4b shows the observed dynamics as a function of f , and A rel 0 . At stationary inflow concentrations Turing patterns in the range of 0.43 ≤ a h0 ≤ 0.46, that corresponds to a range of H + concentration in tank B 7.43 × 10 −6 ≤ h B ≤ 1.24 × 10 −5 . This defines domain of ±0.25 in relative change around the h B = 9.81 × 10 −6 at a h0 = 0.5 value. Accordingly, when A rel 0 is larger than 0.25, then the H + concentration in tank B during the forced oscillations exceeds the stability range of Turing patterns under constant flow.
At low values of f , as A rel 0 exceeds 0.25, mixed mode patterns first develop, then Turing patterns disappear. In that domain, the gel content shows oscillations, where a stripe with high H + concentration appears and disappears by following the forcing. At high values of f pulsating spots form as A rel 0 is smaller than 0.25. Interestingly, pulsating spots can be stable at medium values of f , although A rel 0 exceeds the stability range of Turing patterns under constant flow.
The amplitude of the response of patterns to forcing is presented in Figure 6a. We define the amplitude of patterns as the is the corresponding H + concentration at the unperturbed case (h st (x * ,y * ) = h (x * ,y * ) at A f = 0). These curves have a maximum at low frequencies and the pulsation of spots decays strongly with the forcing frequency.
The forcing of the patterns is performed indirectly, as the sinusoidal forcing is applied on the input of tank B. The variation of concentration of the activator in tank B sets the boundary condition for the gel and is characterized by A rel 0 (Equation 4). When A rel p is scaled with the amplitude of the oscillations in tank B (A rel 0 ), the response of the pulsating spots follows the same curve almost independently to the A f (Figure 6b). This is not the case when the forcing results in spatiotemporal oscillations (red points in Figure 6) or mixed mode patterns.
Topaz and Catllá studied the sinusoidal perturbation of Turing patterns near the Hopf-Turing bifurcation point, where the forcing drives the Hopf mode, which affects the Turing pattern modes. [21] They found a quadratic scaling between the forcing strengths and Turing pattern enhancement (E T ), as the square of the Hopf amplitude appears in the term that couples the Hopf mode to the Turing mode in the corresponding amplitude equations. According to their result the dependence of E T on the difference between the forcing ( f ) and the Hopf ( H ) frequency, is described by the following equation: H ) and c 1 , c 2 and c 3 are parameters. [21] At the conditions we used, our system is far from the Hopf-Turing bifurcation point, but the scaling between the amplitude of forcing (A rel 0 ) and the amplitude of pulsation (A rel p ) is almost quadratic (Figure 7a). In our system, oscillations are suppressed by the effects of the low mobility complexing agent (S) in the absence of forcing. The unforced system is far from the Hopf-Turing bifurcation point at the applied conditions. However, forcing supports the www.advancedsciencenews.com www.advtheorysimul.com reappearance of oscillations at low f above a critical A f (Figure 4). The quadratic scaling presented in Figure 7a indicates that the interaction of the reappeared forced oscillations and Turing patterns is similar to the one observed theoretically by Topaz and Catllá.
The response of spot pulsations also fits well to the theoretical prediction, as it follows the equation A rel (Figure 7b).
The pulsation frequency of the spots matches that of the forcing, but the phase is different. The phase shift (Δ ) between the oscillations of the forced tank B and the pulsation of a spot at (x * , y * ) is constant in time and increases with f (Figure 8). The pulsation of the spots is governed by the wave initiated by the time-periodic boundary condition. According to Equation (6) the phase shift of this wave increases with 1∕2 f and the distance from the boundary. In agreement with this equation, from the position of tank B (x = 4) toward tank A (x = 0), an increase of Δ can be seen in Figure 9a.
At the region of patterns, the chemical reactions strongly determine Δ . At the zone of patterns (around x = 1.5), a sharp decrease and an increase of Δ can be seen as x decreases. The spot pulsation drives the oscillations on the other side of the pat-  terns (0 < x < 1). Along a line of the maxima of the spots ((x * , y)) Δ changes periodically. The minima are at the center, while the maxima are at the border of the spots. In between two spots, where almost no reactions take place, Δ relaxes back to the

Conclusion
Sustained reaction-diffusion patterns can be observed experimentally in open reactors, where the continuous feed and removal of chemicals are typically provided by diffusion from the surfaces. The design of these reactors generates concentration gradients, and the control of the input of the chemicals is less direct than assumed. The effect of a sinusoidal forcing on the input feed composition was recently explored experimentally and numerically using simplified models. [20] In this study, we used a chemically realistic mechanism that includes a reaction responsible for differential diffusion. In the 2D description of the TSFR, we only neglected the feedback of the gel composition on that of the tank. This latter assumption is acceptable when the volume of the gel is small compared to that of the feeding tanks. We focused on the effect of a sinusoidal forcing on the input feed composition on the oscillations in the tank, on the spread of periodic signals in the gel, and principally on the dynamics of Turing patterns. The damping caused by the tank can be considered an essential but technical issue that must be calculated during the design of an experiment. We showed that the periodic changes of concentrations result in a wave that interacts with the Turing patterns in the gel. Besides oscillations and pulsating spots, a tiny range of mixed-mode behavior was also found that is similar to the complex Turing-Hopf mode behavior that have been observed theoretically and experimentally. [28,29] According to our numerical results, the basic properties of the forced Turing patterns far from the Turing-Hopf bifurcation point are similar to the theoretical predictions made at the vicinity of that. [21] The phase shift between the pulsation and the forcing is strongly determined by the local reactions.
The use of periodic input feed concentrations in reactiondiffusion systems is a promising way to control the dynamics. Furthermore, it might open a new direction in nonequilibrium synthesis, where the aim is to perform a reaction under conditions that allow producing materials with desired properties, which are challenging to make using standard synthetic methods. [30][31][32] The input feed concentration of B in CST-B and HA in CST-A was set to zero. The sum of the input feed concentration of A and HA in both tanks was fixed. The input feed of C is equal on both tanks and [S − ] tot is fixed. Equations (7)-(24) were used to calculate the boundary conditions.
The partial differential equations (Equations 7-12) were discretized with a standard second-order finite difference scheme on an equidistant grid having 100 × 100 gridpoints with a grid spacing Δx = Δy = 0.04. The resulting set of the ordinary differential equations was solved by the SUNDIALS CVODE [33] solver using the backward differentiation formula method. The absolute and the relative error tolerances were 10 −10 and 10 −5 , respectively.

Supporting Information
Supporting Information is available from the Wiley Online Library or from the author.