E Beauvier

S Bodea

A Pocheau

Front propagation in a regular vortex lattice : dependence on the vortex structure

Keywords: numbers: 47.70.Fw, 82.40.Ck, 45.10.Db

We investigate the dependence on the vortex structure of the propagation of fronts in stirred flows. For this, we consider a regular set of vortices whose structure is changed by varying both their boundary conditions and their aspect ratios. These configurations are investigated experimentally in autocatalytic solutions stirred by electroconvective flows and numerically from kinematic simulations based on the determination of the dominant Fourier mode of the vortex stream function in each of them. For free lateral boundary conditions, i.e. in an extended vortex lattice, it is found that both the flow structure and the front propagation negligibly depend on vortex aspect ratios. For rigid lateral boundary conditions, i.e. in a vortex chain, vortices involve a slight dependence on their aspect ratios which surprisingly yields a noticeable decrease of the enhancement of front velocity by flow advection. These different behaviors reveal a sensitivity of the mean front velocity on the flow sub-scales. It emphasizes the intrinsic multi-scale nature of front propagation in stirred flows and the need to take into account not only the intensity of vortex flows but also their inner structure to determine front propagation at a large scale. Differences between experiments and simulations suggest the occurrence of secondary flows in vortex chains at large velocity and large aspect ratios.

I. INTRODUCTION

In many reactive systems, reaction does not occur homogeneously but within thin fronts which propagate at some proper velocity V 0 in a still medium. However, when the medium is stirred, the advection and the distortion of reaction fronts by vortices considerably enhance the global reaction rate with important consequences regarding technology (e.g. combustion [START_REF] Williams | Combustion Theory[END_REF], chemistry [START_REF] Field | Oscillations and travelling wave in chemical systems[END_REF][START_REF] Epstein | [END_REF]), ecology [4][5][6] or geophysics (e.g. porous medium [7]). A large attention has then been devoted to understand the mechanisms of enhancement of reaction efficiency by a stirring flow in both a laminar [8][9][10][11][12][13][14][15][16][17][18] or a turbulent [START_REF] Williams | Combustion Theory[END_REF][19][20][21][22][23][24] context.

In this issue, it is essential to appropriately identify the parts of the flow that play a dominant role in the way reactive fronts invade a stirred medium. In particular, in flows composed of single-scale vortices, as in a vortex chain or a vortex array, many models have only relied on flow features at the vortex scale and overlooked the inner structure of vortices. They then considered the same prescribed vortex structure [8][9][10][11][12][13][14] and, drawing on the transport mechanisms of passive tracer in a vortex chain [25][26][27][28][29], focussed analytical derivations on the enhancement of the diffusion of contaminants across vortex separatrices [8][9][10][11], taken as the limiting phenomena for large scale front propagation. Although separatrices play an essential role in front propagation, these models nevertheless implicitly overlooked the possible implication of other parts of vortices. They then provided universal relationships between the mean front velocity and the flow intensity that make no difference between the types of vortices and no case of their inner structure.

Here, we revisit this issue by experimentally changing the inner vortex structure and further comparing the resulting effects on the front propagation. We find a large sensitivity of the mean front velocity to the vortex structure which questions the modeling of front propagation in stirred flows. In particular, this sensitivity stresses that the extremely small thickness of the fronts as compared to the vortex scale enables them to probe the details of the flows with benefit for their propagation [16,17]. This ability, which makes difference with passive tracers, calls for taking into account the inner structure of vortices to appropriately model front propagation in stirred media and further address the effects of flow instabilities or of multi-scale flows.

Our study involves a single scale flow made of a regular lattice of planar vortices involving same sizes and same depth and in which a reaction front can propagate. Front propagation is achieved experimentally by using an autocatalytic reaction in a solution stirred by electroconvective flows in a Hele-Shaw cell. In order to modify the inner structure of vortices in a controlled way, we use two external means which both refer to their boundary conditions : i) a change of lateral boundary conditions from free to rigid, ii) a change of the vortex depth.

In each configuration, a systematic experimental investigation is performed and the evolution of the mean front velocity with the vortex magnitude is determined. It reveals a variability of the mean front velocity with the boundary conditions and, for rigid boundary conditions, with the vortex depth. To check the role of flows in these outcomes, we analytically derive the vortex structure in each configuration in the Stokes regime and we determine the resulting mean front velocity by kinematic simulation based on a lattice dynamics algorithm. The agreement between experiment and simulation strengthens the experimental evidence of the responsibility of the vortex structure in the variability of the mean front velocity in vortex lattices. Altogether, our study thus provides original experimental, analytical and numerical results which converge to evidencing that front propagation in stirred media not only depends on the magnitude of vortices but also on their inner structure.

We denote U the vortex intensity (i.e. the maximum velocity of the flow field), V 0 the proper front velocity (i.e. the velocity of the front with respect to the fluid), V f its effective velocity (i.e. the mean front velocity in the mean fluid frame), L the vortex width, d the cell depth and λ the front thickness. In this framework, the goal of previous modelings was to relate these variables by determining the relationship between their non-dimensional forms : V f /V 0 = f (U/V 0 , d/L, λ/L). However, the main issue addressed in this study will be to question whether a single function f (.) is sufficient for representing the various configurations of vortex flow or whether other flow parameters referring to the inner vortex structure will prove to be also relevant.

Section II describes the experiment and specifies the different lateral boundary conditions applied on vortices. Section III addresses the analytical determination of planar electroconvective flows in the Hele-Shaw regime together with the implications of different boundary conditions on the vortex structure. Section IV reports a kinematic simulation of front propagation in the different flows determined in section III. Section V discusses the results obtained on front propagation enhancement for various vortex aspect ratio d/L and boundary conditions as well as the origin and the implications of the sensitivity of front propagation on the vortex structure. A conclusion on the study is brought about in section VI.

II. EXPERIMENT

The experiment aims at achieving front propagation in a periodic array of steady vortex flows in well controlled conditions suitable for a quantitative study of its dependence on the vortex structure. For this, it is essential to avoid any feed-back of the reaction onto the flows, as displayed for instance in combustion due to gas expansion following heat release. This is achieved by considering an athermal oxydo-reduction reaction in solution. Next, the flows must be steady, which discourages to generate them from primary instabilities (as thermoconvective or Rayleigh-Taylor instabilities), since they are known to induce phase versatility in the vortex pattern and secondary instabilities. Instead, flows have to be forced so as to be under control. This led us to use electroconvective flows in a sufficiently dissipative hydrodynamical regime to avoid flow dynamics. Restricting the study to a Stokes regime in a Hele-Shaw cell, we thus used steady electroconvective flows, directly monitored by magnetic intensity and electric current.

A. Set-up

The chemical reaction used is the chlorite-iodure oxydo-reduction reaction in aqueous solution :

ClO - 2 + 4I -+ 4H + → Cl -+ 2I 2 + 2H 2 O (1) 5ClO - 2 + 2I 2 + 2H 2 O → 5Cl -+ 4IO - 3 + 4H + (2) IO - 3 + 5I -+ 6H + → 3I 2 + 3H 2 O (3)
It is autocatalytic as the primary consumption of iodide I -in ( 1) is enhanced in (3) by the resulting production of iodate IO - 3 via iodine I 2 in (2). The reactants are mixed in a solution which thus starts from a chemically homogeneous state. It is initially blue due to the formation of a complex between iodide, iodine and a Starch indicator introduced for this purpose, and turns uncolored after reaction following the depletion of iodide.

Different modes of reaction may occur : i) a global mode in which the whole medium evolve the same way on any of its parts with, therefore, a homogeneous reaction rate, ii) a local mode in which the reaction rate is localized in a thin interface, called reaction-diffusion front, which propagates at a definite velocity V 0 in a still medium. In the latter case, as the medium is at rest, the phase change induced by the front corresponds to a wave. Its propagation mechanism consists in the diffusion of a contaminant to neighbor parts of the front that then react due to a consecutive enhancement of their reaction rate. A popular analog corresponds to the burning by a fire front where heat diffusion yields the fresh parts close to the front to burn, yielding fire advance and propagation. The simplest, one-dimensional, modeling of the reaction-diffusion wave is provided by a contaminant balance :

∂θ ∂t = D θ + f (θ) (4)
where θ denotes the contaminant, D its diffusivity and f (.) the reaction rate. This dynamics includes propagation of a wave θ(x, t) ≡ θ(x -V 0 t) at a velocity V 0 that is solution of the non-linear equation :

-V 0 ∂θ ∂x = D θ + f (θ) (5)
with appropriate conditions at infinity. It thus corresponds to a non-linear eigenvalue [30][31][32]. The alternative between the two kinds of solutions depends on initial conditions. In practice, the front solution needs a germ to develop. Here, we force it by initiating, at an extremity of the set-up, a localized depletion of iodide by oxydo-reduction at a metallic plate. We note that the velocity V 0 [START_REF] Saul | Oscillations and travelling Waves in chemical system[END_REF]34] slightly depends on the temperature of the solution due to Arrhenius factors controlling the chemical kinetics and on the electric current density following its implication on ionic transport [START_REF] Saul | Oscillations and travelling Waves in chemical system[END_REF]34]. However, this latter dependence is weaker when the front propagates on a direction opposite to the electric current [17]. Therefore, to minimize the variability of V 0 to less than few percents, the temperature has been stabilized to 20 • C and the front has been initiated so as to propagate in average in a direction opposite to that of the electric current. We then obtained values of V 0 varying between 1.0 and 1.4 mm/mn depending on the current density.

