We used the wild-type algal strain Chlamydomonas reinhardtii (CC125) for our studies. Single colonies of these algae were picked from slant cultures and inoculated into Tris-minimal growth media (Supplemental Material Section S1). These media are based on the standard TAP medium [37], but omit acetic acid and HCl is used to titrate to pH  7 , and were chosen to ensure purely photosynthetic growth of cells so that they would be synchronised to light-dark cycles, following [38]. Indeed, liquid cultures of the microalgae were then grown in a 14:10 h light-dark cycle on a rotary shaker at 100 rpm and continuously bubbled with air, as in [38]. The shaking incubator (Infors Minitron) was maintained at a temperature of 25°C, and provided photosynthetically active radiation (PAR) at 315 − 325  μmol / m 2 s , as measured with a PAR meter (Skye SKP200). It took around 7 − 10  days for a culture to reach a concentration of 1 − 2  million cells / mL . Thereafter, it was sub-cultured by mixing 10  mL of grown algae into 140  mL of fresh Tris-minimal media until the cell count, measured with a Z2 Coulter counter (Beckman Coulter, Brea, CA), reached 1.5   million / mL ; this took about a week. Subsequently, algae were diluted everyday by replacing 50  ml of the culture with fresh media. This protocol maintains the algal count between 1.2 − 1.5   million / mL with mean diameter of 4.5 − 5  μm (estimated by using the Coulter counter); the subcultures can be used for 10 − 15  days . Care was taken to do experiments with microalgae harvested during the light phase of the growing cycle to avoid variations in the swimming parameters, and in particular the swimming speed, which have been observed at the onset the dark phase [38]. All the experiments were carried out in square plastic cuvettes of external dimensions 12.5 × 12.5 × 45   mm 3 (Sigma-Aldrich, filling volume 2.5 ml) filled with 2   ml of algal suspension. The imaging was performed using a monochrome CMOS camera (Pointgrey, Grasshopper3 GS3-U3-23S6M) fitted with a macro lens (Sigma 17–70 mm f2.8-4). The cuvette was illuminated from the side by a red ( 660   nm ) square 100 × 100 mm LED array (Advanced Illumination BL1960, Rochester, VT, United States), as shown in Figure 1A. This illumination was used as it allowed to image the suspension laterally without triggering a phototactic response [3]. The concentration of microalgae in the cuvette was estimated from the transmitted light intensity across the short dimension of the cuvette by applying the Lambert-Beer law: the intensity recorded by the camera (measured in arbitrary units, a. u.) can be converted into algal concentration ( million / mL ) from the calibration curve shown in Supplementary Figure S1. The intensity decays as I = I o exp ( − A . C ) , where I o = 179 is the intensity in arbitrary units in the presence of cuvette containing just Tris-min medium, A = 0.22 is the attenuation coefficient and C is the algal concentration in million / mL . This exponential decay provides a mapping to concentration, with an excellent fit for intensity data higher than 50   a . u . , and an R-squared value of 0.99 for a fit across the range of values (see Supplementary Figure S1). In the experiments described below, the swimmer concentration was then quantified from images by first measuring integrated pixel intensity in selected regions (see e.g., Supplementary Figure S2) of the cuvette using ImageJ, and then mapping to actual concentration values using the calibration curve just described.

For the phototaxis experiments, a blue LED (Thorlabs M470L2, nominal wavelength 470 nm) is mounted above the cuvette at a distance of 47   mm from its base. Using a PAR meter, the light intensity at the base of cuvette containing only the media was 16 − 18   μmol / m 2 s , whereas the intensity immediately below the LED is 150 − 160   μmol / m 2 s . The cuvette is separated into an upper “harvest” and a lower “reservoir” region by a porous layer of glass beads. The latter was achieved by folding a rectangular wire mesh so that it attaches to a cuvette, and overlaying it with glass beads of diameter 425 − 600  μm and various weights, as shown in Figure 1B. Glass beads of diameter 425 − 600   μ m , provide a porous medium with spacings large (∼ 50 − 200   μ m ) compared to the size of individual algae ( < 10   μ m ), so that the latter could easily move through the pores. To initialize experiments, first an empty cuvette was filled with an algal suspension approximately up to the mesh height and thereafter the mesh was installed. To make a porous layer of various thicknesses, beads of appropriate weight were placed over the mesh. Finally, more algal suspension was poured from the top to create a harvest region of height ≈ 0.5  cm . The experiments reported below also considered the case of a bare mesh with no beads.

