Hierarchical microphase separation in non-conserved active mixtures.

Non-equilibrium phase separating systems with reactions, such as biomolecular condensates and bacteria colonies, can break time-reversal symmetry (TRS) in two distinct ways. Firstly, the conservative and non-conservative sectors of the dynamics can be governed by incompatible free energies; when both sectors are present, this is the leading-order TRS violation, captured in its simplest form by 'Model AB'. Second, the diffusive dynamics can break TRS in its own right. This happens only at higher order in the gradient expansion (but is the leading behaviour without reactions present) and is captured by 'Active Model B+' (AMB+). Each of the two mechanisms can lead to microphase separation, by quite different routes. Here we introduce Model AB+, for which both mechanisms are simultaneously present, and show that for slow reaction rates the system can undergo a new type of hierarchical microphase separation, which we call 'bubbly microphase separation'. In this state, small droplets of one fluid are continuously created and absorbed into large droplets, whose length-scales are controlled by the competing reactive and diffusive dynamics.

Introduction

Recently, it has emerged that the non-membrane compartments of cells, also known as biomolecular condensates, can be fruitfully viewed as phase-separated liquid–liquid mixtures undergoing chemical reactions [1-7]. In contrast to passive binary fluids, cells continuously consume fuel, driving the system away from equilibrium. Put differently, the liquid mixtures within cells are examples of active matter. The characteristics of phase separation in these biomolecular condensates are distinct from equilibrium counterparts, in that microphase separation is often observed without any energetic cause [4]. (In equilibrium, the latter requires long-range interactions mediated by charged species or block copolymers, for example. [8-10]) This echoes studies of bacteria colonies whose phase separation (coarsening) is arrested by birth-and-death dynamics to give patterns on a finite length scale [11]. This suggests that, like the latter, the physics of biomolecular condensates can be captured schematically by ‘active field theories’ based on minimal, -type ingredients [12-17].

To construct top-down active field theories for particular systems, one can start with the desired symmetry and conservation laws, and write down the lowest-order theory. This follows lines similar to Hohenberg and Halperin’s classification of equilibrium models [18]. In the present work, we focus on theories for a single scalar field , which could be (for bacteria or swimming organisms) a rescaled density, or (for the binary mixtures of interest here) a chemical composition variable. Within the resulting class of minimal scalar models, several mechanisms have been found to arrest phase separation. In our previous work, we explored the lamellar phase and the droplet suspension phase in Model AB [16]. The latter is a canonical representation of phase-separating conservative dynamics coupled to non-conservative reactions, governed by Model B and Model A, respectively, in such a way that the underlying free energy structures in each sector are incompatible. This breaks detailed balance and hence time reversal symmetry (TRS). The finite length scale of the resulting microphase separation (between regions rich in species 1 and 2, respectively) arises from the balance between the reactions (which convert the majority species into the minority) within the two separated phases, and the conserved, diffusive currents that transport species across the interfaces from high to low chemical potential. In steady state, the interfaces cannot be too far apart, otherwise the buildup of minority species within each phase would cause a spinodal instability and create new interfaces. Microphase separation, rather than bulk phase separation, occurs whenever the stable fixed point composition of the homogeneous reaction dynamics lies between the binodals of the conservative phase separation, so that the uniform bulk phases are chemically unstable [11, 16, 19–22].

Meanwhile, phase-separated states have also been reported in fully conservative active systems described by Active Model B+ (AMB+). This model is constructed by adding a complete set of lowest-order non-equilibrium current terms to the standard (passive) Model B [18] (with differing outcomes, in more than one dimension, from an earlier, less complete, Active Model B [12]) [14]. The active terms in AMB+ can send into reverse the Ostwald process for droplets, which would normally lead to full-phase separation by causing small droplets to shrink and disappear while large ones grow [17]. The reverse Ostwald process, which formally

emerges through a negative effective interfacial tension in the expression for the Laplace pressure at a droplet surface, instead equalizes the droplet sizes. (Note that other definitions of interfacial tension do not become negative, so that droplets remain stably spherical for example.) In the absence of noise, this mechanism creates a set of uniformly sized droplets whose number is governed by the initial condition. However, when noise is present, this ultimately selects the length scale for microphase separation by a subtle interplay of various effects. (Those effects include noise-induced nucleation of small droplets and coalescence of large ones, coupled to a collective homogenization of droplet size caused by reverse Ostwald [14].)