The homogeneous solution is poured into a long glass recipient whose depth d is controlled by spacers. This channel is then placed above an array of magnets whose magnetic field B alternates direction between neighbors (Fig. 1-a). An electric current density j, parallel to the channel axis, is then induced in the solution by electrodes placed at its extremities. Following the presence of a magnetic field, it generates a density of Laplace forces in the solution which drives alternate vortices. This electromagnetic forcing has been used in many different contexts to study for instance boundary layer effects in salty water [35,36], mixing [29,37], single scale [38][39][40] multiscale [41] or quasi-two dimensional turbulent or chaotic flows [42][43][44][45], front propagation [15][16][17], fiber deformation [46] or sedimentation [47]. In the Stokes regime in which our experiment stands, the vortex flow intensity is linearly related to the density of Laplace force and thus to both the current density and the magnetic field. Both can be tuned either by the distance of the solution to the magnet array or by the tension applied to the electrodes. This yields vortex intensities up to 50mm/mn for arrays placed beneath the solution, a channel depth of d = 2 mm and current intensities of about 30 mA.

Initiation of a front is triggered by placing a small metallic piece at a channel end. A slight dependence of front velocity on current density arises following the resulting transport of charged species [16][17][18] and is taken into account in the following. As the front propagates, the solution, observed from above, splits into blue unreacted zones and uncolored reacted zones, the front corresponding to their thin frontier (Figs. 2 and3). This provides the opportunity of finely localizing it by a nonintrusive optical observation. Front propagation is finally recorded by a camera onto a computer so as to provide the front shapes and the front velocity, both locally at each time and globally in average over several vortices (Figs. 2 and3). Following fluid stirring, the front turns distorted and enrolled around the vortices. In addition, its velocity along the channel axis appears to be largely increased by the flow, in proportions that we shall analyze with respect to the boundary conditions.

B. Flow boundary conditions

The channel is 400 mm long, 20, 100 or 120 mm large and involves depths of either 2, 4, or 6 mm. A single kind of square magnets is used with dimensions 20*20*5 mm, yielding the channel widths to extend over 1, 5 or 6 magnets.

All flows experience rigid (no-slip) boundaries at the top and bottom faces of the channel, following which their profile in the depth direction is close to a Poiseuille profile.

However, on the lateral directions, two kinds of situations are encountered depending on the channel width. When the channel width is large compared to the magnet dimensions, most vortices are surrounded by vortices. They then experience free (free-slip) boundary conditions (hereafter denoted b.c.) on their all frontiers, except for the two extreme vortices which touch the channel sides (Figs. 1-b and2). As the leading parts of the front actually travel into the central vortices only (Fig. 2), this configuration thus refers to effective free lateral b.c.. On the opposite, when the channel has the same extension than a magnet, a single vortex takes place in between its lateral sides. Vortices then experience both free b.c. at their frontier with neighbor vortices and rigid b.c. on the lateral sides (Figs. 1-c and3). As front will be advected along these sides, this configuration thus refers to effective rigid lateral b.c.. Figures 4 and5 display the velocity field measured by particle image velocimetry (PIV) on a field extending over four vortices. They all show the regular flow structure obtained for both boundary conditions in the midheight plane of the channel. The vanishing of velocity at the lateral sides of the channel for effective rigid lateral b.c. is noticeable in figure 4. In contrast, for effective free lateral b.c., no vanishing occur in figure 5 between vortices.

Following the narrow depth of the channels, the Reynolds number is low (cf. section III) in which case steady stable flows are expected. The steadiness of the flows is confirmed by PIV measurements together with the absence of symmetry breaking in the vortex chain (the slight asymmetry of the flow field in figure 3 is a visual artifact induced by a phase mismatch between the PIV lattice and the vortex axes of symmetry). On the analytical ground, the onset of stability of planar vortex lattices has been studied in configurations involving the two types of lateral boundary conditions [48]. Two parameters arise: the inverse Re L for d = 6mm. Then, stability analysis confirms that, in our velocity range, the studied vortex flows are actually stable [48]. In particular, the onset of instability stands above a velocity of 1000mm/mn (resp. 9500mm/mn) for d = 6mm (resp. d = 2mm), far from the maximum of our velocity range, 50mm/mn. Although this analysis considered only an even number of vortices in the channel width, the onset of instability found is so large that no instabiltiy may be expected in our parameter range.

Both free and rigid lateral b.c. could therefore be achieved within the same set-up, simply by tuning the channel width. In both cases, the channel depth d was varied between a short value 2 mm, an intermediate value 4 mm and a large value 6 mm. The objective of the experiment is then to address the implications of these different b.c. and aspect ratios on front propagation.

C. Front propagation

Propagation regime

The propagation regime depends on two nondimensional numbers P e = U L/D and Da = L/U τ where τ denotes the characteristic reaction time close to the fresh state θ = 0, θ being the progress variable of the reaction and D the molecular diffusivity of the relevant species.

Then, for a concave reaction rate f (θ) (4), theories select a front velocity V 0 yielding τ = 4D/V 2 0 , Da = (LV 0 /4D)V 0 /U and P eDa = (LV 0 /2D) 2 [30][31][32]. With V 0 = 1mm/mn, L = 20mm, D = 2.10 -3 mm 2 /s and 2 < U < 50mm/mn here, this yields the experimental ranges 300 < P e < 8500, 0.8 < Da < 21 with a constant product P eDa ≈ 7000. This corresponds to the thin front regime (P e >> 1) and within it, either to the flamelet regime for large Da where the front is slightly distorted by the flow or to the distributed regime for Da about unity or lower, where the front, engulfed by many vortices, displays a long wake of burning vortices behind it [13,14,17,[START_REF] Borghi | Theory of flames[END_REF]. In both cases, the front thickness λ is small compared to the vortex size L. In particular, here, λ should be about 1mm for a parabolic reaction rate f (.) and is likely less from direct observation, so that λ/L < 5.10 -2 .

In this thin front regime, the front then satisfies an eïkonal dynamics where, apart from curvature effects, it may be considered as a line advancing at a proper velocity V 0 with respect to the fluid. In addition, in the present Hele-Shaw cell, the effect of the large velocity gradient in the channel depth direction may be determined with respect to a Péclet number η = dV 0 /2D based on the channel depth d and the planar front velocity V 0 [START_REF] Edwards | [END_REF]. For low η or low d, the front displays a large curvature on the depth direction following which it globally advances at V 0 plus the depth averaged component of the flow velocity on its normal. In contrast, for large η or large d, curvature is weak enough for making the front advance at velocity V 0 plus the maximum on the channel depth of the velocity component on its normal n, i.e. V 0 + V(x, y, 0).n for a Poiseuille flow, z = 0 denoting the middepth [START_REF] Edwards | [END_REF]51]. Here, as η > 30 at the lowest channel depth, the experiment stands in this large gap regime. The front then advances in the whole channel depth at a single velocity corresponding to a propagation at velocity V 0 in a fluid moving at velocity V(x, y, 0). Whereas, the experiment is actually 3-dimensional, front propagation in a Poiseuille flow is thus equivalent to a two-dimensional propagation in the mid-depth plane.

Effective velocity and vortex intensity

The propagation of the front has been studied at various current densities j for both free and rigid lateral b.c. and various channel depths. In all cases, it led a periodic motion of the front along the channel direction. The mean velocity of the front displacement along this axis provided its effective front velocity V f . This effective velocity results from a balance between two phases : a quick transport from one vortex side to the other for which the front benefits from vortex advection at the large speed U ; a slow crossing of separatrices to enter the next vortex at the low speed V 0 since advection no longer helps the front to propagate on the direction of the channel axis at these locations. The overall result is a mean velocity standing in between the large U and the weak V 0 velocity to a proportion that we seek to study here.

The vortex intensity U is defined as the maximum velocity of the flow field V. It is not measured directly but deduced from the measurement of the current density j to which U is proportional in a given configuration. The proportionality constant between U and j was obtained from the detailed study of front dynamics in a vortex. As shown in figures 2 and 3, the front gets advected from one separatrix to the next once it has entered a vortex. This advection occurs far from the vortex center and close to the quickest vortex streamline, as depicted in a previous study [17]. During this motion, the displacement of the front by propagation in the co-moving medium is negligible, since U is quite large compared to V 0 . Accordingly, the mean velocity during this advection phase stands as a relevant indicator of the vortex intensity. As expected, it is actually found to be proportional to the current density j. streamline, this mean velocity is slightly lower than the maximum flow velocity U . To take into account this bias, we determined both of them on different runs of the simulations of front propagation reported in section IV. They were found proportional one to the other, the measured mean velocity being lower than the vortex intensity U by a factor 0.73 (resp. 0.79) for free (resp. rigid) lateral b.c.. Together with experimental measurements, this provided us with the proportionality constant between the current density j and the vortex intensity U , in any configuration.

However, as the flow velocity weakly varies on a flow

To facilitate the comparison between the experimental results and the literature on simulations or modelings, we shall adopt the vortex intensity U as the scalar indicator of the vortex magnitude. Its linear difference with the mean velocity explains the data differences with reference [18].

The measured mean velocity over front advection from one separatrix to the next is given by πr/t where r is the mean radius of the corresponding streamline and t the advection time. Given a picture length of 760 pixels englobing 4 vortices and the approached value r ≈ 2L/3 found here, one gets with an uncertainty of ±1pixel at each extremities of r and of ±1s for t, the relative uncertainty δU/U ≈ 3.10 -2 . Similarly, measuring the mean effective velocity over 4 vortices and the same picture length, one obtains V f = 4L/T where T denotes the front travel time over these vortices. For the same basic uncertainties, one then gets the relative uncertainty δV f /V f < 1.4.10 -2 .

