Towards personalized medicine : chemosensitivity assays of patient lung cancer cell spheroids in a perfused microfluidic platform †

3076 | Lab Chip, 2015, 15, 3076–3085 This journal is © The R a ARTORG Lung Regeneration Technologies Lab, University of Bern, Murtenstrasse 50, CH-3010 Bern, Switzerland. E-mail: olivier.guenat@artorg.unibe.ch; Fax: +41 31 632 7576; Tel: +41 31 632 7608 Graduate School for Cellular and Biomedical Sciences, University of Bern, Switzerland c Division of General Thoracic Surgery, University Hospital of Bern, Bern, Switzerland Department of Clinical Research, University of Bern, Bern, Switzerland e Division of Pulmonary Medicine, University Hospital of Bern, Bern, Switzerland † Electronic supplementary information (ESI) available: Video showing spheroid formation. See DOI: 10.1039/c5lc00454c Cite this: Lab Chip, 2015, 15, 3076


Introduction
Microbubbles are powerful active elements, 1 and have been used in a variety of applied setups, including particle transport, 2 manipulation 3 and for mixing purposes. [4][5][6][7] These applications harness the powerful flows generated by remotely excited 4,5,[7][8][9][10][11] microbubbles. However, one weakness of microbubbles is their short life span. Indeed a free micron-sized air bubble dissolves in water at room temperature in tens of seconds. In a previous study, 7 we overcame this issue by proposing a new type of acoustic micropropulsor, namely armoured microbubbles (AMBs). An AMB is composed of a capsule, which is a 3D printed hollow partial sphere fabricated on top of a short pole, inside which a microbubble can be caught when immersed in water (see sketch in Fig. 1a and experimental image in Fig. 2a). Recently, stereolithographic 3D printers (also called SLA printers) were used to build channels, valves and pumps, 12,13 and to fabricate electrospray sources in millifluidics. 14 In order to print the capsules, our setup uses 2-photon absorption lithography to reach submicron fabrication precision, allowing us to manufacture capsules on the 10 μm scale. We fabricate capsules with opening diameter at most half of the capsule's overall diameter, due to shorter bubble lifespan at larger opening diameters. Under ultrasound, the free surface of the bubble near the opening of the capsule oscillates, generating a powerful flow via acoustic streaming. We can predict this flow with a detailed 3D hydrodynamic model and also predict the resonance frequencies of the bubbles. 7 Having previously performed a detailed analysis for an isolated AMB, 7 we consider in this article the variety of flows achievable by multiple AMBs. We first consider multipropulsor compounds, where multiple AMBs are physically attached together, and then address active surfaces, where multiple AMBs are placed in a pattern short distances apart. Virtually any arrangement of AMBs may be fabricated using our method and we demonstrate the variety of flows achievable. We consider below specific AMB compounds such as double propulsors, which produces a wider jet than a single AMB, and triple propulsors, which produce multi-directional flows. Moving on to active surfaces, we see that the flows they produce can display large vortices which, through Taylor-Aris dispersion, 15,16 are able to enhance significantly effective diffusion. This feature, combined with the remote actuation of an active surface, its high adaptability and minimal impact on the base flow when the ultrasound is switched off, make active surfaces a good candidate for a microfluidics mixer.