As just described, Model AB and AMB+ each involve a single scalar field on its own. More elaborate active field theories can be constructed by coupling a compositional scalar to a momentum-preserving fluid flow, along lines established in passive systems as ‘Model H’ [18]. It is known that this coupling can give further distinctive routes to microphase separation, for example when activity in effect reverses the sign of the interfacial tension that governs the stress exerted on the fluid at the boundary of the droplet (but without doing the same for the Ostwald currents). This happens when the droplet is made of so-called contractile material [13, 15]. In this paper, we omit these further complications and consider only how to combine the types of physics already captured by Model AB and AMB+, in which there is no separate momentum conservation, and the reaction-diffusion dynamics of give a complete description of the pertinent physics.

In our previous work on Model AB [16], we broke TRS only at the lowest possible order in . This entails choosing chemical potentials for the Model B and Model A sectors that could each give an equilibrium system in the absence of the other, so that , but also choosing . This leads to phase separation via the mechanism described previously. However, in a given physical system, once TRS is broken in this fashion it is likely to be broken elsewhere as well. In particular, we know that the active currents represented by the -conserving AMB+ (which do not derive from any local chemical potential, of the form or otherwise), can independently cause phase separation via the reverse Ostwald process [14]. Even if these are formally subdominant when expanding in , such terms can in principle interact in a strongly non-additive way with the microphase separation in Model AB, particularly if parameters are chosen so that the two mechanisms act on well-separated length- and time-scales.

Accordingly, in this paper we construct Model AB+ which includes both types of TRS breaking at once, and use it to give a first account of the interplay between the two mechanisms of microphase separation. Our aim in this brief and exploratory study is not to comprehensively explore the parameter space, but rather to show proof of principle that new physics can indeed emerge from this combination of mechanisms. In particular, in the case where the chemical reactions are slow, we will give numerical evidence for a hierarchical phase separation in which the conservative and non-conservative mechanisms act simultaneously at small and large length-scales, respectively. This results in a new, dynamical steady-state structure that we call ‘bubbly microphase separation’.

The paper is organized as follows. In Sect.2, we summarize the original Model AB as presented in our previous work [16] and recall its various stationary patterns before introducing Model AB+. Section 3 then explores some steady states of Model AB+ with comparisons to the corresponding Model AB pattern without the active current terms from AMB+. These include droplet emulsions alongside the bubbly microphase-separated state. We then present in Sect. 4 a mean-field Ostwald-type calculation for growth of a (liquid) droplet in a bath with (vapour). This explains a transition (or crossover) between statistically different droplet emulsion states upon increasing the target density of the chemical reactions above the dilute binodal of the conservative phase separation. Section 5 briefly summarizes these findings and suggests directions for future work.

Model AB+

We will first briefly introduce both Model AB and Active Model B+, before combining their two dynamics to construct Model AB+. Model AB, constructed as a combination of Model B and Model A, is a minimal field theory for a scalar composition variable , in reaction–diffusion systems in which the conservative part of the dynamics, taken alone, would give phase separation [16]:We assume that . Hence, is the absorbing state composition of the reaction diffusion dynamics (this would be the zero density state for a birth-death system), which marks the physical lower limit for , whereas is the ‘target density’ of the chemical reactions, such that there is a stable fixed point of the dynamics at , if the diffusive tendency to phase separate is switched off. In addition, are mobility constants, are spatiotemporal white noises and are positive constants. The parameter is temperature, or an equivalent noise parameter in cases where the primary noise source is not thermal.

Given the form of chemical potential chosen here, if we had chosen instead , the conservative and nonconservative dynamics would share a free energy so that the model would be an equilibrium one. Much less obviously, but as discussed in our previous paper [16], if is not equal to but differs from it only by terms linear in , the system maps to an equilibrium model with long-range attractions. Hence, by virtue of the quadratic term in (alongside the absence of the cubic term), Eq. (1) admits the lowest-order mismatch in and that cannot be incorporated into an equilibrium model by judicious matching of higher-order terms.