We present here the details of a simplified model of the concentration of swimming microalgae into the upper “harvest” region by light. The model describes the case of a suspension of microalgae with a porous layer near the top (mesh + beads), as shown in Figure 1, and we shall also apply it below to consider the case of a bare mesh. For the mesh + beads case, the suspension of microalgae is divided into three regions, an upper harvest region (u), a porous layer region (p) and a lower “reservoir” region (l). Photogyrotactic migration delivers microalgae to the upper region from the lower region through the porous region. As evident from our results and discussion (see Section 3.4, 4 below), the dynamics underpinning the concentration are complex; the challenge of describing them with a full photogyrotactic model is beyond the scope of this paper. Instead, we seek here to formulate a model to capture the essential features of the concentration process into the upper harvest region. We make the reasonable simplifying assumption that: 1) the average concentrations in the upper, porous and lower regions evolve slowly compared to the observed photogyrotactic dynamics; we consider here spatial and temporal variations separately, and assume a steady state for the fast dynamics in the upper region. We further assume that 2) phototactic migration is the dominant process and brings cells to from the lower region to the porous region, with swimmers migrating straight upwards toward the light (there is no dependence on light gradients, only light direction) at the maximum phototactic speed, equal to the mean swimming speed of the population, V _s . In the lower region, we assume that 3) the mean concentration is representative of the concentration of cells swimming into the porous layer. In the porous layer region, we assume that: 4) the speed of the swimmers is slowed down by collisions with the porous medium, but the swimming direction continues on average to be upwardly directed by phototaxis. In the upper harvest region, as well as the average concentration dynamics, we also consider a 0^th order spatial model of phototactic concentration. To set this up, as assumed above, we posit that there is a separation of timescales between the migration of cells from the lower region (slow) and the redistribution of cells in the upper region (fast). We further assume that: 5) the effect of flow is negligible prior to the formation of the plume from the upper surface (see Figure 4B); 6) upward phototactic swimming at the maximum speed V_s and diffusion dominate the fast suspension dynamics (gyrotactic effects are negligible); 7) diffusion is assumed approximately isotropic; 8) the meniscus at the top of the suspension is flat (any effects of curvature are neglected).

With the assumptions above, denoting by   c ¯ i   the average concentrations in regions i = u ,   p ,   l for the upper, porous and lower regions, respectively, and with A the cross-sectional area of the cuvette, the average cell numbers N ¯ i   in the three regions evolve according to the following balance equations: d N ¯ l d t = −   A   V s   c ¯ l , (1a) d N ¯ p d t =   A   V s   c ¯ l −   A   V e f f   c ¯ p , (1b) d N ¯ u d t =   A   V e f f   c ¯ p . (1c)