Fabrication of capsules and setup
We use a 2-photon absorption microscopic device (Microlight 17 ) to fabricate hollow capsules of inner diameter 2r with an opening of diameter 2a (see Fig. 1a). The resin is OrmoComp®, 18 a biocompatible hybrid polymer used for ultraviolet (UV) imprint and moulding. We prepare a solution of 6 mg 1,3,5-trisĲ2-(9-ethylcabazyl-3)ethylene)benzene photoinitiator dissolved in 1 ml dichloromethane, and mixed with 1 g OrmoComp®. The resin is only polymerised in the vicinity of the beamwaist of the laser, and the precision is on the order of the 500 μW Nd:YAG microchip laser wavelength, 532  nm. 19 The laser writing speed is estimated at 50 μm s −1 , with 2 ms laser pulses. The setup is mounted on the epifluorescence port of an inverted microscope.
In the measurements shown below, all capsules have an internal radius r = 9 μm, and a double wall construction. The distance between the two walls is 1 μm so the external radius is 10 μm. The opening radius is a = 5 μm unless otherwise stated. Indeed, we have discovered that the optimum ratio for increasing the bubble lifespan is a/r ≈ 0.5. We build AMBs with 2 μm radius poles of height H = 10 μm, which provide more stability during fabrication than the H = 30 μm poles considered in earlier work. 7 This is especially important for our heavier multipropulsor compounds. To overcome the difficulty inherent to printing multiple objects, we added an autofocus capability to the setup. The tower and capsule are printed onto a glass coverslip used as a substrate, but at the micron scale, this surface will not be completely flat. The object must be printed on the substrate, and not above, else it will be washed away during rinsing. Therefore at the start of fabrication of a new object, the beamwaist of the laser should be positioned at or just below the surface of the glass coverslip. A lateral movement of the piezoelectric stage of 100 μm may result in a vertical shift of approximately 0.5 to 1 μm or a tilt angle of the coverslip of 0.3 to 0.6°. We automatically correct for any vertical discrepancy using our autofocus script. It analyses at low laser power the light intensity transmitted at various altitudes, in order to find the exact position of the laser beamwaist relative to the coverslip.
In order to allow us to print a large number of objects in predetermined patterns such as spirals (Fig. 1c), circles (Fig. 1d), and triangles (Fig. 3a), we program a long-range translation stage for a list of positions using a Python script. Since changing the direction of an AMB only requires rotating it around the z axis (see notation in Fig. 1a), our patterns include capsules with openings pointing in all possible directions.
The fabrication process is autonomous once the printing area is defined. A propulsor takes approximately 15 min to print, so 100 propulsors can be fabricated in 28 hours. The only time limitation is due to the resin, which naturally cures over about 48 hours. Scanning Electron Microscope (SEM) observations of the fabricated objects confirmed the reproducibility of the process (see Fig. 1b-d).

Acoustic streaming flow prediction around individual propulsors and AMB arrays
We have previously derived a theoretical model for calculating the streaming flow around a single AMB both in free space and near a wall. 7,20 In order to approximate the streaming flow induced by double and triple propulsors, we linearly superpose the streaming flows each of the AMBs produce individually when isolated and close to a wall. Since we are superposing individual solutions which only obey the boundary conditions on their own AMB, this introduces error and the resulting streaming flow approximation is not valid close to where the AMBs meet. Additional errors are added close to the AMBs due to the importance of inertia in the boundary layer close to the AMBs. But outside the boundary layer the streaming flow is a Stokesflow and so a short distance from the AMBs our approximation is valid.
We use a different approximation to calculate the streaming flow around AMB arrays. Indeed, for AMB arrays the length scales of configurations are O(100 μm) as opposed to O(10 μm) for individual AMBs. Therefore we may approximate the flow field induced by each AMB in free space as its leading-order solution in the far field, 7 namely a Stokeslet (point force), plus its next-order correction of a Stresslet. The Stokeslet strength for an AMB f t was chosen so the flow velocity ≈10 μm from the Stokeslet was of the order of 1-100 mm s −1 , the range of velocities observed in front of an AMB experimentally. The flow field of a Stresslet is determined by a tensor, which for an AMB pointing in the x direction modeled at leading order by a Stokeslet of strength f t in the x direction, is given by (1) This form and strength of the Stresslet is given by the full solution of the flow around the AMB. 7,20 Both the Stokeslet and Stresslet are singular solutions; in order to make the flow field finite everywhere we then mathematically regularise the Stokeslet and Stresslet. Additionally, to model the AMBs near one wall we add the images of the Stokeslet and Stresslet. More explicitly, near a flat plate each AMB is modelled using the regularised version of the well known Blake solution for a Stokeslet above a no-slip wall 21-23 plus the Stresslet with its known solution near one wall. 22,23 We then linearly superpose the flow fields generated by each AMB individually near one wall to obtain the flow field induced by the whole AMB array. If there is a base flow present, we linearly add that contribution to the AMB array flow field.
To model AMB arrays in a channel (i.e. sandwiched between two close walls and with two side walls), we use a slightly different approach. We approximate each AMB by its leading order Stokeslet of strength f t only (so we ignore the Stresslet correction). Between two walls a Stokeslet has a well known solution. 24 To obtain the flow field in a vertical cross section of the channel, we linearly superpose this individual solution for each AMB between two walls.
We then use the fact that the solution for a Stokeslet parallel to two confining walls has a vertical component of velocity which decays exponentially, whereas the planar component of the velocity decays as a power law. Therefore in the far field, the vertical velocity component is negligible. The horizontal velocity is given by (2) where z = 0, h is the bottom and top plate respectively, and is a 2D Stokes Doublet which, for an AMB pointing in the x direction, is where j = 1, 2, μ is the dynamic viscosity of the fluid, ( y 1 , y 2 , H) is the position of the Stokeslet, is the distance in the plane from the Stokeslet, and r 1 = (x − y 1 ), r 2 = (y − y 2 ) is the straight line distances from the Stokeslet in the x and y direction. We see that in the far field the flow is equal to a twodimensional Stokes Doublet as in eqn (3) with a magnitude which varies quadratically across the channel as in eqn (2). 24 We can then regularise this Stokes Doublet to remove singularities using the technique described in ref. 22 and 25. We then linearly superpose this individual solution for each AMB to obtain the total flow field around the AMB array. Note that the result in eqn (2) shows that the choice of horizontal plane (i.e. the value of z) only affects the total magnitude of the flow field thus the flow field is structurally the same in every plane (at leading order). When there is a base flow of strength 4v 0 /h 2 (for constant v 0 ), we assume it is a traditional 2D Poiseuille flow with velocity field , for pointing down the channel, and thus has the same quadratic dependence across the channel as the flow generated by the AMB array, so when we linear superpose it the flow is still structurally the same in each plane. We note that this model ignores the effects of the side walls since they are 1 mm apart as opposed to 80 μm for the channel height, an accurate assumption when the AMB arrays are placed away from these side walls and the AMBs are weak enough for their effects to not extend across the channel width.