Free boundary conditions

Figure 2 displays a propagation sequence over a period corresponding to the crossing of two counter-rotating vortices. Five vortices stand within the channel width. However, the leading part of the front, which sets the mean velocity, stands on the three central vortices for which free b.c. are in order. The front motion gets quite periodic soon after the front initiation, typically till the second vortex crossing, so that an asymptotic dynamics is reached well before the channel end. See the experimental movie [57] performed at d = 2, L = 20, U/V 0 = 23.1 and free BC.

Figure 6 reports images of front propagation for three reduced vortex intensities U/V 0 and three channel depths. It shows similar front shapes displaying thin tongues around separatrices and wakes extending behind the fronts until all vortices are burnt. However, tongues get thinner and wakes get longer as the relative flow velocity U/V 0 increases. The latter effect results from a decreasing propagation efficiency in a vortex in comparison to advection. In particular, as contaminated vortices burn mainly by radial propagation of a front towards the vortex center with no benefice from advection, they take a fixed time t ≈ L/V 0 to disappear. Meanwhile the leading front advances over a distance V f t ≈ LV f /V 0 which is found to be proportional to U/V 0 (Fig. 8-a). Except close to boundaries, front wakes are thus all the longer than vortex intensity is larger.

Figure 8-a reports the enhancement of the mean front velocity V f by the vortex flow U . It is expressed with the non-dimensional variables V f /V 0 , U/V 0 . Data show a linear raise with no difference regarding the channel depth d. The change of vortex structure with d, if any, is thus insufficient to noticeably modify front propagation. In agreement with the expectation of section II C 1, this supports the picture of a leading part of the front moving in the mid-height plane where the flow is the largest and of a planar flow whose profile on the depth direction does not affect front propagation.

Rigid boundary conditions

Figure 3 displays a propagation sequence somewhat similar to that reported for free b.c., except that a single vortex stands within the channel width. In particular, the quick phase of front advection between successive separatrices takes place at a distance of about half the channel depth from the lateral boundaries [16,17]. As this is the expected width of the lateral boundary layer in a Hele-Shaw cell, the front trajectory thus presumably fits within this boundary layer. See the experimental movie [58] performed at d = 2, L = 20, U/V 0 = 17.7 and rigid BC.

Figure 7 reports images of front propagation on conditions similar to figure 6 but for rigid lateral b.c.. Here too, the scenario of front propagation over the vortex chain looks quite similar, independently of flow inten-sity and of channel depth. However, the frontier between fresh and burnt domains seems to enlarge here with the channel depth. As images are taken from above, this reveals that fronts are less superposed on the channel depth direction than for free lateral b.c., this effect increasing with the cell depth. Finally, some patterns within burning vortices begin to appear at the large channel depth d and the largest flow velocity U/V 0 . We attribute this to secondary flows raising with both d and U/V 0 .

Figure 8-b reports the enhancement of the reduced mean front velocity V f /V 0 by the reduced vortex flow U/V 0 . Different enhancement curves are obtained depending on the channel depth d. In particular, they exhibit at d = 4 or 6 mm a concavity which differs from the linearity displayed for free b.c. and which noticeably increases with d. For instance, at U/V 0 ≈ 34, increasing the channel depth from d = 2mm (V f /V 0 = 11.7) to d = 6mm (V f /V 0 = 6.8) decreases the mean front velocity V f of about 40%. This reveals that different vortex structures are in order depending on the aspect ratio d/L and that their differences are sufficient for noticeably affecting the enhancement of front propagation, even at low vortex intensity.

This observation means that a tiny variation of the vortex structure may yield noticeable implications on front propagation. This sensitivity reveals that inner details of the flows, i.e. their sub-scale structure, is important in the large-scale behavior of front dynamics, a fact that is overlooked in most large-scale theories or modeling of front propagation [8][9][10][11][12][13][14]. The remaining of this paper is devoted to analyze this situation first by investigating a model capable of quantifying the implication of b.c. on the flows and then by determining the implications of flows on front propagation by numerical simulation.

III. STREAM FUNCTION MODELING

We call V the flow velocity, B the magnetic field, j the electric current density, p the fluid pressure, ρ its volumic mass and ν its kinematic viscosity. We denote (x, y, z) the directions along the channel length, width and depth (Fig. 1) and (e x , e y , e z ) the corresponding unit vectors. The vortex intensity is still labelled U .

Assuming a planar flow in the bulk, the advective term V.∇V may be estimated to about U 2 /L since the flow variation scale in the flow plane is L [59]. However, in case of rigid lateral b.c., this scale decreases to about the channel depth d near the boundaries, yielding another estimate of about U 2 /d there. On the other hand, the dissipative term ν V may be estimated to νU/d 2 for a near Poiseuille profile. Accordingly, the Reynolds number reads Re = U d/ν close to a boundary and Re = U d 2 /νL, i.e. Re = Re d/L, in the bulk. As the water viscosity is about 1 mm 2 /s and L = 20 mm, one gets at d = 2 mm, Re = 1 for U = 30 mm/mn and Re = Re/10. Accordingly, in our experimental range where U ≤ 50mm/mn and 2 ≤ d ≤ 6mm, we may assume a Stokes regime (Re 1), except near the boundaries at large velocities (Re > 1). Neglecting advection terms provides a linear relationship between the flow and its sources, that we shall solve mode by mode.

Taking into account the density of Laplace forces f l = ρj ∧ B, we thus obtain, at a distance of order d from the boundaries :

- ∇p ρ + ν V + j ∧ B = 0 (6) ∇B = 0 ; ∇ ∧ B = 0 (7) ∇j = 0 (8) ∇V = 0 (9)
Here B stands for the magnetic field produced by the magnets, the one that is induced by the electric current being quite negligible in comparison. Relations ( 8) and ( 9) stand for charge conservation and flow incompressibility.

Taking the curl of equation ( 6) and using the divergenceless character of j and B yields the vorticity equation :

Ω = ν -1 (j.∇)B = S ( 10
)
where Ω = ∇ ∧ V denotes the vorticity and S the density of vorticity sources. On the other hand, relations (7) show that B is Laplacian free : B = 0. Assuming a uniform current density j = je x , in agreement with parallel electrodes placed at the extremities of the channel in the present set-up, then yields S = 0, ∇ ∧ S = 0 and V = 0. In the following, we consider a channel of width l containing an odd number of square magnets of length L, placed so that the channel boundaries fit with the magnet boundaries. We place the origin at the middle of the channel width, on a magnet boundary on the x axis and at half the channel depth. This configuration forces the magnetic field to be antisymmetric with respect to the y-axis and symmetric with respect to the x-axis. As it is periodic with a period 2L on both the x and y axis and symmetric with respect to all the magnet centers, this selects its Fourier modes to be B k,p (z) sin[(2k +1)πx/L] cos[(2p+1)πy/L] ; (k, p) ∈ N 2 . Following (10), this yields the vorticity sources density S to read :

S = ∞ k=0 ∞ p=0 S k,p (z) cos[(2k + 1)π x L ] cos[(2p + 1)π y L ]
(11) As S = 0, the Fourier components S k,p (z) must satisfy :

d 2 S k,p dz 2 -a 2 k,p π 2 L 2 S k,p = 0 (12)
with

a k,p = [(2k + 1) 2 + (2p + 1) 2 ] 1/2 (13)
They thus express as exponentials. However, having placed the magnets below the fluid layer, we shall select the sole exponential that does not diverge at z = +∞ :

S k,p (z) = s k,p e -a k,p πz/L (14)
We now assume that the flow is planar, V.e z = 0. This allows us to seek it in term of a stream function ψ(x, y, z) : V = ∇ ∧ (ψe z ) = ∇ψ ∧ e z . The vertical vorticity Ω z = Ω.e z then writes Ω z =h ψ where h = ∂ 2 /∂x 2 + ∂ 2 /∂y 2 denotes the horizontal Laplacian. Labelling S z the z-component of S, one thus obtains from (10) :

h ψ = -S z (15)
A. Boundary conditions and general solution

As equation ( 15) is linear in ψ, its solution to definite b.c. corresponds to the sum of a particular solution ψ p and of complementary functions that are solutions of the homogeneous equation :

h ψ = 0 ( 16
)
Our goal now consists in determining for each mode (k, p) of B, the resulting mode induced on the stream function ψ. This determination will be found below to depend on boundary conditions. We note that it a priori requires to explicitly express the modes of B from the distribution of its sources. However, we shall show in section III D that, among the modes of ψ, the fundamental mode will be overdominant, leaving the determination of the other harmonics and thus of B itself, useless.

Boundary conditions

Rigid b.c. take place on the solid plates of the channel, i.e. the top and bottom plates z = ±d/2 and the lateral faces y = ±l/2. They impose a vanishing of the velocity there, V = 0, and thus a constant ψ on the formers and a uniform ψ on the x-axis on the latters.

In all horizontal planes, we fix at 0 the uniform value of ψ for no flow. We then notice that, in each of them, as vortex separatrices correspond to streamlines that are common to counter-rotating vortices, they must keep unchanged under the symmetry ψ → -ψ that reverses velocity and exchanges the vortex types. Accordingly, they correspond to ψ = 0. One of these streamlines is imposed by the channel plates at y = ±l/2. Those that occur in the bulk are however a priori unprescribed and will result from both the distribution of the vorticity sources and the b.c.. Finally, the intersection between vortex separatrices and the top, bottom, and lateral plates yields, by continuity, ψ = 0 on them.

In term of stream function, the b.c. thus read :

• on the top and bottom faces z = ±d/2 : ψ = 0 ψ(x, y, ±d/2) = 0 (17)