Equation 1a describes the loss of cells from the lower region due to the phototactic flux of cells, of concentration   c ¯ l , swimming into the porous region at speed V s . Correspondingly, the porous region, as described by Equation 1b, gains an equal and opposite flux. This region also has a loss term due to cells, of concentration   c ¯ p , swimming at a speed V e f f into the upper region. The upper harvest region, as shown in Equation 1c, has an equal and opposite gain. The speed V e f f is the effective swimming speed of the microalgae within the porous layer, which is given by [39]. V e f f = V s τ c τ + l o b s τ , (2) where V s is the “free” mean swimming speed of the microalgae, τ c = λ V s is the time between collisions with the beads in the porous layer, and λ is the swimmer mean free path. The timescale τ = τ c + τ R is the total porous travel time, including the residence time τ R that a swimmer spends at an obstacle. These parameters were recently measured experimentally for C. reinhardtii (see Table 1). To express system (1) in terms of concentrations only, we note that the mean number of cells in regions i = u ,   p ,   l can be written as N ¯ i = Ah i   c ¯ i , where, as above, A is the cross-sectional area of the cuvette, and h i ,   c ¯ i   are the height and mean concentration in region i , respectively. Substituting into (1), we thus obtain, dividing both sides by the respective h i , d c ¯ l d t = −   α   c ¯ l , (3a) d c ¯ p d t =   β   c ¯ l − γ   c ¯ p , (3b) d c ¯ u d t =   δ   c ¯ p (3c) where we have defined the upswimming rate constants α = V s / h l , β = V s / h p , γ = V e f f / h p and δ = V e f f / h u . Equation 3a has immediate solution   c ¯ l = k 0 e   −   α  t , where k 0 is a constant. The system of Eq. 3 can then be solved analytically by substituting this solution into (3b), and the resulting solution (e.g. by using the integrating factor e β  t ) into (3c). Applying the initial conditions   c ¯ l ( 0 ) = c ¯ l 0 ,   c ¯ p ( 0 ) = c ¯ p 0   and     c ¯ u ( 0 ) = c ¯ u 0 , where c ¯ i 0 represent the initial average concentrations in the three regions, we find:   c ¯ l ( t ) = c ¯ i 0 e −   α  t , (4a)   c ¯ p ( t ) = c ¯ l 0 β γ −   α e −   α  t + ( c ¯ p 0 − c ¯ l 0 β γ −   α )   e −   γ  t , (4b)   c ¯ u ( t ) = c ¯ u ∞ − c ¯ l 0 βδ α ( γ − α ) e −   α  t − ( c ¯ p 0 − c ¯ l 0 β γ −   α ) δ γ e −   γ  t , (4c) where we have defined the long-time concentration in the upper region as c ¯ u ∞ = c ¯ u 0 + c ¯ p 0 δ γ + c ¯ l 0 βδ α γ = c ¯ u 0 + c ¯ p 0 h p h u + c ¯ l 0 h l h u   , (5) and where, recalling the definitions of the constants  α , β , γ   and δ , we have re-written c ¯ u α in terms of the heights of the regions. Thus, it is clear from Equation 5 that, in this simple model, the long-time (maximum) concentration in the upper region occurs when all swimmers from the porous and lower regions have concentrated themselves into the upper region.

We also consider the “mesh-only” case (without a porous layer of beads). The derivation, shown in Supplementary Appendix SB, is similar and provides the temporal evolution of the mean concentrations as c ¯ l m ( t ) =   c ¯ l m 0 e −   α  t , (6a)   c ¯ μ m ( t ) =   c ¯ l m ∞ −   c ¯ l m 0 η α e −   α  t , (6b) where the superscript ‘m’ denotes concentrations in the mesh-only case, and we have defined the rate constants α = V s / h l , which is as in the porous layer model (but takes a slightly different value because of the different value of h l , see Supplementary Table S2), and η = V s / h u . For the mesh-only case the concentration in the upper region at long times is given by c ¯ u m ∞ =   c ¯ u m 0 + c ¯ l m 0 η α =   c ¯ u m 0 +   c ¯ l m 0 h l h u   . (7)

This corresponds to the concentration in the upper region occurring when all microalgae have swum into it from the lower region.