Theoretical mixing measure of AMB arrays
Most theoretical studies of mixing use detailed numerical simulations 26,27 or scaling arguments. 28,29 Here we present an alternative approach, namely a numerical method applied to a simplified mixer-channel setup. The advantage of our method is that, due to its simplicity, it can analyse hundreds of mixers a day, and is designed to identify groups of good mixers which can then undergo further detailed numerical and experimental study. We make use of the far-field flow field outlined above for the mixer between two plane walls. Since every plane is structurally the same we focus our analysis on one plane. The regularisation of the Stokes Doublet breaks into a source and sink, but by choosing a regularisation length of half the radius of the AMB these are contained well within the AMB. To prevent fluid unphysically becoming trapped in a sink, we add a Ae −B/(r−r 0 ) 2 velocity contribution in the direction of the Stokes Doublet within the AMB radius r 0 (where A, B are constants) to push flow from the sink to the source. The maximum velocity within the AMB is maintained at the same order. This contribution is continuous at all derivatives on the boundary of the AMB and is added only within the AMB. The AMB strength was taken as 0.45 nN with a base flow in the AMB plane of 1.2 mm s −1 .
Diffusion will ultimately cause the mixing but its effectiveness will be dependent on the size of the surface area between the two fluids which the mixer will increase. In our plane we consider a section of line 3.2 mm in length which divides the two fluids. We measure mixing capability with the ratio of the length of line after it has partially passed through the mixer to its initial length. We take our stretch of line starting away from the mixer where the base flow dominates, and thus we can assume that each point on the line follows the same path, only with a time delay. The analysis is terminated when the end of this line element reaches the starting position of the front of the line element, at which point the front of the line element will have well passed out the other side of the mixer (if not caught in a mixer circulation).

Individual propulsors: single, doublets or triplets
Bubbles resonate at a wavelength much greater than their radius. 30 In the case studied here, the r = 9 μm AMB with an a = 5 μm opening resonates at a frequency of f r = 320 kHz, 7 which corresponds to an acoustic wavelength of λ = 4.7 mm in water. The resonance frequency for a = 3 μm increases to f r = 510 kHz and for a = 7.5 μm it decreases to f r = 160 kHz. We verified experimentally by measuring the resonance frequency that, at constant a and increasing volume, f r is indeed constant. Since λ ≫ r, we are able to excite a large number of propulsors with the same focused ultrasonic wave (55 propulsors on Fig. 3), the acoustic beam width being on the order of λ. In Fig. 2a, the flow generated by a H = 10 μm capsule is illustrated, showing two pairs of counter-rotating vortices. The backward, slower pair of vortices is due to a boundary effect absent when H = 30 μm. 7 Such a propulsor is able to generate forward flows with mean velocity v mean = 1-100 mm s −1 .
We next take advantage of the combined flow from different AMBs and design multipropulsor compounds. In order to extend the lateral range of the forward flow of a single AMB, we first developed a double propulsor. It consists of 2 solidly bound, identical, adjacent capsules, mounted on the same 10 μm pole. The width of the generated forward liquid jet is doubled as a result, see Fig. 2b (with similar jet widening also observed at a = 7.5).
We then considered what multipropulsor compound could generate vortices in all directions in the x-y plane, similarly to a free 2D bubble in a microchannel, 8,31 but protected from dissolution by the capsules. To this end, we developed a triple propulsor consisting of 3 capsules back to back, with openings oriented at 120°from one another generating 3 pairs of vortices (see Fig. 2c). Note that these multi-propulsor compounds have the advantage of having the same resonant frequency as each of the individual propulsors composing it. Both the double and triple propulsor produce flows that can be predicted accurately with our streaming model as shown on Fig. 2.
The AMB and multipropulsor compounds are single units that can be duplicated on a surface, in any pattern, using our new automated fabrication setup, as we now explore.