• at the lateral faces y = ±l/2 :

ψ = 0, ∇ h ψ = 0 ψ(x, ±l/2, z) = 0 (18) ∂ψ ∂y (x, ±l/2, z) = 0 (19)
• at vortex separatrices : Vortex separatrices stand at abscissa x b = ±L/2 mod(L) following the forcing by the magnet lattice and, in case of free b.c., on ordinates y b = ±L/2 mod(L), by reason of symmetry. Then : ψ(x b , y, z) = 0 ; ψ(x, y b , z) = 0 (20) with, because of effective free b.c. on these separatrices, the continuity of stresses

σ i,j = ν/2(∂ i V j + ∂ j V i ), or, equivalently : continuity of ∂ 2 ψ ∂x∂y and of ( ∂ 2 ψ ∂y 2 - ∂ 2 ψ ∂x 2 ) ( 21
)
2. Particular solution Following (11), a particular solution ψ p of equation ( 15) can be sought as a Fourier series :

ψ p = ∞ k=0 ∞ p=0 ψ k,p (z) cos[(2k + 1)π x L ] cos[(2p + 1)π y L ] (22)
with, following (14) :

d 2 ψ k,p dz 2 -(a k,p π/L) 2 ψ k,p = s z k,p
(a k,p π/L) 2 e -a k,p πz/L (23

)
where s z k,p denotes the z-component of the Fourier coefficient s k,p .

As the right member of equation ( 23) belongs to the kernel of the operator of its left member, a secular solution is generated. It writes :

ψ k,p = - s z k,p
2(a k,p π/L) 3 z e -a k,p πz/L (24)

Complementary functions

Complementary functions are solutions of the homogeneous equation (16). As they belong to the kernel of a composition of operators, secular solutions are expected. Added to the particular solution (22) (24), they will enable the b.c. on the z and y directions to be satisfied.

The linearity of equation ( 16) allows us to seek its base of solutions in terms of functions of separate variables. In addition, as we disregard b.c. on the x axis, the expected stream function ψ as well as its complementary function ψ c = ψ -ψ p will thus be, as the particular solution ψ p , periodic in x and with the same period 2L. Accordingly, the relevant base of solutions for complementary functions will involve functions of the form :

ψ c = f (z)g(y) cos(γx + ξ) ; γ = (2k + 1)π/L k ∈ N (25)
where no phase shift ξ will be required in the following. In addition, following the symmetry of the addressed configuration with respect to the plane (x, z), attention will be restricted to even functions g(.).

The corresponding complementary functions are identified in appendix A by a systematic procedure. Only six kinds of real complementary functions, even in x and y, are selected :