If the reactions are switched off (), the diffusive sector of the dynamics drives the system towards bulk phase separation with the two coexisting phases at the binodal densities, , as long as the mean density ( ) is between the binodals. With reactions, broadly speaking, if the target density of the non-conservative sector is between the binodals, the system exhibits micro-phase separation, where material is created in a dilute region, transported across the interface by the Model B current, and eventually destroyed in a dense region. The morphology depends on the value of , as shown in the upper panels of Fig. (1). For close to zero, the system shows lamellar patterns; otherwise, we either see a droplet phase or its inverse, a bubble phase, depending on which binodal is closer to. In this paragraph, and below, we refer to the two bulk densities as though they were liquid () and vapour () although of course for biomolecular condensates is composition variable describing whether a fluid mixture of species 1 and 2 is composed mainly of species 1 () or species 2 (). In this context, a rightward current of ‘material’ means a combination of a rightward current of species 1 and a compensating negative current of species 2 at constant total mass density.

Model AB breaks time reversal symmetry solely by having mismatched free energies in the nonconserved (Model A) and conserved (Model B) sectors. But of course, each sector can, in principle, break detailed balance on its own. In the Model A sector, this requires gradient terms in that do not stem from a free energy. At leading order in the gradient expansion, the resulting ‘Model A+’ will not lead to microphase separation under steady-state conditions because, without a conservation law, there are no interfaces in the steady state at which such terms could be large. Instead we expect relaxation towards a uniform state at zero , or . These active ‘A+’ terms can of course modify the microphase separation caused by the competing conserved and non-conserved dynamics in Model AB, but we don’t expect them to introduce a second type of microphase separation.

In what follows, we therefore focus on enhancing Model AB by adding active gradient terms to the conservative sector where such terms can lead to microphase separation in their own right. The terms in question are those of Active Model B+ [14],where is the same as before. The and terms are the lowest-order Landau–Ginzburg terms that break TRS in systems with conservation. (One other term at this order can be absorbed into a -dependent square-gradient coefficient without breaking the free energy structure.) It has been shown that the and terms not only shift the binodals, but that the latter can also lead to effectively negative surface tensions for either liquid droplets () or vapour bubbles (), albeit not both at once, depending on the sign of [14]. The tension in question is the one governing Laplace pressures and hence Ostwald ripening, not other physics (such as fluctuations at the interface between phases which remains stable). The resulting reverse Ostwald process causes microphase separation into an emulsion of finite droplets that do not coarsen, which is of course not seen in equilibrium Model B. In renormalization-group terms, this microphase separation appears to be connected with a strong coupling regime at large enough values of . As such, these variables remain important even though, being higher order in the Landau–Ginzburg expansion, they are formally irrelevant in the neighbourhood of the Wilson–Fisher fixed point which controls the critical onset of bulk phase separation [23].

This suggests that equally important physics could also be lost by ignoring these non-integrable gradient terms in the conservative sector for systems with both conserved and non-conserved dynamics, whereas we did ignore them when constructing Model AB as a canonical model for that case [16]. In light of this, we now construct Model AB+ by adding the and terms to the conservative sector of Model AB, so that (1) is replaced by:The rest of this paper addresses Model AB+ as written in the form (3). All numerical simulations are performed in d=2 with periodic boundary conditions. Spatial derivatives are computed in real space following the methods detailed in “Appendix” A of [14]. Time integration is performed using the explicit Euler–Maruyama method [24]. Throughout the paper, unless otherwise stated, we set . Note that we make no attempt at a systematic exploration of parameter space, but focus on selected regimes where new physics can be expected by adding nonzero terms to parameter sets that we used previously to study Model AB. (In particular, the somewhat arbitrary choice for u was adopted in our previous paper [16] and is kept here for consistency. In principle, one can simply set and the results will be qualitatively the same.)

Steady states

Model AB+ has at most one stationary state of uniform density, which lies at the target density . (Any other uniform state has from the reactions, with no diffusive currents that could balance this.) Linear stability analysis finds this uniform state to be unstable to spatial perturbations ifwhere , . The two sides of the inequality represent characteristic relaxation rates via diffusion and reactions, respectively, on length scales set by . In many biological situations, relatively rapid thermal diffusion of chemical species is accompanied by a relatively slow reaction-driven turnover time, so that this ‘slow reaction’ condition is easily obeyed. Hence, in this paper we focus on the unstable regime far away from the threshold of linear instability. We further assume that so the reaction dynamics is approximately linear in for values lying between the binodals. This means that we have an approximate symmetry in our system, similar to Active Model B+ where this symmetry is exact [14].