Arrays of propulsors: active surface
As a first test of active surfaces, we consider a triangular array of 55 AMBs in 10 lines, staggered so that no propulsor is being blocked by another one in front (see Fig. 3a). A polydimethylsiloxane (PDMS) circular ring (2 cm inner diameter) is glued on the coverslip on which the AMBs are printed. It is filled with phosphate-buffered saline (PBS) and 2 μm tracer beads before being closed with a PDMS disk. PBS was used since it increases the lifespan of the bubbles over that in pure water. 7 We spaced propulsors 100 μm apart (from those adjacent on the same line and the closest in the line in front), thus producing an equilateral triangle array of side length 1 mm. We can activate all the propulsors in the triangle simultaneously using the same ultrasound wave since its wavelength λ is approximately 5 times the size of the triangle, as shown by the particle tracking velocimetry (PTV) analysis in Fig. 3b (we used 2 μm beads and FAST software 32 for this).
We can see a net flow generated by this triangular array in Fig. 3a. Each propulsor produces a pair of counter-rotating vortices similar to those seen on Fig. 2a. However, there are unexpected non-uniformities in the flow. These are not due to the fabrication process but may originate from the tracer beads. Indeed we observed that the tracer beads tend to aggregate on the air/liquid interface around the opening, which affects the efficiency of the acoustic streaming and can almost cancel out the flow generation in some instances. Using smaller beads lead to no discernible differences. Another cause may be the non-uniformity of the acoustic field. Although the wavelength is significantly greater than the array size, there may be some modification of the acoustic field due to it crossing the PDMS. Therefore each AMB may not experience exactly the same acoustic excitation. This has been quantified on the histogram in Fig. 3c which shows the disparity of average velocities in front of the AMBs.
We are able to print 0.43 mm 2 of propulsors, all containing bubbles that can be activated by the same acoustic wave. We now use these active surfaces, with AMBs placed closer than 100 μm, in order to generate collective flows for mixing inside microchannels.
4 Active surface application: microchannel mixing

Mixing at small scales
Mixing in microchannels is known to be a challenging problem. At low Reynolds number, turbulence is not available and thus diffusion is solely due to molecular diffusion. 33 There has traditionally been two distinct types of micromixers: passive vs. active mixers, which exploit different strategies to generate mixing. Passive mixers use the shape of the microchannel, 34 either through carved grooves in the floor, 35 or through sophisticated shape designs. 36 In contrast, active mixers use external forces such as a magnetic field to activate a moving part inside the channel. 37 Our active surfaces fall in the category of active mixers, since the acoustic field generates flow through the streaming of AMBs. 4,5,7 The advantage of using individual micropropulsors, as we are now going to demonstrate, is the versatility of the arrangements able to mix. We now explore the effectiveness of different array shapes for mixing two fluid streams in a channel. One of our objectives is to use the minimum number of capsules for the maximum mixing efficiency.