ψ 1 = f 1 (z) cosh(

General solution

The general solution to equation (15) writes as the sum of ψ p and of complementary functions ψ c , up to coefficients to determine according to the b.c. to satisfy. We use below the allowable complementary functions to satisfy sequentially the b.c.. Solutions to free or rigid b.c. will be labelled with a superscript "f" or "r" respectively. Considering a Fourier component (k, p) of ψ p asks for complementary functions involving sinusoidal factors cos(βy) cos(γx), with γ = (2k + 1)π/L and β = (2p + 1)π/L. Among the six kinds of complementary functions, only one, ψ 5 , satisfies this criterion. The corresponding functions involve an exponential factor e αz with α 2 = β 2 + γ 2 with, therefore, two opposite values of α, yielding two independent complementary functions. These are sufficient to satisfy the two required b.c. at z = ±d/2 while preserving the free lateral b.c.. The resulting combination writes :

ψ f = ∞ k=0 ∞ p=0 ψ f k,p (z) cos[(2k + 1)π x L ] cos[(2p + 1)π y L ]
(26) with :

ψ f k,p = s z k,p d 4 32 b 3 k,p -ze -b k,p z (27)
+ e b k,p sinh[b k,p (z -1)] + e -b k,p sinh[b k,p (z + 1)] sinh(2b k,p ) b k,p = a k,p π 2 d L
a k,p being given by ( 13) and z = 2z/d with -1 ≤ z ≤ 1.

The difference between this determination and other modelings [52] traces back to the variation of the eigenvalue a k,p with the mode indexes (k, p), which is explicitly considered here in order to satisfy the properties (7) of the magnetic field.

It appears that the functions ψ f k,p (z) are slightly asymmetric following the asymmetry of the vorticity source density (14). They thus differ from both a sinusoidal mode (Fig. 9-a) and from a parabolic profile (Fig. 9-b). This therefore reveals a difference between the present flow induced by a vorticity source and the Poiseuille flows induced by a pressure gradient. In practice however, the difference is negligible on the fundamental mode k = 0, p = 0 (Fig. 9-b).

C. Rigid lateral boundary conditions

In order to satisfy, for each Fourier mode k ∈ N , the additional constraints (18) and (19)

k,p (z) that is C 1 in [-d, d] (Appendix B
). This yields :

ψ f n,k,p (x, y, z) = ψ f n,k,p cos[nπz/d + ϕ n ] cos[(2k + 1)πx/L] cos[(2p + 1)πy/L] (28)
where, following the odd symmetry of the function ψf k,p (z) at z = ±d/2, the phase ϕ n is π/2 (resp. 0) for odd (resp. even) n and where no mode n = 0 is in order following the zero average of ψf k,p (z) in [-d, d] . As equation ( 15) is linear, we may seek its solution for each modes (n, k, p) separately and independently. For this, we look for complementary functions ψ c n,k,p (.) suitable for making the resulting stream function ψ r n,k,p (.) = ψ f n,k,p (.) +

yielding the Fourier amplitude :

ψ f n,k,p = -S z n,k,p (L/π) 4 [(2k + 1) 2 + (2p + 1) 2 + n 2 ] -1 [(2k + 1) 2 + (2p + 1) 2 ] -1 (31)

D. Fundamental mode

Returning to the expression of the Fourier amplitudes ψ f k,p (z) (27) relevant to our issue, we notice that it quickly decreases with k and p. In particular, ψ f k,p (0) = s z k,p b -3 k,p tanh(b k,p ) d 4 /32 with s z k,p usually decreasing with k and p (in particular, we notice that it vanishes at infinite (k, p) when S z is a C 1 function). Accordingly ψ f k,p (0) decreases at least as the cube of k and p with already ψ f 1,0 (0)/ψ f 0,0 (0) < 5 -3/2 tanh √ 10/ tanh √ 2 < 0.1. This allows to restrict attention in the following to the fundamental mode (k, p) = (0, 0) and to neglect harmonics. In addition, we notice that the Fourier expansion of ψ f 0,0 (z) is dominated by the fundamental mode n = 1 :

ψ f 0,0 (z) = 1.03 cos[πz/d] -0.04 cos[3πz/d] -0.03 sin[2πz/d] + 0. sin[3πz/d] + 0.004 sin[4πz/d] + • • • (32)
the amplitude of the mode n = 3 being about 4% of that of the fundamental while that of the mode n = 2, about 3% of the fundamental, reflects the slight asymmetry of ψ f 0,0 (z). Here too, this legitimizes to neglect higher harmonics n > 1. Altogether, this yields us to consider the fundamental mode (n, k, p) = (1, 0, 0) as a relevant approximate solution with the corresponding expressions :

• for free b.c. on the y axis :

ψ f 1,0,0 (x, y, z) = ψ f 1,0,0 cos[πz/d] cos[πx/L] cos[πy/L] (33)
This is the solution usually adopted in theories or modeling [9-13, 46, 47]. Its streamlines are depicted in figure 10.

ψ r n,k,p (x, y, z) = ψ f 1,0,0 cos[πz/d] cos[πx/L] cos[πy/L]+(-1) q cosh[βπ µ/2] cosh[π y/L] -cosh[βπ y/L] cosh[π µ/2] cosh[βπ µ/2] sinh[π µ/2] -β sinh[βπ µ/2] cosh[π µ/2] (34) with µ = l/L = 1 + 2q, q ∈ N and β 2 = 1 + (L/d) 2 .
In particular, for a channel involving a single vortex between its lateral boundaries, as in the present experiment, one has µ = 1 and q = 0. The corresponding streamlines for the fundamental mode (n, k, p) = (1, 0, 0) are then depicted in figure 11 for different aspect ratios d/L.

Interestingly, one notices that, for free b.c., the stream function ( 33) is independent of the aspect ratio d/L whereas for rigid b.c. it actually depends on it, following the dependence of the parameter β on d/L (34). Accordingly, one may expect, for free b.c., a constant vortex structure independent of the aspect ratio d/L and, for rigid b.c., a vortex structure varying with this aspect ratio.

The dependence of the stream function with the aspect ratio d/L is noticeable in figure 11 in the streamline ellipticity which increases with it. Also, one notices the striking similarity between streamlines for free b.c. (Fig. 10) and for rigid b.c. at small aspect ratio d/L = 2/20 (Fig. 11-a). As boundary layers in a Hele-Shaw cell scale with d at fixed L, rigid b.c. yield an effective free b.c. at the end of boundary layer but for a reduced vortex size. For small d/l, this reduction is minor so that the vortex looks similar to that obtained for free b.c.. In contrast, for larger d/L, the reduction is large enough for making the vortex streamlines noticeably elliptic.

IV. SIMULATION OF FRONT PROPAGATION

In agreement with the absence of noticeable 3dimensional effects in the experiment and with the usual assumption of planar flows in Hele-Shaw cells, we intend to simulate front propagation in a plane stirred by a 2dimensional flow satisfying the above flow modeling. The main objective of the simulation will be to recover the effects of b.c. on front propagation so as to point out their crucial influence on velocity enhancement. As the actual front propagation is driven in the plane where the flow is maximal, i.e. the mid-depth plane z = 0 up to the slight asymmetry pointed out in (27), the stream functions at z = 0, ψ f (x, y, 0) and ψ r (x, y, 0) will be considered as those of the 2-dimensional flow used in the following.

To simulate front propagation, we do not need to detail the sharp transition between burnt and fresh domains but only its spatio-temporal evolution. We may thus skip the simulation of the concentration fields by reactiondiffusion-advection to focus on a simulation of the advance of the transition between the two domains. This corresponds to turning to a kinematic simulation of front propagation. In addition, to avoid the geometrical complexity inherent to a moving line, we will not simulate the line dynamics itself but the dynamics of the domains that it separates. These dynamics will be based on two transport phenomena : an advective process modeled by a Lagrangian transport at velocity V and a contamination process at the proper front velocity V 0 modeled by simple lattice dynamics.

The advantage of the simulation will be to provide front propagation at a sufficient accuracy to discuss boundary effects while escaping geometric complexity. The disadvantage will be to overlook the front thickness and the effect of large curvature. However, the former being smaller than any other characteristic length scale (λ d) and the relevant curvature radii ρ being large compared to the front thickness in the experiment (ρ > ∼ O(d) λ), this bias will induce no prejudice on the modeling of front propagation here.

As the transition between burnt and fresh domains is sharp, they both refer to uniform concentrations to which discrete values, 1 and 0 respectively, will be assigned.

The objective of the simulation is then to implement both their advection by the flow field V and the contamination of one by the other at speed V 0 . Simulating this evolution looks similar to implementing the so-called G-equation [START_REF] Williams | Combustion Theory[END_REF] in which a continuous field G relating a reacted domain where G = 1 to an unreacted one where G = 0 undergoes the dynamics DG/Dt = V 0 |∇G|, D/Dt denoting the Lagrangian derivative attached to the flow field V. However, the difference is that the field considered here is discrete, so that its evolution turns out to be kinematic.

We describe below the mapping used to implement this kinematic model on a lattice.

A. Mapping

A lattice is designed to simulate a chain of n vortices of width L extended over distances L x = nL, L y = L in the x and y direction. As the lowest velocity is V 0 , the elementary time and spatial steps ∆t and ∆x are linked by ∆x = V 0 ∆t. The simulation of the vortex chain therefore requires (L x /∆x, L y /∆x) pixels. Choosing ∆t = 1s, this corresponds with V 0 = 1mm/mn and L = 20mm to (1200n, 1200) pixels. The kinematic evolution of domains during an elementary time step ∆t is achieved by successively applying two mappings, one corresponding to advection and the other to contamination.

Advection of a pixel by the flow V on a time step ∆t is achieved by translating its position along a streamline at speed V . A forward Euler method was preferred to a Runge-Kutta method which appeared more time costly without noticeable improvement. However, to better take into account streamlines curvature, the time step ∆t has been decomposed into V /V o sub-steps on which the streamline direction, the flow intensity and the resulting advection have been computed. The final position hence obtained is finally projected on the nearest pixel of the discrete lattice.

Altogether, it appeared that the advection mapping alone conserved trajectories on streamlines at a quite satisfactory accuracy. For instance, at the vortex medium streamline ψ = 1/2, the variation of ψ was only δψ = 7.5 10 -4 after a vortex turn-over time T = πL/U and even smaller on more external streamlines. Meanwhile, as V f ≈ U/3 for U > 5V 0 (Figs. 13-a,b), the front head then already reached the next vortex : V f T ≈ π/3L > L. This ensures that the accuracy in the advection scheme is quite sufficient for adequately modeling front propagation in the vortex chain. On the other hand, care has been taken to address sufficiently large flow velocities to make the advection distance crossed over the time step ∆t large enough compared to the pixelization so as to ensure a weak enough relative uncertainty. In practice, for the lowest bound considered, U/V 0 = 10, the distance crossed at each time step, about 10 pixels, implied a relative accuracy on the trajectory of one half pixel over ten, i.e. of 5%.

Contamination by reaction around a given burnt pixel on a time step ∆t is achieved by changing its neighbors into burnt pixels, whatever their actual state. This procedure is of course redundant in most of the burnt domain where neighbor pixels are already burnt and only efficient at its boundary where it determines the front advance into the fresh domain. It is however simple to handle, time cheaply and provides the advantage of not having to consider the domain boundary, i.e. the front, specifically.

As the spatial and temporal steps are linked by ∆x = V 0 ∆t, the propagation velocity is a priori V 0 . However, as the mapping involves a square lattice, it does not respect the intrinsic isotropy of the actual contamination process by reaction-diffusion. To better approach isotropy on a sufficiently short delay to accurately simulate front propagation, we mix contamination along the lattice axis directions with contamination on the lattice diagonal directions. We thus apply two kinds of mappings, one M a along the axis and the other M d along the diagonal (Fig. 12), each of them occurring on the elementary time step ∆t. One has nevertheless to take into account that contamination along diagonals yields a propagation √ 2 times quicker than along the lattice axis. Therefore, over a period of 3∆t, the mapping M a is applied twice and the mapping M d once, in the following order : M a ⊗ M d ⊗ M a . This way, the difference of propagation with respect to directions is reduced while preserving the propagation velocity V 0 . Altogether, without flow, the iteration of the three mappings yields the contamination from a pixel to be much more isotropic than would be obtained with a single mapping (Fig. 12).

The combined kinematics given by advection and contamination is obtained by applying successively the advection mapping and the contamination mapping, care being taken to alternate the elementary mappings M a , M d , as reported above. The natural change of directions of domain boundaries provided by advection is then found to also improve the isotropy of contamination at large time. The initial condition applied was a line of burnt pixels on the width of the vortex chain, at one of its extremities. It however proved to be unessential for front propagation on distances of several vortices (Fig. 15).

The mean front velocity is accurately measured as the mean velocity of progression of the burnt domain on the chain direction over an integer number of vortices. In practice, the numerical model has been applied from U/V 0 = 10 to 100 for different kinds of flow fields reflecting different b.c. and different vortex aspect ratios d/L. The lowest bound U/V 0 = 10 ensures a sufficient accuracy of the advection scheme and the largest bound U/V 0 = 100 allows the system to remain close to the Hele-Shaw regime (Re = 3 for d = 6 mm and U/V 0 = 100 here).

t 0 ! t 3 ! t 1 ! t 2 !
FIG. 12: (Color online) Propagation mapping starting from a single burnt pixel P at time t0. To approach isotropy on the square lattice, propagation on axes (mapping Ma) and on diagonals (mapping M d ) are alternated at each time steps ti = t0 + i∆t. Mappings at t1 : MaP ; at t2 :

M d ⊗ MaP ; at t3 : Ma ⊗ M d ⊗ MaP .

B. Front propagation

Numerical simulations are worked out for flows corresponding to the fundamental modes for free [START_REF] Saul | Oscillations and travelling Waves in chemical system[END_REF] or rigid (34) lateral b.c.. They qualitatively reproduce the front evolution found in experiment and apparently show quite similar front trajectories for both b.c. and various depths (Fig. 15). See the movie [60] reporting a simulation performed at U/V 0 = 60 and free BC and the movie [61] reporting a simulation performed at U/V 0 = 60, d = 6, L = 20 and rigid BC. We address below the quantitative outcomes of these simulations regarding front velocity. The full points refer to the kinematic simulation worked out here. No dependence with respect to the channel depth d is in order since the stream function [START_REF] Saul | Oscillations and travelling Waves in chemical system[END_REF] of the flow does not involve it. No dependence with respect to the vortex width L is in order either, since the stream function involves spatial variables nondimensionalized by L.

Free lateral boundary conditions

The full line interpolates the data provided by the numerical simulations of the advection-reaction-diffusion dynamics in the same flow and for the same b.c. [12]. Interestingly, both kinds of simulation agree remarkably. This confirms that, in the present regime of advectionreaction-diffusion, the front thickness is not a sensitive parameter and may be overlooked. This also validates the basement of the kinematic simulation, in particular its advection and contamination schemes, as well as the relevance of its basic time scale and resolution. In particular, figure 14 shows the implication of a decrease of the resolution on front propagation below the 1200 pixels per length L used here. Only below a much smaller resolution of 120 pixels per length L, do we notice a significant effect on velocity enhancement.

Finally, one notices the large linearity of the relationship between the front effective velocity V f and the flow velocity U . (34) relevant to rigid b.c.. As the stream function then depends on the ratio d/L through the parameter β, considering different channel depths d for a given vortex width L makes sense. We then observe that the larger d/L, the weaker the effective front velocity V f . This echoes the similar dependence on d found in the experiment for rigid b.c.. However, the effect is weaker here than in the experiment (Fig. 8-b). In particular, at the boundary of the experimental range, U/V 0 = 40, the relative variation of V f with d is about 10% whereas it is about 55% in experiment.

Rigid lateral boundary conditions

The linear fit of data obtained for free b.c. has been reproduced on figure 13-b for comparison. It appears that instead of the linear trend found for free b.c., the enhancement curve is concave : the rigid b.c. lower the effective front velocity V f significantly. Here too, this behavior is similar to that found in experiment, but weaker.

V. DISCUSSION

We analyze the variations of velocity enhancement with the channel depth and the boundary conditions (b.c.) and then comment on the sensitivity of the front propagation to the vortex structure.

A. Variations with the channel depth and the b.c. Figure 16 gathers the evolution of front velocity V f with vortex intensity U obtained in experiment and lattice dynamics simulation for free (Fig. 16-a) and rigid (Fig. 16-b) b.c.. Experiment and simulation agree on the following : for free b.c, on the same linear raise of V f with U , independently of the channel depth d ; for rigid b.c. on a concavity of the relation between V f and U and its dependence on the channel depth d. This concavity means that the change of vortex structure induced by turning the b.c. from free to rigid weakens the relative enhancement of front velocity by flow advection. The independence (resp. dependence) with respect to the channel depth d echoes the same independency (resp. dependency) found analytically for the dominant mode of their stream function (33) (resp. (34)).

This overall agreement between experiment and simulation on the salient features of the velocity enhancement by vortex advection supports the modeling of vortices by planar flows and the analytical solution for their dominant modes derived in section III.

However, for rigid b.c., front velocities appear weaker in simulation than in experiment for the largest channel depth and the largest vortex intensities U/V 0 > 30. This effect may be due to phenomena related to deviations from the Stokes regime and especially to secondary flows (see. Appendix D). These flows are induced by the inertial term of the Navier-Stokes equation and result from the part of it that is not compensated by a pressure gradient. Their sources are actually localized here close the rigid boundaries, either the lateral or the top and bottom boundaries. They begin to be noticeable when the flow departs from the Stokes approximation, i.e., in our velocity range, mainly for d = 6mm and either U/V 0 = 15 at lateral boundaries or U/V 0 = 30 at top and bottom boundaries. They may then generate secondary flows that would break the planar flow assumption and might explain the additional depletion of velocity observed at d = 6 mm in the experiment.

B. Sensitivity of front propagation to vortex structure

Following the large dissipation inherent to Hele-Shaw cells, the change of b.c. confines flow modifications to the vicinity of boundaries. Yet, these local changes suffice to significantly affect the effective front velocity, making its raise with vortex intensity turn form linear to concave and from an independence to a dependence on the channel depth. This expresses a large sensitivity of front propagation to the flow details. This conclusion is also supported by the possible depletion of front velocity induced by secondary flows at the largest flow intensity and channel depth, since their amplitude is nevertheless faint. We address below the nature and the implications of this sensitivity. The origin of this sensibility lies in the large scale separation between the front thickness λ and the vortex size L. Following it, the front interacts with local parts of the vortex and can then probe its inner structure. In particular, during vortex crossing, the front visits some definite zones. Should the flow be modified in them, the front propagation could be noticeably altered whereas the overall flow modification would nevertheless seem tiny. Accordingly, it is the local nature of the front with respect to the vortex that makes it sensitive to localized changes of the flow. This gives rise to the sensitivity of front propagation to the flow structure.

In comparison, the diffusion of passive tracers in the same vortex flows appears to refer to a different nature since diffusion spreads tracer concentrations over distances much longer than the vortex scale. This yields a large homogeneity of tracers within each vortex following which the effective tracer diffusion from vortex to vortex depends on their concentration level in vortices and on the structure of vortex separatrices. Therefore, local changes of the vortex flow elsewhere than on these separatrices would negligibly modify the concentration field in a vortex and thus the diffusion dynamics. The absence of scale separation between the diffusion field of passive tracers and vortices thus prevents their diffusion over the vortex chain to be sensitive to the flow structure.

The implication of the local nature of front propagation on front velocity has been explored in a model of trajectory optimization [16,17]. It has been found that the quickest way for a front to cross a vortex was to be advected on a streamline belonging to the boundary layer but close to the quickest streamline [16,17]. Then, the front could take benefit from both a quick advection from one side of the vortex to the other and a short distance to the next separatrix which has to be crossed at its slow proper velocity V 0 . Of course, should the boundary layer be modified by a change of b.c. or by a raise of channel depth, the location of the quickest streamline and thus of the optimal trajectory would be affected, as well as the resulting effective front velocity.

The features of this quickest trajectory help in understanding that, among the different configurations of rigid b.c., the one which provides the closest velocity enhancement to that of free b.c. corresponds to the thinner channel depth, d = 2 (Fig. 16-b). At first sight, this statement might appear somewhat paradoxical since the most dis- When open flows are added to the vortex array [53], we finally notice that a link between spatial structure and sensitivity may be guessed. In some regimes, the front actually appear frozen in the vortex frame. Its geometry then determines the effective net burning rate and thus appears as the analog of the mean front velocity in closed flows. The particular structure of this so-called burning invariant manifold [54] results from a fine equilibrium be-tween propagation by reaction-diffusion and advection by the flow. Its complex geometry may then be guessed to sensitively depend on the vortex structure, thus providing an analog of the dynamical sensitivity evidence here.

Implications

The sensitivity of large scale front propagation to the small scale structure of vortices questions the objective of previous modelings that aimed at relating the effective front velocity V f to variables at the vortex scale only, e.g. U . This goal consists in providing a universal relationship for describing propagation in vortex arrays, independently of the vortex structure, whereas variations of V f of 40% are displayed here by simply changing the lateral boundary conditions (Fig. 16). In contrast,our study shows that, following the large scale separation between front and vortex, all sub-scales cannot be removed by coarse-graining or by renormalization. Accordingly, models of large scale front propagation should obtain universal relationships, but restricted to classes of vortex flows.

On the practical side, this sensitivity might open the way to strategies for enhancing or lowering front propagation by acting at definite parts of a vortex instead than on its whole. In particular, imposing lateral boundary conditions at a channel depth of d = 6mm and at large vortex amplitude, U/V 0 > 30, proves to be equivalent here to lowering the vortex intensity by nearly a factor 2 (Fig. 16-b). This might find applications on epidemics, fire burning or chemistry in microfluidics for instance.

Regarding multi-scale flows, the sensitivity of large scale propagation to sub-scale flows yields to reconsider the relevance of renormalization procedures. One may expect their validity to be restricted to some classes of imbricated flows, yielding for instance a difference between turbulent flows in which sub-scales permanently evolve and laminar flows where they do not. We seek to address this issue in the near future in the same experiment set-up.

VI. CONCLUSION

We have studied front propagation in a periodic lattice of vortex flows and addressed the dependence on the vortex structure of the effective front velocity at a large scale. For this, we considered three channel depths d and two kinds of boundary conditions (b.c.) : free b.c. in an extended lattice and rigid b.c. in a vortex chain. Front propagation was studied both experimentally using an athermal autocatalytic reaction in a solution stirred by electroconvective flows, and by kinematic simulation based on flow solutions derived in the Stokes regime. Altogether, our combined experimental, analytical and numerical study thus provides the first systematic and comparative analysis of front propagation for different boundary conditions.

Solutions to planar electroconvective flows in the Stokes regime were obtained in term of Fourier series whose fundamental mode appeared to be dominant. It showed no modification of the planar flow structure with the channel depth d for free b.c. but a boundary layer of order d for rigid ones. Accordingly, the flow modifications when changing the b.c. were only local, close to the channel boundaries, and tiny at the vortex scale. Nevertheless, they induced a noticeable reduction of front velocity both in experiment and in simulation according to which the vortex structure is as important as the vortex intensity regarding the effective front propagation. This reveals a sensitivity of front propagation to the fine details of vortices and thus to their sub-scales.

This sensitivity proved to be larger in experiment than in simulation. As the front velocity for free b.c. was well recovered by simulation, we attribute this difference to the development of secondary flows for rigid b.c., follow-ing the raise of inertial terms. The fact that front propagation enabled the detection of secondary flows at low Reynolds numbers, Re < 4, prior to any other noticeable manifestation, emphasizes its sensitivity to flow details.

On a physical ground, this sensitivity relies on the scale difference between the front and the flow structure, the front thickness λ being quite smaller than the vortex scale L. This allows the front to take advantage of tiny flow features at its proper scale λ to enhance its propagation. On the practical side, this sensitivity might be used to diagnose some fine features of flows from the propagation of fronts in them or to enhance or weaken front propagation by tiny modifications of the vortex flow in which fronts propagate. On the fundamental side, it makes a definite difference with passive tracers for which the scale of their diffusion field is comparable to that of vortices.

Regarding scale analysis, the sensitivity of front propagation to sub-scales revealed in this single scale study denies the emergence of a universal solution for front propagation in stirred flows, except when the flow structure is fixed, as found here for free b.c.. On a multi-scale context, where different vortex scales might be imbricated, this feature raises the question as to whether the sensitivity to sub-scales will persist or whether universality will be recovered by statistical effects in a scale cascade.

The linear differential equation (A1) admits a base of solution involving products of complex exponentials and, because of possible secularity, of monomials : ψ c ≡ e k.r z e y f x g , k ∈ C 3 , r = (x, y, z), (e, f, g) ∈ R 3 . However, the sinusoidal dependence on x involved here requires k x = iγ, γ ∈ R and g = 0. Calling for real solutions then yields the corresponding basic modes to be products of circular functions, real exponentials or monomials. We enlarge this quest below by seeking basic solutions involving separate variables :

ψ c = f (z)g(y) cos(γx), ; γ = (2k + 1)π/L k ∈ N (A2)
As equation (A1) involves a composition of operators, secular terms in f (.) and g(.) are expected. Finally, following the symmetry of our set-up with respect to the (x, z) plane, attention will be restricted to even functions g(.).

Application of the Laplacian operator yields :

h ψ = f (z) cos(γx)G (A3) G = g -γ 2 g (A4) h ψ = [(f -γ 2 f )G + f G ] cos(γx) (A5)
Equation (A1) then provides the following constraint on functions f (.) and G(.) :

(f -γ 2 f ) G = -f G (A6)
which yields the following alternatives :

1. G = 0, ∀f Equation (A4) then selects g(y) = cosh γy

G = δ 2 G f = 2 f with 2 = γ 2 -δ 2 .
The solutions read :

i. = 0 Then f (.) = a + bz and G(y) = cosh(γy).

Accordingly, G appears as a resonant forcing term in relation (A4). It thus generates a secular term g(y) = y sinh(γy).

ii. = 0 One obtains f (z) = e ± z and G(y) = cosh(δy). As δ = γ, G appears as a nonresonant forcing term in relation (A4).

It then generates a normal mode g(y) = cosh(δy). Depending on the sign of δ 2 = ±β 2 and 2 = ±α 2 , either circular or exponential functions are generated. Altogether, one finally obtains the following base of complementary functions :

ψ 1 = f 1 (z) cosh(γy) cos(γx) ψ 2 = (a + bz) y sinh(γy) cos(γx) ψ 3 = e ±γz (a + by) cos(γx) ψ 4 = e αz cosh(βy) cos(γx) with α 2 + β 2 = γ 2 ψ 5 = e αz cos(βy) cos(γx) with α 2 = β 2 + γ 2 ψ 6 = cos(αz + ϕ) cosh(βy) cos(γx) with β 2 = α 2 + γ 2
where the wavenumbers (α, β, γ), the coefficients (a, b), the phase ϕ and the function f 1 (.) are real.

Following (15) and (A1), all these stream functions refer to no vorticity sources. In addition, as the stream function ψ 1 satisfies h ψ 1 = 0, it refers no vertical vorticity. On the opposite, as the remaining stream functions ψ i , 2 ≤ i ≤ 6, satisfy h ψ i = 0, they generate some vertical vorticity which corresponds to that induced by the rigid lateral boundaries. We seek a Fourier expansion of the Fourier amplitude ψ f k,p (z) (27) which is defined in the interval [-d/2, d/2]. For this, we wish to prolongate it so as to make it periodic prior to expanding it in a Fourier series.

As ψ f k,p (z) vanishes at z = ±d/2, d may be considered as a multiple of its half-period. This allows us to prolongate it over an interval of length 2d and then to duplicate it so as to make it periodic. This is obtained by considering the function ψf k,p (z) that is odd both in [-d, 0] with respect to -d/2 and in [0, d] with respect to d/2 and that is equal to ψ f k,p (z) in [-d/2, d/2]. The vanishing and the odd symmetry of function ψf k,p (z) at z = ±d/2 ensures both its continuity and that of its derivative, as well as the vanishing of its average on [-d, d]. Being C 1 , ψf k,p (z) will then be equal to its Fourier series in ] -d, d[ and thus in [-d/2, d/2], where it reduces to ψ f k,p (z). The corresponding Fourier modes ψ f n,k,p (x, y, z) are made explicit in relation (28). They then yield the Fourier modes S z n,k,p (.) (30) of the vorticity source S z (11) (14). We then notice that, as the modes ψ f n,k,p (.) or ψ r n,k,p (.), all the modes S z n,k,p (.) vanish at z = ±d/2 whereas the vorticity source S z does not. This discrepancy comes from the fact that the Fourier series based on the modes S z n,k,p (.) corresponds to that of the vertical vorticity source S z prolongated in the interval [-d, d] by the same method of odd symmetry than that applied to the stream function. However, as S z does not vanish at z = ±d/2, its prolongation involves a discontinuity at these points. This forbids equality between S z and its Fourier series at these discontinuity points. where µ = l/L. Considering, in agreement with the experimental configuration, a channel width l containing an odd number of vortices (or equivalently, to a phase shift, an odd number of magnets), we get µ = 1 + 2q, q ∈ N .

Noticing that sin[(2p + 1)π/2 µ] = (-1) p+q , the solution to the algebraic system (C3) yields the stream function ψ Notice that, in the combination yielding this stream function, ψ f n,k,p (.) stands as the solution to the vorticity source for free lateral b.c. whereas ψ 1 and ψ 6 refer to solutions to no vorticity source. Among them, ψ 1 induces no vertical vorticity since h ψ 1 = 0. In contrast, ψ 6 which satisfies ψ 6 = 0, generates some vertical vorticity equal toh ψ 6 = ∂ 2 ψ 6 /∂z 2 and thus to -n 2 (π/d) 2 ψ 6 . This corresponds to the vorticity induced by the rigid lateral boundaries.

APPENDIX D: SECONDARY FLOWS

Secondary flows are induced by the irrotational part of the inertial term of the Navier-Stokes equation, i.e. by the part that can therefore not be compensated by the pressure gradient. Their sources appear close to boundaries since the velocity shrinks to zero rapidly whereas the pressure gradient remains constant, thus yielding a sudden disequilibrium in the Navier-Stokes equation. Their intensity may thus be expected to be noticeable when inertial terms get no longer negligible with respect to dissipative ones, i.e. when the Stokes approximation begins to be invalid. This criterion corresponds to a Reynolds number reaching a value of order unity, the appropriate length and velocity scales used to built it depending on the flow zones that are addressed. As evoked in section III, two kind of zones need to be distinguished : the vicinity of lateral boundaries and the bulk. Close to lateral boundaries, i.e. at a distance of the order of the channel depth d, all characteristic flow scales are of order d. Then the Reynolds number writes Re = U d/ν so that, with ν = 1 mm 2 /s, it turns out to be unity for : U = 30 mm/mn at d = 2 mm, U = 15 mm/mn at d = 4 mm and U = 10 mm/mn at d = 6 mm. In our experimental range U ≤ 50 mm/mn, a departure from the Stokes regime may thus be feared close to boundaries. However, to our knowledge, no detailed study of its implication has been achieved in our configuration made of straight lateral boundaries and circular streamlines. Nevertheless, a detailed study in an analogous configuration involving a uniform flow encountering a bumpy boundary in a Hele-Shaw cell revealed the generation of secondary flows localized in the vicinity of the bump [55]. By analogy, one may expect here secondary flows to be induced within the boundary layers of the lateral boundaries, especially at d = 6mm and for values of U/V 0 of several tenths.

In the bulk, the Stokes regime ends when the boundary layers induced by the horizontal plates separate. This corresponds to a transition from a Poiseuille profile to a flatter profile. As the scale of variation is L in the horizontal plane (x, y) and d on the direction z, one gets |V.∇V| ≈ U 2 /L, |ν V| ≈ νU/d 2 and the Reynolds number Re = Re d/L = U d 2 /Lν. It then reaches unity for : U = 300mm/mn for d = 2mm, U = 75mm/mn for d = 4mm and U = 33mm/mn for d = 6mm. Then, an inviscid zone progressively emerges around the mid-depth of the channel. Inertial terms are then no longer equili-brated by dissipative terms but by a pressure gradient. However, in the vicinity of top and bottom boundaries, the sharp decrease of velocity and thus of dissipation, together with the quasi constancy of pressure, yields an imbalance which generates a secondary flow. This mechanism, analogous to the Ekman pumping [56], generates vertical flows which may perturb front propagation. In our velocity range, these secondary flows are too weak at small or moderate channel depths, d = 2 or 4 mm, but may be noticeable at the largest channel depth d = 6mm above about U = 30mm/mn. They may thus partly explain the specificity of this channel depth and the noticeable departure of experimental data from numerical ones beyond U/V 0 ≈ 30.

FIG. 1 :

 1 FIG. 1: Sketch of the set-up showing the chemical solution stirred by Laplace forces induced by an electric current flowing over alternate magnets (a) and the resulting array (b) or chain (c) of counter-rotating vortices. The channel depth, the vortex size and the channel width are denoted d, L and l respectively. The magnetic field and the electric current are labelled B and I.

FIG. 2 :

 2 FIG. 2: (Color online) Front propagation with free lateral boundary conditions, except for the side vortices. Times regularly increases from left to right while fronts propagate upwards. The channel depth is d = 2mm, the vortex width L is 20mm and the channel width is labelled l : l/L = 5, U/V0 = 33.0. Time interval between pictures : 49s.

- 1 L

 1 = ν/LU of the Reynolds number Re L based on the vortex width L, ν denoting the fluid kinematic viscosity ; an effective damping coefficient µ that takes into account the viscous friction implied on the planar flows by the top and bottom plates of the channel. For a Poiseuille profile in the channel depth direction, one gets µ = 8Re -1 L (L/d) 2 which is quite larger than Re -1 L here : µ = 800Re -1 L for d = 2mm and µ ≈ 90Re -1

FIG. 4 :

 4 FIG. 4: (Color online) Flow field measured by PIV in a channel with rigid lateral boundary conditions. The channel depth is d = 4mm, the vortex width L is 20mm and the channel width l equals the vortex width L : l = L = 20 mm. (a) Velocity field. (b) Superimposed intensity of the velocity field intensity obtained by extrapolation techniques. The vanishing of the flow field at the lateral boundaries makes a difference of vortex structure with free boundary conditions (Fig.5).

FIG. 5 :

 5 FIG. 5: (Color online) Flow field measured by PIV in a channel with free lateral boundary conditions. The channel depth is d = 4mm, the vortex width L is 20mm and the channel width l amounts to 6L. The slight distortion at the left bottom originates from an uncompensated defect in the camera frame. (a) Velocity field. (b) Superimposed intensity of the velocity field intensity obtained by extrapolation techniques.

FIG. 6 :

 6 FIG. 6: (Color online) Front propagation for free b.c.. The vortex width is L = 20mm. Columns correspond to channel depths d = 2, 4, 6 mm. Lines correspond to reduced vortex intensities U/V0 weak (≈ 10), moderate (≈ 20) and large (≈ 30).

FIG. 7 :

 7 FIG. 7: (Color online) Front propagation for rigid b.c.. The vortex width is L = 20mm. Columns correspond to channel depths d = 2, 4, 6 mm. Lines correspond to reduced vortex intensities U/V0 weak (≈ 10), moderate (≈ 20) and large (≈ 30).

FIG. 8 :

 8 FIG. 8: (Color online) Reduced front velocity V f /V0 as a function of the reduced vortex intensity U/V0. (a) Free b.c. : a single linear relationship is displayed whatever the channel depth d. (b) Rigid b.c. : concave relationships are displayed with an increasing concavity with the channel depth d.

  γy) cos(γx) ψ 2 = (a + bz) y sinh(γy) cos(γx) ψ 3 = e ±γz (a + by) cos(γx) ψ 4 = e αz cosh(βy) cos(γx) with α 2 + β 2 = γ 2 ψ 5 = e αz cos(βy) cos(γx) with α 2 = β 2 + γ 2 ψ 6 = cos(αz + ϕ) cosh(βy) cos(γx) with β 2 = α 2 + γ 2 where the wavenumbers (α, β, γ), the coefficients (a, b), the phase ϕ and the function f 1 (.) are real.

  As the free lateral b.c.(20) (21) are naturally satisfied by sinusoidal factors, the particular solution ψ p satisfies them both on the x and y directions. It thus only remains to satisfy the rigid b.c.(17) on the z direction : ψ(x, y, ±d/2) = 0.

  imposed on the y axes by rigid b.c., we now seek to complement the solution (26) by complementary functions. We thus start from the stream function ψ f (26) (27) which satisfies rigid b.c. on the top and bottom plates z = ±d/2 but free lateral b.c. at y = ±l/2 and we seek to expand its Fourier amplitude ψ f k,p (z) in a Fourier expansion in [-d/2, d/2] based on the modes ψ f n,k,p (x, y, z), n ≥ 1, k ≥ 0, p ≥ 0. This requires prolongating the function ψ f k,p (z) by a periodic function ψf