In the upper region we observe that cells accumulate strongly at the surface. To describe this, we can use a simplification of the Williams and Bees model [26]. By virtue of assumptions 5)–8) above, as shown in Supplementary Appendix SA, the full swimmer conservation equation in the Williams and Bees model simplifies to: ∂ c u ∂ t =   − ∇ ⋅ [   V s   c u  k   −   D   ∇ c u   ] , (8) where k is a unit vector pointing upwards and D is the diffusivity, approximated as isotropic, by assumption 7) (see Supplementary Appendix SA for more details). By assumption 1), we have a steady state, so that (8) implies V s     c u  k   −     D   ∇ c u =   c o n s t . , (9) where k is a unit vector pointing upwards. Imposing a no flux condition at the upper boundary (flat for simplicity, assumption viii) requires (   V s     c u   k   −     D  ∇ c u ) ⋅ k = 0 on z = h , so that Equation 9 becomes d c u d z =   V s D c u   , (10) which integrates to c u     =   k 1   e z l p   , (11) where we have defined a characteristic phototactic accumulation lengthscale l p = D   V s , and where k 1 is a constant. To find the latter, we use the fact that the average background concentration is given by c ¯ u , that is, taking z = 0 at the bottom of the upper region and z = h u at its top, c ¯ u = 1 h u ∫ 0 h u c u d z . Thus, integrating Equation 11 gives k 1   = c ¯ u     ( e   h u / l p − 1 ) − 1   h u / l p , so that finally the distribution in the upper region is given by c u ( z ,  t )   =   c ¯ u ( t )     h u l p     e z l p e   h u l p − 1       , (12) where that the mean concentration as a function of time, c ¯ u ( t ) , is provided by Equation 4c.

Prior to considering the effect of light on a suspension of C. reinhardtii placed in the cuvette, we will consider the initial condition of the suspension in the lower reservoir region, which will be the same starting point for all subsequent experiments. With the blue LED light off, microalgae were mixed into the cuvette and the suspension was allowed to stabilize in the presence of only red illumination from the side (see Figure 1A), which does not elicit a phototactic response (Supplementary Video S1) [3, 40]. The suspension images and profiles are shown as a time-series in Figure 2: over a few minutes, the suspension (initial concentration ∼   1   million / mL ) settles into a distribution where the majority of cells reside at the bottom of the cuvette; a steady distribution is observable beyond 4 min. The gradient in concentration already visible for the concentration profile at t = 0 is due to a lag in transferring the cuvette to the imaging setup after mixing: some settling has already occurred at the first instance of imaging. The steady distribution observed beyond 4 min would, for other species of gyrotactic swimming algae such as Chlamydomonas augustae [2], also display features known as “bottom-standing plumes” (see also Figure 1A in [7]). For C. reinhardtii, however, the bottom-standing plume structures are not discernible in the images of the cuvette (though plumes can be seen to descend from its sides, see Supplementary Video S1), though they can be seen in a larger container (see Supplementary Figure S3). It is possible that, for this species, the cuvette width is too narrow to give rise to central bottom-standing plume structures.

We consider here the effect of light on a suspension of microalgae in a cuvette when the surface of the suspension is free (the metal mesh applied in the next section has been raised above the surface). This experimental scenario can be seen in Supplementary Video S2, stills of which are shown as the sequence in Figure 3. Initially the blue LED illumination is switched off and the suspension is distributed with the majority of cells at the bottom, as described in the previous section. Then the LED is switched on, and the cells in suspension phototactically respond to the light, migrating upwards toward the surface (Figure 3, t = 3  min ). Concomitantly, instabilities arise throughout the suspension, resulting in meandering plumes (Figure 3, t = 3 − 5  min ). These are of photogyrotactic origin, as discussed below. In the span of ∼   6 minutes phototactic migration appears to have delivered many swimmers to the surface, leaving the bulk of the suspension depleted. This surface accumulation is gravitationally unstable because of the negative buoyancy of surface-accumulated cells: it results in the formation of a plume instability seen to originate from the middle of the meniscus of the suspension surface (Figure 3, t = 9  min ). The plume structure wiggles around but once formed, is dynamically stable (Figure 3, t = 18  min ), delivering cells to the bottom of the container. Once they reach this, the microalgae migrate back up to the surface to join the plume, and so forth. When the light is switched off (Figure 3, t = LED Off + 17  sec , +   3 ,   6  min ), the phototactic migration toward the surface stops and the surface accumulation sinks as a broader, non-meandering plume. This takes the cells to the bottom of the cuvette, where they once more settle into a distribution similar to the initial one.