Setup
We use an 80 μm high microchannel of lateral width e = 1 mm (see Fig. 4a and b), with two entrances: one is fed with a pure PBS solution while the other is fed with PBS + blue dye solution (see Fig. 4a). The flow rate is controlled through two syringe pumps and was in the range Q = 2-4 μl min −1 .
Experiments were performed with a 350 kHz focused transducer (Olympus), in the acoustic pressure range of P ac = 100-2000 kPa. Pressure measurements are done using results from a separate experiment where a hydrophone (Onda) is placed behind a ring of PDMS of similar thickness to the microchannel in order to simulate the loss of pressure from transmitting sound from water to PDMS to PBS. The voltage applied to the transducer during the mixing experiments is then converted to the acoustic pressure evaluation deduced from the hydrophone measurements.
In order to assemble the microchannel with the glass coverslip on which the capsules are printed, marker lines are drawn on the glass on each side of the array of capsules. The microchannel is carefully placed on top after plasma treat-ment, with the pattern chosen from the possibilities of Fig. 4c in the middle of the microchannel. Note that the plasma treatment temporarily makes the surface of the propulsors hydrophilic so no bubbles are trapped. As such, a 12 h waiting period is thus required for the plasma treatment to wear off and before the microchannel is ready for use. In the case where the pattern is not perfectly positioned in the middle of the channel, we can move the demarcation line between the two liquids by adjusting the flow rate of one pump. The line of demarcation stabilises in a matter of seconds at the flow rates that we use.

Protocol
Experiments are conducted with cycles of 2 second ultrasound followed by two seconds at rest. Two seconds is sufficient to restore the flow to an unmixed state and prevent additive effects between two successive mixing cycles. It is also much greater than the mixing time that we will calculate below in section 5.1. The acoustic pressure is increased over the course of each experiment to investigate its effects on mixing and there are typically 20 to 50 cycles in each experiment.
In order to quantify mixing, we first define an 'analysis box' downstream of the pattern of AMBs (see Fig. 5). The width of the box is equal to the width of the channel while its height has been set to 100 pixels (287 μm). The distance from the last capsule to the beginning of the box is between 700 and 1100 μm. The box is always positioned sufficiently downstream so as not to include any of the vortices generated by the AMBs. We then measure the RMI mixing