FIG. 9 :

 9 FIG. 9: (Color online) Profile of the fundamental mode ψ f 0,0 on the channel depth direction z (blue full line) and comparison with a sinusoidal shape (a) (red dashed line) and a parabolic shape (b) (red dashed line).

  ψ cn,k,p (.) satisfy the rigid lateral b.c. at y = ±l/2. Among the set of complementary functions involving separate variables, only those exhibiting the same dependance on z than ψ f n,k,p (.), i.e. a sinusoidal dependance, can be relevant to our issue. Following appendix A and section III A 3, only two kinds of them satisfy this constraint, the functions ψ 1 and ψ 6 , with the respective dependance on y : cosh[(2k + 1)πy/L] and cosh[βπy/L] withβ 2 = (2k + 1) 2 + n 2 (L/d) 2 .Considering the linear combination ψ r n,k,p = ψ f n,k,p + a ψ 1 + b ψ 6 with real coefficients (a, b), provides two degrees of freedom to satisfy the rigid lateral b.c. (18) (19) : ψ r n,k,p (x, ±l/2, z) = 0 ; ∂ψ r n,k,p /∂y(x, ±l/2, z) = 0. This yields a linear algebraic system whose solution is reported in appendix C. Labelling µ = l/L, it provides the following stream function : ψ r n,k,p (x, y, z) = ψ f n,k,p cos[nπz/d + ϕ n ] cos[(2k + 1)πx/L] × cos[(2p + 1)πy/L] + (-1) p (2p + 1) ψn,k,µ (y/L) (29) The difference with free b.c. is thus conveyed by the uni-variate function ψn,k,µ (y/L) whose expression is reported in appendix C. When the channel width is taken large compared to the magnet width, i.e. when µ 1, the function ψn,k,µ (.) decreases as e -(2k+1)π|l/2-y|/L . As expected, it thus only serves close to the lateral faces y = ±l/2 to satisfy the rigid b.c. and leaves in the bulk a stream function close to the expression found for free b.c.. The connexion of Fourier modes to the vorticity source S is made by the Stokes vorticity equation (15). Both the modes ψ f n,k,p (.) and ψ r n,k,p (.) are excited by the same Fourier mode S z n,k,p (.) of the vertical component S z of S : S z n,k,p (x, y, z) = S z n,k,p cos[nπz/d + ϕ n ] cos[(2k + 1)πx/L] cos[(2p + 1)πy/L]