In this section we consider the case of a metal mesh immersed at the top of the microalgal suspension. As described in the methods, the pore size of the mesh is 350   μm , so individual microalgae (∼ 5  μm in diameter) easily swim through it. A typical experiment is shown in Supplementary Video S3, stills of which are presented as a sequence in Figure 4A (top row). As for the free surface case, blue LED illumination is initially off, and the suspension is distributed with most cells at the bottom (Figure 4A, first still). The LED is then switched on and microalgae migrate upward in response to the light (Figure 4A, t = 3  min ). The response is broadly similar to the free surface case, but there are some interesting differences. One such striking difference is that the mesh creates a pattern of light and shadow to which the microalgae visibly respond photogyrotactically, forming accumulations (‘phototactic curtains’), see Figure 4C. The average width of phototactic curtain feature is 630 ± 53   μm , nearly twice the mesh pore size, showing that the curtains are not the result of shadowing by the mesh, but genuine phototactic structures originating from the response of the microalgae to the local light profile. As in the case of a free surface, when the density of cells phototactically accumulated at the surface becomes too high, a plume of dense cells forms and sinks. However, viscous resistance caused by the mesh pores prevents the plume from completely sinking beyond the mesh, and instead a cloud-like plume structure is seen to be trapped, hovering above the mesh (Figure 4A, t = 9 ,   18  min ). Not all the plume-cloud is trapped, negative buoyancy is sufficient to cause some of it to escape through the mesh forming a meandering secondary plume, similar in appearance to the one observed in the free surface case (Figure 4A, t = 9 ,   18  min ). While the light is on, these structures appear dynamically stable. As the light is switched off, however, the curtains and cloud structure disperse, cells sink through the mesh, and the escaped plume sinks down straight, again similarly to free surface case (Figure 4A, t = LED Off + 19  sec , +   3 ,   6  min ). This emphasizes the stabilizing influence of phototaxis: none of the observed structures could be possible in the absence of the light. Both mesh and light are critical for supporting the plume-cloud.

We next turn to the case where a porous layer is placed on top of the suspension. As described in the Methods, the porous layer consists of glass beads overlaid onto a metal mesh (the same as was used in the previous section). The beads are around 425 − 600   μm in diameter, which results in interparticle spacings ∼   50 − 200  μm (from microscopic observation). Thus, individual algae ∼ 5   μm in diameter can swim through the porous layer. We studied the effect of light on suspensions of microalgae overlaid with porous layers, quantified by the weight of the beads placed on the mesh. A typical experiment with a layer weighing 400 mg is shown in Supplementary Video S4, and stills from this video are presented in Figure 4B. As in previous cases, the LED light is initially off and the suspension is distributed with the majority of cells at the bottom (Figure 4B, first still). When the LED is switched on, the initial suspension dynamics are similar to the mesh-only case (Figure 4B, t = 3 − 6  min ), displaying instabilities as the microalgae respond to the light (but with no curtains visible). However, for this case, we were also able to observe clusters of cells swimming upwards as waves in response to the light, see Figure 5A for an example. Averaging over five such waves, we found them to have a mean speed of 190 ± 60   μm / s . This is faster than mean swimming speed of individual algal cells, 80  μm / s [38], possibly as a result of advection by upwelling fluid in the lower region of the cuvette generated by the photogyrotactic suspension dynamics. The large deviation in the speed of the waves could also be due to the interaction of the waves with other photogyrotactic structures and up/downwelling flows in the suspension. Once cells have had time to accumulate in the harvest region and on the surface of the suspension, a plume-cloud structure originating at the low point of the meniscus forms above the porous layer (Figure 4B, t = 9 ,   18  min ). The plume-cloud appears more diffuse than in the mesh case. The time taken for the plume-cloud to arise in 10 out of 12 experiments used for the analysis is between 6 − 10  min from when the LED light is switched on, as shown in Supplementary Figure S4. Unlike the case of the mesh, the plume does not leak through the porous layer into the suspension: in the presence of light, the viscous resistance offered by the porous layer is sufficient to stabilize the plume-cloud. Instead of sinking the plume-cloud is observed to gradually expand into the upper region. Figure 5B charts this expansion. The lateral extent of the plume structure increases the most between 7 and 9 min after the LED has been switched on, when the plume begins to drop and propagate along the porous layer. After that the plume-cloud achieves a steady structure, probably as a result of balance between influx of cells from the surface, where the plume originated at the low point of the meniscus, and loss to the edges of the harvest region (and resorption to the suspension surface by upswimming).