Micropropulsor array optimization
AMB arrays are active mixers but their physical presence in the channel could also contribute to passive mixing. Separate experiments show that when the base flow passes near nonexcited capsules it is deflected on a length scale close to the diameter of the capsules before resuming its linear trajectory (see Fig. 4c). Passive mixing is thus negligible so we will attribute mixing effectiveness in what follows to the flows generated by the active capsules. 4.4.1 Line mixers: effect of orientation. In order to make the microchannel assembly process easier, we first designed an active surface consisting of a straight line of propulsors, each 50 μm apart (centre to centre), stretching across the channel and pointing against the input flow. This is easy to assemble since a line of propulsors longer than the channel width can be printed, and thus the centre of the array does not need to be aligned with the centre of the channel when attaching the microchannel and substrate, unlike with other patterns. The line has 19 AMBs within the channel, each containing an 18 μm diameter bubble.
The flow generated is shown in Fig. 6a and takes the form of a large vortex pair, taking up the whole width of the channel, with each vortex ∼500 μm wide. This is a large collective flow generated by individual small elements. Mixing-wise, it only reaches RMI = 61% (see Fig. 6b) and the efficiency decreases for higher acoustic pressures. The effect is clearer when the peak mixing efficiency measured during each pressure cycle is plotted against the acoustic pressure (see    . 6c). This decrease in efficiency is due to the capsules gradually filling with water. Our observations suggest therefore that when submitted to the permanent flow in the microchannel and forced with ultrasound, bubbles tend to disappear much quicker than when placed in a liquid bulk. As a comparison, the propulsors used in the triangle experiments of subsection 3.2 can be used for several hours. Note that in the triangle experiments, silanisation even further increases the bubbles' lifespan as it makes the capsules hydrophobic 39 but silanisation is not an option in microchannels where a silanised surface leads to poor wetting and the formation of very large bubbles. Nevertheless, this experiment does confirm that closely-placed AMBs in a microchannel can generate large flows.
Next we tilt the line shape by 45°and point half of the AMBs down/left and half up/right (see Fig. 7a). The aim was to both lengthen the lifespan of the bubbles by having them not pointing directly into the permanent flow, and to increase the mixing efficiency. Both goals were achieved; indeed, for 150 seconds almost all propulsors kept air inside themselves and a mixing efficiency of 90% (Fig. 7b) was reached. However this shape requires a large number of capsules in order to fill the channel and to achieve good mixing.
This then motivated a study of a left/right line consisting of only 6 propulsors (see Fig. 8a), three times fewer than the tilted left/right mixer. All propulsors there stay active for the whole 250 seconds of the experiment and reach a maximum mixing efficiency of 98% (Fig. 8c). Therefore we are able to fully mix a wide channel (e = 1 mm) with only 6 propulsors, arranged in a shape that takes as little as 90 minutes to print. We propose that the high efficiency of this pattern is due to the two vortices communicating with each other as seen in Fig. 8a, exchanging fluid from the two sides of the channel.
In order to investigate whether we can predict this behaviour numerically we used the Stokeslet model described in subsection 2.2. The numerical flow field shown in Fig. 8b contains similar vortex patterns to the experiments for a base flow velocity of 1.2 mm s −1 and quantitatively similar parameters as the experiments. More importantly, it predicts the exchange of liquid from one side of the channel to the other.
4.4.2 V shape mixers. Motivated by trying to separate the flow in order to maximise mixing, we developed a V-shapedmixer (see Fig. 9), aligning the point in the V shape with the line of demarcation between the two liquids. We examined the effect of changing the distance between propulsors on the mixing efficiency. Gaps of 100 μm and above between propulsors gave insignificant mixing (RMI < 20%). Smaller gaps of 50 μm performed far better but, similarly to the line mixers, the lifespan of the bubbles is greatly diminished when propulsors are facing into the flow. In the case of the V mixer, the bubbles can disappear in 7 or 8 mixing cycles, as seen in Fig. 9b. This explains why the maximum mixing is reached (for 50 μm gaps) at P ac = 806 kPa and decreases at higher pressures (Fig. 9c). However, this pattern is the most efficient at very low acoustic pressure, where it reaches RMI = 59% for P ac = 224 kPa.
4.4.3 L shaped mixers: formation of a vortex pair and mixing time. We next designed a L mixer, inspired by the natural phenomenon of fire whirls which forms when trees in the shape of an L burn inducing a wind vortex which fire is channelled into. 40 The L mixer was designed to break the symmetry of the flow and generate a large mixing vortex. We now look in detail at the formation of this mixing vortex.
Firstly, we place the L propulsor arrangement in an open pool and investigates if a mixing vortex is generated. We use a circular pool without a ceiling, made in a PDMS cell filled with PBS. Unlike with the triangle arrangement of section 3.2, this active surface (Fig. 10a) induces a collective flow whose shape agrees with our numerical model (Fig. 10b). A collective flow is induced by the L and not the triangular arrangement because of the smaller distance between propulsors (50 μm for the L vs. 100 μm for the triangle). However in this open pool we do not observe the formation of a large vortex but only small individual vortices generated at both extremities of the L (both experimentally and numerically).
We next confine the L in an 80 μm high microchannel in order to investigate the role played by the top wall. The channel has the shape of a circular pool of diameter 5 mm. The specific microchannel used for Fig. 11a and b had a thick PDMS layer, so attenuates ultrasound from the transducer more than usual since higher acoustic pressures were used for these experiments than the ones shown in Fig. 10. However the streaming velocities generated by the propulsors (v = 10-30 mm s −1 ) were comparable.
In this confined channel, a large leftward vortex is observed, together with a rightward vortex (see Fig. 11a), in agreement with numerical simulations (Fig. 11b). From this result, we deduce that the non-constrained flow in Fig. 10 is three-dimensional, flowing upward, and that the addition of a ceiling constraining that upward flow is compensated for by the formation of a pair of vortices. This hypothesis is confirmed by simulations. A vertical cut alongside the length of the L shows that with no ceiling, the liquid flows upward. The ceiling blocks that flow (see Fig. 11c) and conservation of mass must result in the formation of two vortices.
We next confine the L further, placing it in an e = 1 mm wide and 80 μm high microchannel. We will analyse the formation of its large vortex as P ac is increased, see Fig. 11d-f. At lower pressures, we see in Fig. 11e individual vortices, including the backward flows characteristic of individual propulsors on 10 μm poles (see Fig. 2a). As the pressure increases, the backward flow becomes small compared to the strong forward flow, and a large leftward vortex takes shape (see Fig. 11f). There is a threshold pressure at which the vortices form (specifically, P ac = 216 kPa).  Considering the flow velocities in these vortices in more detail, we again use the L in a 80 μm high circular pool, the setup which produced Fig. 11a and b. Particles alongside the vortices follow near-ellipsoidal trajectories, but their average velocities change as they move around the vortex. The vorticity analysis in Fig. 12a (using the freeware FAST 32 ) confirms that it is a rotational vortex. The mean velocities of the particles is shown in Fig. 12b. The particles reach a maximum of 27 mm s −1 inside the L and decrease as they follow the elliptical trajectory further away from the L. One of the leftward trajectories is analysed in detail in Fig. 12c. Starting inside the L the velocity is at a maximum. The velocity along the vortex v t then decreases until the attracting force from the AMBs affects the particles again and the velocity increases. There is thus no simple relationship between the particle velocity and the radial distance from the centre of the vortices. From this analysis, we can also calculate the time τ m taken for particles to complete a full loop of the vortex. This is shown on the insert of Fig. 12c, which reveals that it takes between 5 and 400 ms to complete a loop of the vortex, depending on the distance from the vortex centre.
We finally test the mixing efficiency of the L pattern inside our Y-shaped microchannel. Our measurements show that the L shape gives highly efficient and sustained mixing as shown in Fig. 13. The generated vortices produce nearly perfect mixing in our microchannel for P ac ≥ 1 MPa (see Fig. 13b), which is a significantly lower pressure threshold value for nearly perfect mixing than the other array shapes. This arrangement is the most efficient of all the patterns we considered. Moreover, the propulsors suffer significantly less from loss of air than any other patterns since they do not face into the permanent flow. This pattern is also robust, since experiments showed the level of mixing is unaffected whether the L is placed in the centre of the channel or slightly closer to either side wall, or even if the L is tilted. This makes the L array an excellent candidate for microfluidic applications.