Some typical steady states are shown in Fig. (1) for , where the Ostwald tension for droplets of the dense phase () are negative, as shown in Tjhung et. al. [14]. These are directly compared with the pure Model AB case (). For low values of , where pure AB is in a uniform state we see an emulsion of these dense (‘liquid’) droplets in the dilute phase (‘vapour’, ), stabilized by reverse Ostwald. At high , chosen so that the Model AB mechanism leads to a microphase-separated state with liquid phase in the majority (hence an emulsion of vapour bubbles in liquid), the active current terms cause phase inversion of this state so that the system is liquid-in-vapour. This is inevitable for sufficiently negative where the conservative dynamics causes microphase separation on a relatively short length scale since, importantly, these terms only stabilize liquid in vapour droplets and not vice versa (unless their signs are reversed, in which case so is the whole phase diagram).

Because even at these parameter values the binodals are only modestly shifted from the equilibrium values , the global phase volumes of the liquid and vapour phases remain primarily under the control of the reaction dynamics and are hence set by . At large enough values of this quantity, the majority liquid phase unavoidably percolates, but the minority phase of disconnected large vapour bubbles contain within them small liquid droplets that are stabilized by reverse Ostwald. (The interior of such a bubble thus resembles a piece of the liquid-in-vapour emulsion seen at smaller .) These small liquid droplets are continuously produced within the vapour bubbles but then grow, diffuse and merge into the surrounding liquid phase. This behaviour closely resembles the ‘bubbly phase separation’ reported previously for AMB+, except that in the latter case, there is only a single vapour domain in the system which can effectively be viewed as a bulk phase separation between the liquid-in-vapour emulsion and excess liquid [14]. (Note also that [14] mainly addresses the case with for which the identities of the ‘liquid’ and ‘vapour’ phases are, trivially, interchanged from those in the present discussion.) In the presence of the Model AB mechanism, this type of emulsion/liquid bulk phase separation is itself unstable since the bulk liquid phase is not at the target density. Accordingly, for slow reaction dynamics the system must homogenize once again at some larger scale so that the chemical conversion of the majority into the minority species in each neighbourhood can be balanced by diffusive mass transport. Accordingly we see a finite density of the large liquid-in-vapour emulsion bubbles whose size cannot grow further. Echoing the language of [14], we refer to this state as ‘bubbly microphase separation’.

A snapshot of the life-cycle of droplets is shown in Fig. (2) showing the growth of a small nucleated droplet, its growth, and coalescence into the surrounding liquid phase. As stated previously, this echoes the findings of [14] for the life-cycle of liquid droplets within a bulk phase separation between an emulsion of such droplets and excess liquid.

Three snapshots of the life-cycle of a dense droplet created inside the dilute phase, as indicated by the arrow. The simulation parameters are

Ostwald dynamics in model AB+

In the mean field limit, the binodals of AMB+, defined as the steady-state coexisting densities for a flat interface, can be computed analytically [13]. We will not show the calculation here but merely quote the results for : the dense binodal shifts to and the dilute binodal shifts to . These replace the values for the passive model (with as chosen above).

In addition, around circular two-dimensional droplets, the density close to the interface is further modified by an offset depending on the radius of the droplet R. In equilibrium Model B, this curvature-induced shift around a liquid droplet is proportional to , where is the surface tension. As shown in [14] for AMB+, the surface tension defined this way can be negative for sufficiently negative of and , thus reversing the Ostwald ripening for droplets, arresting their growth at a finite size. If we now switch on the chemical reactions to give Model AB+, then so long as these reactions are slow (which is the regime addressed in this paper) the dynamics on the length-scale of the interface width is dominated by the conservative current. This means that we can carry across from AMB+ the results for the curvature-induced offset to the binodals.