In view of quantifying microalgal concentration in the upper harvest region above the mesh or the mesh + beads porous layer, it is instructive to chart the evolution of the average concentration of the suspension in this region (Supplementary Figure S2). To identify a porous layer thickness that would not leak into the suspension below, we considered layers of several weights in trial experiments presented in Supplementary Figure S5. We found a general qualitative trend that was similar for all cases: the concentration grows as the light is switched on, and then saturates to a constant value. We focus here quantitatively on the mesh-only case and the “minimally-leaky” mesh + beads ( 400   mg ) case, shown for three repeats in Figures 6A,B, respectively. The averaged profiles are shown in Figure 6C. This makes it clear that the concentration in the upper harvest region of both the mesh and mesh + beads cases, following a dip in concentration due to phototactic accumulation of cells to the upper surface, grows after the LED is switched on and then tends to saturate. The mesh case, however, saturates earlier, probably because of losses to the lower region, such as the plume visible in Figure 4A (9 min). Another interesting quantitative difference between the two cases is the initial rate of concentration, which appears slightly larger for the mesh case. This indicates the concentration process is initially slower when a porous layer is present, than in its absence. As discussed below using the essential model, this makes sense in terms of the microalgae having to make their way through the porous layer, which reduces the swimming speed that sets the concentration rate. The difference in swimming speed will also affect the average time it takes to form the plume, which was measured to be 7.3 ± 0.6   min for the mesh case, while it is 8.7 ± 1.5   min for mesh + beads (Supplementary Figure S4). When the LED is switched off (Figures 6A–C inset), the concentration in the upper harvest region is seen to rise briefly before steadily falling. This is because, with the light off, the concentrated algal suspension in the harvest region no longer responds phototactically and cells accumulated to the surface are released, sinking down as dense fluid. The increase in concentration due to the cells coming off the surface shows that our measurements likely underestimate the concentration in the upper harvest region because of cells “hidden” at the surface. This could account, at least in part, for discrepancies with model predictions discussed below. For our setup, the time after switching the LED off is optimal for harvesting the suspension, yielding a harvest concentration ≈ 5  million / mL (gain ≈ 4.2 compared to the initial concentration) for the mesh + beads case, as compared to ≈ 4   million / mL cells for mesh-only (gain ≈ 3.6). This highlights the advantage of concentrating using a porous layer. The latter also slows down the rate at which the cells sink back through to the lower region, which depends on the layer thickness.

Also shown in Figure 6 are profiles charting the temporal evolution of concentration in the lower reservoir region. As for the upper region, we have measured triplicate repeat profiles for the mesh (Figure 6D) and mesh + beads (Figure 6E), and also evaluated averaged profiles (Figure 6F). We see that, after the LED is switched on, the concentration for the mesh and mesh + beads falls, as phototactic swimming into the upper regions depletes the lower region of cells. However, the depletion appears to saturate, and to a higher concentration in the case of mesh-only, reflecting the greater leakiness of the mesh, as discussed below. Insets in Figures 6D–F display how, with the LED off, the concentration in the lower region rises due the influx of cells sinking from the upper regions.