Mixing process and mixing time
Microfluidic mixers aim, essentially, to accelerate molecular diffusion. One method, used by grooved microchannels, is to generate chaotic mixing by folding the flow onto itself. Another method, exploited by our AMB arrays, becomes apparent when a blob of dye is placed inside a thin channel and a laminar flow is imposed. The dye is stretched, taking on the parabolic shape of the flow. At the front of the parabola, the dye diffuses towards the walls, while at the back of the parabola it diffuses inward towards the channel centre. This classical process is called Taylor-Aris dispersion. 15,16 Stretching the liquid with high shear rates can increase mixing efficiency as used by circular mixers. 41 The Péclet number for our experiments is given by where v = 1-100 mm s −1 is the flow velocity generated by the propulsors, h = 80 μm, and D ≈ 10 −5 cm 2 s −1 so that Pe ≈ 10 2 -10 4 . Due to the very strong flows generated by our propulsors, the Péclet number of our experiments is 1 to 2 orders of magnitude larger than is common in microfluidic experiments. This confirms that the convection-diffusion at play in Taylor-Aris dispersion will be significant. Assuming a steady incompressible pressure-driven laminar flow in the channel, the dispersion is governed by the equation: (4) where c is the concentration and v the velocity field. 41 For long time scales and a narrow channel cross-section (h ≪ e), the effective diffusion coefficient D eff is D eff = D(1 + κPe 2 ), where κ is a shape coefficient. With Pe ≫ 1, the relative significance of Taylor dispersion to molecular diffusion is κPe. 2,42 If we assume v ≈ 20 mm s −1 such as in Fig. 10a and   11a, then . This requires an approximate value for κ which, in general, is difficult to compute exactly. But we can estimate κ as being between its lower limit of two infinite parallel plates (κ = 0.0048) and its higher limit of a trapezoidal channel (κ = 0.0329). 42 Thus, in our case, Taylor dispersion contributes at least 99.99% to the total dispersion. In addition, we note that the Taylor-Aris diffusion time scale τ Taylor is . 41 Mixing time is typically defined as the time taken for a particle to travel from an unmixed to a perfectly mixed region. 4 However, for our mixers, the border between the unmixed and mixed regions is dynamic, making it difficult to exploit this definition to measure the mixing time. This is corroborated when analysing the grey profile alongside the segment defined by the dashed white line on Fig. 13a as using this definition leads to an unrealistic mixing time of 10 ms.
For the L-shaped mixer we instead hypothesise that the mixing is of the same order of magnitude as τ m , the time taken for a fluid particle or tracer to complete a loop inside the vortex. For P ac = 1.53 MPa, a bead takes τ m = 10-400 ms to perform a complete loop of a vortex (Fig. 12c), depending on its distance from the centre of the vortex. Therefore, a good approximation for the mixing time in our e = 1 mm microchannel, is to average τ m over all tracers within 500 μm of the vortex centre (assuming the vertical branch of the L is in the middle of the channel). This gives a value of τ average = 111 ms, which is much lower than the 2 second cycle period we use in our experiments. Unsurprisingly, τ average is also on the same order of magnitude as the Taylor diffusion time scale τ Taylor . The mixing time will be on the same order of magnitude as τ average , since we expect mixing will take only a few loops. On the insert in Fig. 13b the RMI measurements are plotted against time, focusing on the start of a new cycle when the driving force is turned on. The acoustic pressure was P ac = 1.48 MPa, similar to the experiment of Fig. 11a. The results show a steep increase in RMI before plateauing. Here the velocity inside the vortex is much greater than the input flow and the mixed fluid reaches the analysis box in only 30 ms as compared to at low acoustic pressure where it takes about 0.3 s for the mixed liquid to be advected to the analysis box downstream (see Fig. 5). To go from RMI = 20% to RMI = 90%, it takes τ exp = 380 ms so just under 4 loops.