Consider N liquid droplets of radius for in a dilute bath in a domain of volume V with periodic boundary conditions. We now single out one droplet of radius and treat the rest as homogeneous background. Let be the composition field inside and outside the singled-out droplet, respectively, and define variables as the deviation from the modified binodal densities:For reasons given in [16], provided that is close to the modified dilute binodal (so that there is a large excess of vapour around the droplet), and so long as the timescale separation between conservative and nonconservative dynamics holds (slow reaction regime), we can linearize the deterministic part of Eq. (3) to obtain,where and we have defined for conciseness. The prefactor accounts for the change in the chemical turnover as a result of the N droplets. In the last term in the equation for , denotes the current at the surface of the jth droplet and thus this term represents the injection of mass by all the droplets. We refer to Sect. 3.3 of our previous paper [16] for a more detailed account of the calculation. Note that we neglect the term here as varies on a length-scale much larger than the interfacial width. The linearized equations need to be solved with the appropriate boundary conditions: at and , we only need to be finite; at the boundary of the droplet, we require , where is the aforementioned offset for the dense and dilute phase, respectively.

We will now sketch the perturbative Ostwald calculation with results presented only at key steps and refer to Sec. 3.2 of [16] for more detailed explanations of the scheme. This is inspired by the standard Ostwald calculation for both passive and active systems [14, 17] and its adaptation to multiple droplets [16, 19], but we implemented it here, for the first time, to include both the nonconservative reaction terms in the bulk phases (following [16]) and the effects of active currents on the matching conditions (following [14]) which cause reversal of the Ostwald dynamics.

For fixed droplet radii , we solve for the stationary state of the linearized equations with , then compute the current at the interface with the quasi-static solutions as a function of the droplet radii ,where are the modified Bessel functions, and . At the surface of the singled-out ith droplet (treating the rest as a homogeneous bath), the mismatch of the diffusive currents on the two sides of the interface shifts the droplet boundary, thus changing the radius of the droplet,So far, we have reduced the multiple-droplet dynamics to N coupled equations for by treating each droplet in isolation. One can then immediately spot a solution where all the droplets have the same radius. Let this radius be R, the currents can then be obtained explicitly by substituting the equation into the equation:where the arguments of the modified Bessel functions have been omitted to ease notation. Substituting Eq. (9) into Eq. (8) yields the growth rate as a function of R. The resulting curve is plotted in Fig. (3a) for one droplet in a box of size for and . Both plots show a stable fixed point: smaller droplets grow towards the stable radius, and larger droplets shrink towards it. These facts, respectively, explain the spontaneous emergence of liquid droplets within the vapour phase (given the presence of noise at short scales) and also the arrested growth of droplets shown in Fig. (1). There are, however, some important differences between the two plots in Fig. (3a). Observe that the curve for the currents inside the droplet is similar in the two panels, whereas the current in the vapour phase is qualitatively different. In the case where , as whereas for large R, leading to diffusive currents away from the droplet on the outside. In contrast, for , we have positive at infinity, thus the current outside the droplet is always toward the droplet. This has profound consequences: for , the zero of the has an upper bound at , which is independent of the reaction rate ; but for , there is no upper bound for the stationary fixed-point radius (at least at the mean field level). Comparisons of the above analysis against numerical simulations with low noise are shown in Fig. (3b): although the perturbative Ostwald calculation does not quantitatively capture the exact time evolution of R, due to the quasi-static assumption, the stable radii between theory and simulations agree reasonably well.

Fig. 3

a Plots of the growth rate and the currents inside and outside for a dense (liquid) droplet of radius R in a dilute (vapour) bath for and (recall ). In both cases, the curve crosses the X-axis from above, indicating a stable fixed point for R. The key difference is that in the left panel the curve is positive as whereas remains negative on the right panel. b Comparisons between theoretical predictions (dotted line) against simulations (solid line) for and . All plots are produced with the same reaction rate . The noise level used in the simulations is

As one increases the number of droplets in the Ostwald calculation (by changing in Eq. (9)), the stable fixed-point radius decreases, also observed in [16]; though qualitatively, the fixed-point radius is still substantially larger for than . This explains why, at fixed noise level, the emulsion of small and rather monodisperse liquid droplets seen in the bottom left panel of Fig. 1 takes on a quite different character in the two bottom centre panels, with relatively large and quite polydisperse droplets. Clearly in this high phase-volume emulsion state the droplets are neither uniform nor well separated, so that any further attempt to rationalize their statistics requires a much better understanding of the role of noise in the steady-state kinetics. This lies beyond our scope; as noted previously, the issue is a surprisingly complicated one even in AMB+ where chemical reactions are absent [14].

Fig. 1