FIG. 10 :FIG. 11 :

 1011 FIG. 10: (Color online) Stream function of the fundamental mode (n, k, p) = (1, 0, 0) for free b.c. : contour plot of ψ f1,0,0 in the mid-depth plane (x, y, 0). Lines correspond to increased level of ψ f of amplitudes i/8, 0 ≤ i ≤ 8.

Figure 13 -

 13 Figure 13-a shows the enhancement of the effective front velocity V f by the flow for free b.c.. The full points refer to the kinematic simulation worked out here. No dependence with respect to the channel depth d is in order since the stream function

Figure 13 -

 13 Figure13-b displays the results of the simulation of front propagation for the fundamental flow mode(34) relevant to rigid b.c.. As the stream function then depends on the ratio d/L through the parameter β, considering different channel depths d for a given vortex width L makes sense. We then observe that the larger d/L, the weaker the effective front velocity V f . This echoes the similar dependence on d found in the experiment for rigid b.c.. However, the effect is weaker here than in the experiment (Fig.8-b). In particular, at the boundary of the experimental range, U/V 0 = 40, the relative variation of V f with d is about 10% whereas it is about 55% in experiment.

FIG. 13 :

 13 FIG. 13: (Color online) Kinematic simulation of front propagation for the fundamental modes (33) (34). Front velocity enhancement for : (a) Free b.c.. Full points correspond the present kinematic simulation of front propagation. The full line interpolates data obtained from simulation of advection-reaction-diffusion dynamics for the same flow and the same b.c. [12]. (b) Rigid b.c.. The full line, which refers to a fit of data for free b.c. (a), is reported for comparison.