Theoretical outlook
Using the experimental mixers as templates we finally considered variations on their design and addressed theoretically their impact on mixing using the method described in section 2.3. Specifically, we considered rotations of both the mixers and of their rows of AMBs (or the AMBs themselves). We also varied the starting position of the fluid boundary across the width of the channel, considering the line extension at every 10 μm position across a 1 mm width. Our theoretical results indicate that the position of the boundary between the two fluids to be mixed relative to the AMBs does affect the mixing efficiency. If there is a circulation in the flow, it will benefit mixing if it is able to pull in both fluids rather than just rotating one of them. The position of the boundary affects which local features stretch and contract and by how much (see Fig. 14), although some mixers (like the long left/right mixer) do have wide bands of similar boundary stretching. The two biggest differences in behaviour is that either the boundary is stretched but passes out the other side of the mixer or it is pulled and caught in local vortices.
Our analysis best compares the cases where the boundary passes completely through the mixer, which is useful for a mixer which is permanently 'on' as vortices trap fluid, delaying fluid from progressing down the channel even once it is mixed. In this scenario the best of the experimental mixers was a long left/right followed closely by the tilted left/ right mixers and then the L shape (see Fig. 4). In addition the left/right mixer has a wide band of over 300 μm where a similar high level of stretching is observed, very useful when diffusion starts to smudges the boundary. However, theoretically better mixing can be achieved with mixers such as those shown in Fig. 15. For the L-mixer, generally good mixing is observed theoretically when the two rows of AMBs point along their own length.
We expect vortices to create the best mixing, particularly when the mixer is turned 'on' and 'off' as in the experiments. Local vortices cause the greatest stretching, but the boundary then became tightly wrapped so the diffusion length scale becomes very significant. If the AMB strength was increased, corresponding to higher driving pressures, vortices would become larger in size thus potentially reducing this issue. However the unconstrained large vortices throw flow out sideways rather than capturing it thus the side walls are significant. Future work could include adding side  walls to our model and consider stretching of the boundary already within the mixer.
Both mixing arrangements of AMBs presented in Fig. 15 were studied experimentally. The modified L shape suffered from an unexpected issue of short bubble lifespan in the horizontal line of AMBs pointing across the channel; it was only able to reach RMI = 37% for P ac = 448 kPa before the bubbles disappeared (not shown). This new L shape was already rotated 180°from Fig. 4(a) to Fig. 15 so the AMBs along the vertical line of the L would not point against the flow, but we did not anticipate the AMBs on the horizontal line would experience issues as well. Clearly these lifespan issues would not be identified by the numerical model. The variation of   the left/right mixer was far more promising experimentally, with RMI = 65.5% obtained at a low acoustic pressure of P ac = 537 kPa (see Fig. 16b). However at higher pressures the AMBs pointing to the right lost their bubbles. Future experimental work should tackle the issue of bubble lifespan, for example by making the printed capsules hydrophobic.

Conclusion
In this paper we have demonstrated that, due to the versatility of our fabrication setup, bubble-based acoustic micropropulsors can be arranged in any pattern and so generate a rich variety of flows. Those flows are sustained, powerful, and well predicted by our hydrodynamic model. Moreover, a multitude of objects can be activated by the same acoustic wave. Furthermore, these arrays have shown great promise for use as mixers in microchannels. A pattern of individual micropropulsors can generate large vortices close to the millimetre scale. Those vortices, present for shallow microchannels, mix through Taylor-Aris dispersion in very short time-scales. They can be predicted by modelling the propulsors as individual flow singularities bound by two parallel walls. Our experiments demonstrate the possibility of using relatively few propulsors for mixing a 1 mm wide microchannel: an L pattern containing 12 propulsors is the most efficient from a relatively low acoustic pressure, with long lifespan bubbles, and a left/right pattern with only 6 propulsors is able to fully mix at slightly higher pressures.