Steady-state patterns for for two sets of values. For (as in Model AB), the patterns are stationary modulo noise fluctuations and are shown in the top panel. For , the steady states are dynamical with droplets of the dense phase (red) constantly spawn in the dilute phase (blue), snapshots of which are shown in the bottom panel

For even larger values of , the Ostwald calculation above continues to qualitatively capture the growth of initially small liquid droplets within vapour regions, but as the droplets grow larger, they come into contact and either form or merge into the continuous phase of liquid which is present in the resulting regime of bubbly or hierarchical microphase separation.

Conclusion

In this paper, we have furthered our investigation into scalar field theories for non-equilibrium phase separating systems. Specifically, we have presented an initial study of the effect of having diffusive dynamics that breaks detailed balance on its own, on top of a mismatch of chemical potentials between the conservative dynamics and the non-conservative dynamics. This mismatch is separately a hallmark of broken time reversal symmetry (broken detailed balance) in active phase separation whenever both types of dynamics are present at once. The resulting new model, Model AB+, incorporates the and terms, describing active currents in Active Model B+, into Model AB from our previous paper [16] which focuses on the chemical potential mismatch in isolation. The latter causes microphase separation whenever the chemical reactions drive the system towards a target density that lies between the binodals of a conserved, diffusive phase separation.

This study was motivated by the ability of the active current terms to reverse the classical Ostwald process [14], which is also associated with a distinct strong coupling regime renormalization group calculations [23]. As a result of this, Model AB+ has two distinct mechanistic channels for microphase separation. Our work has shown that these can interact in nontrivial ways, giving two distinct types of emulsion state and also a ‘bubbly’ microphase separation in which the two processes are operative hierarchically on different scales. To be more specific, we focused our studies on parameter regimes where the surface tension for droplets of the dense ‘liquid’ phase is negative. For close to the dilute ‘vapour’ binodal , suspensions of stable droplets of finite size are observed. Ostwald calculations of a single droplet reveal that if there is an upper limit on the droplet size as determined entirely by the diffusive dynamics, whereas for such an upper limit does not exist. The dynamics become more complicated for multiple droplets, though the qualitative conclusion remains that there are two structurally distinct liquid-in-vapour emulsions, as reflected in the phase diagram. On the other hand, for larger values of , the system settles into the bubbly microphase separated state, characterized by large but finite domains with densities at the two binodals, but with droplets of the liquid phase continuously created within the vapour domains before merging with the surrounding fluid.

Thus far, we explored some interesting behaviour of Model AB+, but made no attempt to systematically scan its multidimensional parameter space. There is, accordingly, plenty of room to find new physics beyond that presented here. For example, it would be interesting to see what is the interplay between the reverse Ostwald dynamics and the limit cycles observed in Model AB [16]. In these cycles, there is no stationary phase separation but a cycle where phase separation leads to a slow reduction in global density until the system rehomogenizes, whereupon the global density reverses until the system phase separates again. (This is possible when the local reaction rate is sufficiently nonlinear in density.) It is not clear how this oscillation might interact with a conserved dynamics that favours microphase separation, especially when this is itself a highly dynamical process governed by a nontrivial life-cycle for droplets. Finally, as reported separately for Model AB and AMB+, the noise strength plays an important role in selecting the number and size of droplets in the emulsion state [14, 16] and this aspect deserved further study.

Finally, we have not made any attempt to connect the parameters or mechanisms of Model AB+ directly with any specific microscopic examples of active phase separation, whether within living cells [1-4], bacterial colonies [11], or elsewhere. There have been successes in applying similar Swift–Hohenberg–Halperin-type models to ab-initio simulations and real experiments [25-27], but this is generally not easy because the model is constructed top-down to have a minimal structure. While the results are generic, only a comparison with more microscopic, bottom-up treatments can determine whether individual parameters are large or small in any particular case. For example, in purely conservative phase separation with ABPs (Active Brownian Particles), it was found that hard-core repulsion leads directly to the term in AMB+ but that to recover the term, additional (e.g., soft-core) interactions are also required [14]. Despite these difficulties, we believe that the exploration of a generic model such as Model AB+ can suggest mechanistic explanations for structure formation in binary active systems that might otherwise be quite puzzling. A possible example is when there is emergent structure on more than one length scale, as arises in bubbly microphase separation (bottom right panel of Fig. 1).