FIG. 15 :FIG. 16 :

 1516 FIG. 15: (Color online) Simulation of front propagation for free or rigid b.c.. Burnt (resp. fresh) medium appears in white (resp. in blue (black)). Columns correspond to free b.c. (left), rigid b.c. with d = 2 mm, L = 20 mm (center) and rigid b.c. with d = 6 mm, L = 20 (right). Lines correspond to reduced vortex intensity U/V0 = 20 (above) and U/V0 = 60 (below). The differences are tiny and mainly concern the distance between the front and the lateral boundaries. These distances are similar for free b.c. and rigid b.c. with d = 2mm, but increase for d = 6mm.

2. G = 0 (

 0 a) G = 0 Function g(.) writes g(y) = a + by. Function f (.) then satisfies f = γ 2 f and thus writes f (z) = e ±γz . (b) G = 0 As f (.) and G(.) are functions of different vari-ables, there exists a factor δ 2 , positive or negative, such that :

  APPENDIX B: FOURIER EXPANSION OF ψ f k,p(z)

  APPENDIX C: STOKES SOLUTION FOR RIGID BOUNDARY CONDITIONS We consider the linear combination ψ r n,k,p = ψ f n,k,p + a ψ 1 + b ψ 6 with real coefficients (a, b) and seek which of their values yields the rigid lateral b.c. (18) (19) to be satisfied. These b.c. write : ψ r n,k,p (x, ±l/2, z) = 0 ; ∂ψ r n,k,p ∂y (x, ±l/2, z) = 0 (C1) Satisfying them yields the linear algebraic system : a cosh[(2k + 1)π/2 µ] + b cosh[βπ/2 µ] = 0 (C2) -(2p + 1) sin[(2p + 1)π/2 µ] + a(2k + 1) sinh[(2k + 1)π/2 µ] + bβ sinh[βπ/2 µ] = 0 (C3)

  r n,k,p :ψ r n,k,p (x, y, z) = ψ f n,k,p cos[nπz/d + ϕ n ] cos[(2k + 1)πx/L] × cos[(2p + 1)πy/L] + (-1) p (2p + 1) ψn,k,µ (y/L) (C4)where the difference with the solution for free b.c. is provided by the function ψn,k,µ (y/L) :ψn,k,µ (y/L) = (-1) q cosh[βπ µ/2] cosh[(2k + 1)π y/L] -cosh[βπ y/L] cosh[(2k + 1)π µ/2] (2k + 1) cosh[βπ µ/2] sinh[(2k + 1)π µ/2] -β sinh[βπ µ/2] cosh[(2k + 1)π µ/2] (C5)

Acknowledgments

We acknowledge useful advices from M. Amielh, L. Bergougnoux, B. Metzger, P. Meunier for achieving PIV of the flow field in the Hele-Shaw channels. The project leading to this publication has received funding from Excellence Initiative of Aix-Marseille University -A*MIDEX, a French "Investissements d'Avenir" programme. It has been carried out in the framework of the Labex MEC ("Mechanics and Complexity").

APPENDIX A: COMPLEMENTARY FUNCTIONS

We look to solutions to the unforced Stokes vorticity equation (A1) :

These solutions will be used to complement the particular solution ψ p (22) (24) so as to satisfy rigid b.c. on the directions z and y. However, in the present configuration, the periodic forcing of the magnet lattice in a channel whose ends are overlooked yields the expected stream functions ψ for both free and rigid b.c. to be periodic in x, with period 2L. As the complementary functions ψ c = ψ r -ψ, they will therefore involve a Fourier expansion involving sinusoidal modes on x, with wavenumbers (2k + 1)π/L, k ∈ N . Following the linearity of equation (A1), we shall then determine the complementary function for each of these Fourier modes independently.