We have developed a simple model to capture the essential features of the phototactic concentration dynamics, and evaluate it here using parameters for C. reinhartii concentrated using a mesh + beads setup, as shown in Table 1. In Figure 7A, the model prediction using Equation 4c for the average concentration of cells c ¯ u ( t ) in the upper harvest region is shown as a function of time (the concentration process starts at time t = 0, “LED on”). Qualitatively, the predicted behaviour is as in the experimental curves (Figure 6C), with the concentration initially rising and then saturating. However, quantitatively, the concentration values predicted by the essential model are much larger than those seen experimentally. Indeed, using Equation 5 and the parameters in Table 1, the essential model predicts saturation to a long-time concentration c ¯ u ∞ =   8.2   × 10 6   cells/ml. This is of the same order of magnitude as, but approximately double what we observe experimentally ( ≈ 4 × 10 6   cells/ml). Part of the discrepancy is because, as mentioned above, experimental concentration curves underestimate the concentration in the upper reservoir because of swimmers phototactically accumulated and ‘hidden’ at the surface. Another possible reason is that the essential model unrealistically ignores mechanisms causing losses: as illustrated by Equation 5, c ¯ u ∞ corresponds to the concentration obtained when all the swimmers from the porous and lower regions swim to the upper region and do not leave it thereafter. In reality, swimmer diffusion will cause cells to be transferred from the upper to the porous region, particularly at longer times when concentration gradients between the regions are large. Another possibility not accounted for by the essential model is that, if swimmers respond to gradients of light (as opposed to just its direction, as assumed in the model), the denser suspension of swimmers in the upper region shades the region below, changing the light gradient and reducing the phototactic speed of swimmers below, and thus the rate of accumulation. Figure 7A also shows the model prediction for the concentration in the upper region for the mesh-only case. As in experiment, this is seen to initially rise steeper and saturate at a lower value than the case of mesh + beads; numerically, however, the predicted concentrations ( c ¯ u ∞ =   7.5   × 10 6   cells/ml) are approximately double what we measured experimentally. This is for the same reasons as for the mesh + beads case, and additionally in the mesh-only case there are also losses due to the plume leaking through the mesh, as we have shown (Figure 4A, 9 min). The faster rise in concentration observed for the mesh-only case compared to mesh + beads, is due to the difference in upswimming rates in the two cases: for mesh + beads, swimmers are slowed down when they swim through the porous layer. The essential model also allows the prediction of the concentrations in the porous and lower regions, respectively c ¯ p ( t ) and c ¯ l ( t ) , which are shown in Figure 7A. The lower region concentration is seen to decay exponentially to zero, as swimmers evacuate the lower region by phototactic upswimming. The lower region decay predictions overlap for the mesh and mesh + beads cases, so they are not separately visible in the figure (the prediction equation is the same for these cases, and parameters are practically identical). We can compare these predictions with the experimentally determined concentration profiles. As observed in the previous section, these also decay with time, but not to zero: they saturate to a fixed value (Figure 6F), with the mesh-only case reaching a lower value than mesh + beads due to the greater leakiness of the mesh. The essential model fails to predict this saturation and the important difference between the two cases, demonstrating the need to model diffusive transfer and leaking plumes between the regions, and/or a reduction of the phototactic speed. For the porous region, the essential model predicts that the concentration, c ¯ p ( t ) , initially rises, due to influx from the lower region outpacing losses to the upper region, and eventually decays to zero. It was not possible to optically image the microalgae in the porous region and obtain the concentration there, so we cannot make a comparison with the essential model prediction in this case.

Assuming phototaxis and diffusion processes are dominant in the upper region, and that these occur faster than the accumulation from the porous region, we can also use the spatial extension of the essential model to chart the distribution of swimming algae in the upper region, which is provided by Equation 12. We note that, since this model does not fully account for photogyrotaxis, the predictions are only strictly valid prior to the formation of the plume off the upper surface, which we know from experiment occurs ≈   9 min after turning the light on. In Figure 7B the distribution of swimmers is charted at different points in time (prior to plume formation), predicting that the suspension becomes increasingly top-heavy as time progresses. This accumulation, with concentrations reaching ∼   10 8   cells/ml close to the upper boundary, is unstable against its own negative buoyancy, and eventually results in the formation of the plume we observe experimentally. As it is not possible to accurately image the accumulation of cells around the meniscus, we did not experimentally quantify the spatial concentration distribution in the upper region. However, the increasing accumulation of swimmers at the surface is clearly discernible in our image sequences, see Supplementary Video S4. The model predicts that the cells accumulate strongly at the top of the upper region, with no sizeable concentration below a certain height. Instead, our image sequences reveal that there is also a nonzero concentration in the bottom of the upper region (indeed that is what we have measured to obtain Figures 6A–C). This could be accounted for by losses from the accumulation at the surface to the edge of the cuvette, which are not considered in our model.