Critical transition between intensive and extensive active droplets
Jonathan Bauermann Giacomo Bartolucci Job Boekhoven Frank Jülicher Christoph A. Weber














Department of Physics, Harvard University, Cambridge, MA 02138, USA


Department of Physics Universitat de Barcelona, Carrer de Martí i Franquès 1-11, 08028 Barcelona, Spain



School of Natural Sciences, Department of Bioscience,
Technical University of Munich, Lichtenbergstraße 4, 85748 Garching, Germany




Max Planck Institute for the Physics of Complex Systems,
 Nöthnitzer Straße 38, 01187 Dresden, Germany


Center for Systems Biology Dresden,  Pfotenhauerstraße 108, 01307 Dresden, Germany


Cluster of Excellence Physics of Life, TU Dresden, 01062 Dresden, Germany 






 Faculty of Mathematics, Natural Sciences, and Materials Engineering: Institute of Physics, University of Augsburg, Universitätsstraße 1, 86159 Augsburg, Germany














§ ABSTRACT 
Emulsions ripen with an average droplet size
increasing in time.


In chemically active emulsions, 
coarsening can be absent, leading to a non-equilibrium steady state with mono-disperse droplet sizes. 

By considering a  minimal model for phase separation and chemical reactions maintained away from equilibrium, we show that there is a critical transition in the conserved quantity between two classes of chemically active droplets:
intensive and extensive ones.
Single intensive active droplets reach a stationary size mainly controlled by the reaction-diffusion length scales.
Intensive droplets in an emulsion interact only weakly, and the stationary size of a single droplet approximately sets the size of each droplet.
On the contrary, the size of a single extensive active droplet scales with the system size, similar to passive phases. In an emulsion of many extensive droplets, their sizes become stationary only due to interactions among them. 
We discuss how the critical transition between intensive and extensive active  droplets affects shape instabilities, including the division of active droplets, paving the way for the observation of successive division events in chemically active emulsions. 



                              




Critical transition
between intensive and extensive active droplets




    Christoph A. Weber
    September 9, 2024
===========================================================================








§ INTRODUCTION


Coarsening or ripening refers to the growth of larger domains at the expense of smaller domains that eventually shrink. For passive systems, the kinetics of coarsening stops when the system reaches  equilibrium, corresponding to a single domain in a finite system. Coarsening occurs in various systems, ranging from spin systems <cit.>, liquid emulsions <cit.>, and crystallized precipitates <cit.>. The kinetics of coarsening  is universal and determined by conservation laws, symmetries, and the dimension of the system <cit.>. 

For active systems persistently maintained away from equilibrium <cit.>, the kinetics of coarsening is altered and can be even suppressed <cit.>. The paradigm is reaction-diffusion systems that give rise to non-equilibrium steady state patterns with various spatial morphologies <cit.>. Another example is Model B+, which, in contrast to the classical Model B, also accounts for contributions to the diffusive fluxes that do not arise from free energy <cit.>. These fluxes give rise to anti-coarsening with a condensed phase that stopped growing and a “bubbly”  morphology of material-poor domains <cit.>. Finally, suppressed ripening was also observed in liquid-liquid phase-separated systems with chemical reactions maintained away from equilibrium <cit.>. These systems are also called chemically active emulsions <cit.>.  

The formation of steady-state patterns in reaction-diffusion systems relies on the reaction flux that breaks the detailed balance of the rates. Together with diffusion, this gives rise to various reaction-diffusion length scales that are crucial but not exclusively responsible for pattern morphology. Chemical processes generically come with conservation laws for mass, and if incompressible, also for volume. It has been shown that conservation laws are key determinants for the emerging patterns in mass-conserving reaction-diffusion systems <cit.>. Key implications are that the pattern-forming transition is typically sub-critical, reminiscent of discontinuous phase transitions <cit.>. 

Active emulsions also give rise to reaction-diffusion length scales. These length scales are crucial for various non-equilibrium phenomena in active emulsions, such as dividing droplets <cit.>, formation of steady liquid shells <cit.>, and the suppression of Ostwald-ripening <cit.>. Despite such interesting phenomena, a limitation of some minimal models <cit.> is that two components were considered lacking  conservation laws. Similar to reaction-diffusion systems <cit.>, conservation laws can qualitatively alter the nature of the transition and instabilities in active emulsions. 

In passive, phase-separated systems  with chemical reactions, the reaction-diffusion length scales do not determine the equilibrium state.
In this passive case, conservation laws, i.e., so-called lever rules for quantities conserved by the chemical reactions, determine the volume of the condensed phase(s) at thermodynamic equilibrium. When the chemical reactions are maintained away from equilibrium (active emulsion), the emerging reaction-diffusion length scales compete with the conservation laws. It remains unclear whether such reaction-diffusion scales or the conservation law dominate the selection of pattern length scales in active emulsions.

Here, we study the role of conservation laws in multi-component mixtures for the droplet size and the droplet number in monodisperse active emulsions that do not undergo coarsening.  
Our key finding is that the conserved quantity controls a critical transition between 
intensive and extensive chemically active droplets that differ in their physical behavior when increasing the system size; see Fig. <ref> for an illustration. In the case of intensive droplets, single droplets in large systems are stationary, with droplet sizes independent of system size. 
For extensive droplets, the stationary droplet size increases with the system size. 
We study the consequences of this transition for the collective dynamics of many droplets in emulsions and show how monodispersity can arise for both classes of chemically active droplets.







§ MINIMAL MODEL FOR AN ACTIVE EMULSION WITH A CONSERVED QUANTITY




 §.§ Dynamics of the concentration fields


An incompressible binary mixture with two molecules A and B converting into each other via a chemical reaction has a conserved quantity, the total mass of A and B. However, in a binary mixture, this conserved quantity is trivial as it is constant in space and thus does not change dynamically <cit.>. 
A mixture containing chemical reactions must comprise at least three different molecules to have a conserved quantity that can vary in space. 
In the following, we introduce a minimal model for such a  ternary mixture that is incompressible and composed of a non-reacting solvent S and the molecules of types A, B, that can react via the following reaction scheme

    A k_BA B  ,    
    B k_AB A

For this chemical reaction,
the dynamics of average concentrations ϕ̅_i(t)=V^-1∫_V d^3x  ϕ_i(x,t) (i=A,B,S) are constraint to a specific conserved line in the ϕ̅_A-ϕ̅_B-plane which can be expressed as

    ψ̅  =  (ϕ̅_A+ϕ̅_B)/2

where V denotes the volume of the system and ϕ_i(x,t) are the concentration fields  that depend on position x and time t. The conserved quantity is denoted as ψ̅.

The dynamical equations for the fields ϕ_A(x,t) and ϕ_B(x,t) are
  

    ∂_t ϕ_A    =   ∇·( Γ_A ∇μ_A )+ k_ABϕ_B  - k_BAϕ_A  
    ∂_t ϕ_B    =  ∇·( Γ_B ∇μ_B ) - k_ABϕ_B  + k_BAϕ_A
 


where μ_i =  δ F /δϕ_i
are the chemical potentials of i=A,B, descending from the Helmholtz free energy F=∫_V d^3x ( f_0+ κ_A (∇ϕ_A)^2/2 + κ_B (∇ϕ_B)^2/2 ) <cit.>. Here,  Γ_i denote the mobility coefficients.

The first terms of Eqs. (<ref>)  are diffusive fluxes that are driven by local gradients of the chemical potentials, while the remaining terms, k_ABϕ_B  - k_BAϕ_A, are chemical fluxes. For passive systems, 
the reaction rate coefficients k_ij are concentration dependent such that the free energy F  determines  thermodynamic equilibrium <cit.>. If, however, k_ij are also dependent on some external reservoir or are chosen independently of F, the system is inherently chemically-active <cit.>.

Without chemical reactions, the dynamic system governed by Eqs. (<ref>) relaxes to phase equilibrium, characterized by equal chemical potentials μ_i^I=μ_i^II, and equal osmotic pressures 
Π^I=Π^II, where Π= -f+ ∑_i=A,Bϕ_i μ_i, see Ref. <cit.> for a general introduction. 
Such phase equilibria are also relevant for the dynamic system with chemical reactions since
they locally govern the dynamics of the interface.  
The phase equilibria of the mixture are portrayed by the phase diagram. 
In ph_d(a), we show sketches for three different mixtures where
A always phase separates from S when B is absent. 
The molecular interactions of B with A and S, shape the phase diagram. 
In the left panel of ph_d(a), we show a sketch of a phase diagram where B molecules are similar to S and phase separate from A; in the central panel, B molecules that are neutral and do not interact differently with  A or S. In the right panel B molecules are identical to A molecules and therefore phase separate from the solvent. 

To mimic these interaction characteristics and to avoid tying our conclusions to specific free energy models, in the following, we take a geometrical view on phase diagrams.
Phase equilibria are characterized by the slope of the tie lines in the phase diagram. 
The relevant geometrical property of the phase diagrams that we consider is the local angle α between the tie lines formed with ϕ_A-axes of the phase diagram (ph_d(b)).

This local, geometrical view on a phase diagram allows us to consider free energies that are mathematically simpler even than mean field theories such as the Flory-Huggins theory, allowing us to obtain analytical results for stationary droplet states.
For ease of presentation, we consider cases where tie lines are locally orthogonal to both the dense and dilute branches of the binodal line. 
These geometrical properties can be captured by introducing the shifted concentration field ϕ_i →ϕ_i-ϕ̅_i and the following Ginzburg-Landau-like free energy density f_0:


    f_0(ϕ_A, ϕ_B;α)    = b_1/2 (ϕ_1 +  )^2 (ϕ_1 -  )^2+ b_2/2ϕ_2^2   ,
    [ ϕ_1; ϕ_2 ]   =
    [  cosα  sinα; -sinα  cosα ][ ϕ_A; ϕ_B ]


where b_1 and b_2 are parameters characterizing interactions and entropic contributions. 
The qualitative cases shown in ph_d(a) can be described 
by varying the angle α that parameterizes the rotations of the energy density f_0.
Specifically, different values of α correspond to different  types of interactions of the B molecules with A and S.  
Without loss of generality, we restrict ourselves to the domain -π/4 ≤α≤π/4. Indeed, due to the symmetry of the free energy, the transformation α' = α + π/2 is equivalent to re-labeling A to B and B to A. 
For α = -π/4, B and  A interact equally with the solvent, for α = 0, B does not interact with A and solvent, and for α = π/4, B and S interact equally with A. 
As a consequence, in the case α = -π/4, the two reactants A and B localize in two distinct phases, i.e., they segregate.  In contrast, for α = +π/4, molecules of type A and B   phase separate together from the solvent S. Following Refs. <cit.>, we refer to α = -π/4 as the segregative case and α = +π/4 as associative case. In the literature on mixtures composed of a large number of components, these two cases are often called demixing and condensation, respectively <cit.>. 
In the literature of biomolecular condensates, the shape of the phase diagrams is often related to the underlying molecular interactions.
According to Ref. <cit.>, heterotypic interactions correspond to the segregative case and homotypic interactions to the associative case.
Since the molecular interactions determine the composition of the coexisting phases,  we refer to  α as compositional angle.
Varying α in the range -π/4<α<π/4  interpolates  between the limits of the segregative and associative case. 


In ph_d(b), we show the phase diagram associated with the free energy density f_0 in f. 
Note that the binodal lines (dark blue) separating the mixing from the demixing regimes (white and blue shaded area, respectively) are straight lines tilted by the compositional angle α from the vertical. 
The tie lines (represented as dashed blue lines) are parallel to each other, perpendicular to the binodal, and thus tilted by the compositional angle α from the horizontal.
Finally,  we comment on the role of ϕ_0 in f (ph_d(b)): 
the value ϕ_0 sets the scale of the concentration axis since for cos(α)ϕ_A + sin(α) ϕ_B inside the interval [-ϕ_0, ϕ_0], the passive system can phase separate.


In the following, we consider constant chemical rate coefficients k_ij that are independent of the free energy F. As a result, the system is chemically active since in general no phase equilibrium of coexisting phases exists where the chemical reactions are at a steady state in all phases. 
These homogeneous steady states of chemical reactions are governed by the condition 

    k_ABϕ_B  = k_BAϕ_A

 defining the reaction nullcline shown as a solid green line in ph_d(b).
We refer to the angle β between tie lines and the reaction nullcline with

    β = arctan(k_BA/k_AB) - α

 as the activity parameter. This is because it
 characterizes the strength of non-equilibrium driving.
Indeed, only for the special case β = 0, 
the system is passive and can settle to thermodynamic equilibrium because the coexisting concentrations (phase equilibrium) are also connected by the 
reaction nullcline (Eq. (<ref>)). Therefore, for β=0, the phase equilibrium and chemical equilibrium coincide.  
For β≠ 0, coexisting concentrations connected by a tie line do not lie on the reaction nullcline. 
Thus, no stationary state of reactions (reaction_nullcline) simultaneously fulfils phase equilibrium. Instead, reaction fluxes and diffusive flux between the phases balance each other leading to a non-equilibrium steady state. 
The activity parameter is restricted to the domain 0 ≤β < π/4 because for β = π/4 the reaction nullcline becomes parallel to the tie lines.

In summary, the active emulsions with the chemical reaction (<ref>) 
is characterized, when using a geometrical representation of the phase diagram (f), by the 
compositional angle α, the activity parameter β, and the conserved quantity ψ̅ (Eq. (<ref>)).



 §.§ Dynamics of the conserved and non-conserved fields



To understand the role of conserved quantities in the ripening kinetics, we introduce the conserved field ψ(x,t) and the reaction extent field ξ(x,t), which are defined as: 

    ψ(x,t)   =(ϕ_A(x,t) + ϕ_B(x,t))/2 
    ξ(x,t)   = (ϕ_A(x,t) - ϕ_B(x,t))/2


For the Ginzburg-Landau type of free energy 
(f), we consider a constant mobility for simplicity and also write Γ_i = Γ, leading to: 

  

    ∂_t ψ   = Γ∇^2 μ_ψ  (ψ,ξ; α)   
    ∂_t ξ   = Γ∇^2 μ_ξ  (ψ,ξ; α) 
     
          -K ξ + K (π/4 +α+β) ψ
 

where we have introduced the chemical potentials of the conserved quantity, μ_ψ =(μ_A +μ_B )/2, and the one of the non-conserved quantity, μ_ξ = (μ_A-μ_B)/2. 
When recasting the free energy f_0(ψ,ξ) in terms of the conserved and non-conserved fields, these chemical potentials can be obtained from functional derivatives, μ_ψ=δ F/δψ and μ_ξ=δ F/δξ.
Furthermore, we define the overall reaction rate as

    K= k_AB + k_BA

In Appendix <ref>, we show how 
the rate coefficients k_AB and k_BA can be expressed in terms of the overall rate K, the activity parameter β and the compositional angle α.


We can identify two special cases where the chemical potentials of the conserved field and of the reaction extent field decouple, i.e., the free energy density is devoid of coupling terms proportional to (ψ ξ): 
(i) For the segregative case  α=-π/4, we find f_0(ψ,ξ) =  b_2 ψ^2/4 + b_1(ξ-√(2)ϕ_0)^2(ξ+√(2)ϕ_0)^2/8. Thus, the free energy density contribution in ξ is a classical double-well potential, which leads to a Cahn-Hilliard dynamics with reactions for the reaction extent field ξ( x, t) (dyn_xi). The quadratic term of ψ to the free energy density gives rise to Fick's law of diffusion for the conserved field ψ( x,t) ((dyn_psi)).
(ii) On the contrary, in the associative case, α=π/4 gives a free energy density   f_0(ψ,ξ) = b_1(ψ-√(2)ϕ_0)^2(ψ+√(2)ϕ_0)^2/8 + b_2 ξ^2/4. In this case, the free energy density leads to a Cahn-Hilliard dynamics for the conserved field ψ( x,t) (dyn_psi),  while the reaction extent field ξ( x, t) follows Fick's law of diffusion with linear reactions ((dyn_xi)). 
Though the dynamics of the conserved field ψ evolves independently from the reaction extent field in both cases (i) and (ii), 
the dynamics of the reaction extent is affected by the conserved field ψ( x,t) via a source or a sink term in the reaction dynamics (dyn_xi).

 
The ripening dynamics of these two special cases α=±π/4 is known and has been excessively studied in the past:  
(i) For the segregative case α = -π/4 (A and B segregate in different phases), the conserved variable ψ follows a simple diffusion dynamics that is agnostic to the chemical reactions and the dynamics of the reaction extent ξ.
Thus, the field ψ settles in a homogeneous state  on long-time scales. 
Moreover, the dynamics of ξ is the dynamics of a binary mixture composed of A and B  with active chemical reactions between these two components. This leads to the suppression of Ostwald ripening, i.e., a mono-disperse emulsion with a stationary droplet size <cit.>.
Once droplet size exceeds this stationary value,  
droplets shrink or can undergo shape instabilities, giving rise to droplet division <cit.>, or form stationary, active shells <cit.>.

(ii) For the associative case α = π/4 (A and B phase separate together), the conserved variable ψ follows a classical Cahn-Hilliard dynamics. Thus, droplet-like domains, e.g., rich in ψ, undergo Ostwald ripening in a sea of low ψ, where bigger droplets grow at the expense of smaller shrinking ones that eventually disappear.
On large time scales and in a finite system, the system evolves toward a single droplet with a volume that scales with the system size V.  
The dynamics of the reaction extent ξ 
is affected by the ψ field without disturbing it. 
A similar decoupling occurs in models for reacting diluted “client” particles corresponding to ξ and phase separating scaffold components corresponding to ψ <cit.>.



In the remaining general cases -π/4<α<π/4, the dynamics of the conserved field ψ and the reaction extent field ξ are coupled in the free energy. 
The key property as compared to the special cases (i,ii) is that the dynamics of the conserved field ψ is influenced by the non-conserved field ξ. This coupling is a generic consequence of the different interactions among the molecular constituents of the mixture (see ph_d(a,b) and related discussions). 

Therefore, the active  chemical reactions in the dynamic equation for the reaction extent ξ (dyn_xi) also affect the dynamics of the conserved field ψ (dyn_psi). 
In the next section, we will show that, strikingly, this remains relevant in the thermodynamic limit of large systems (V →∞), i.e., on length scales much larger than the reaction-diffusion length scale. 
Note that solely the non-conserved field ξ gives rise to a reaction-diffusion length scale λ = √(D/K), which can be obtained when linearizing near phase equilibrium. Here, D = Γ∂^2 f_0(ψ^±, ξ^±)^2/ ∂ξ^2 which can be calculated in our model:
D = √((Γ[ b_2+ 4 b_1ϕ_0^2 +(b_2- 4 b_1ϕ_0^2)sin(2α)])/4).
The reaction-diffusion length scale λ is finite and real (for Γ>0, K >0) for all values -π/4 ≤α≤π/4 and independent of β and the conserved quantity ψ̅.



§ STATIONARY SINGLE DROPLETS


To study the effects of the conserved quantity  ψ̅ on the droplet size, we first consider a single spherical droplet in a finite, spherically-symmetric system of radius R_sys=( 3V/(4π))^1/3, where V is the system volume. 
We consider a sharp interface limit, which is  valid when the droplet size is large compared t exceeds the width of the interface <cit.>.
Phase equilibrium is imposed at the position of this sharp interface. This boundary condition couples the reaction-diffusion equations in each phase. These  reaction-diffusion equations can linearized near phase equilibrium at the interface;  for details, see Appendix <ref> and Ref. <cit.>. 

In the following, we discuss the stationary solutions for a single droplet, i.e., ∂_t ψ(x,t) =0 and ∂_t ξ(x,t)=0, and a stationary interface position R_stat. 
For large enough system sizes R_sys, we always find two solutions for R_stat (solid and dashed lines in system_size(a)) at which the total reaction flux of every component in one phase is perfectly balanced by the diffusive flux of the same component over the interface. 
The lower branch (dashed) of these two solutions is the critical nucleation radius corresponding to an unstable fixed point.   The upper branch (solid) of R_stat is a stable fixed point corresponding to a non-equilibrium steady state. 
See Appendix <ref>, analytics_Rstat, for an analytical expression of the stationary radius as a function of the activity parameter β and the compositional angle α for large system sizes.



The conserved quantity ψ̅ (cons_quantity) affects the non-equilibrium steady state and leads to a changed behavior of the stationary radius R_stat as a function of the system size R_sys.
This changed behavior is depicted in system_size(a) which shows  R_stat(R_sys) for two different values of the conserved quantity ψ̅ = -0.61 (blue) and ψ̅=-0.425 (red) in a system with α=0 and β=π/8. 
Interestingly, in the case of ψ̅= -0.61 (blue), the stable stationary solution converges to a finite value in the limit of large systems, i.e., R_sys→∞, while for ψ̅= -0.425 (red), the solution scales linearly with the system size, for large systems. In other words, there are two different cases: 
(i) Intensive active droplets (blue) for which the stationary droplet size R_stat quickly saturates with system size R_sys implying that R_stat 
is set by molecular and kinetic parameters for large system size. 
(ii) Extensive active droplets (red) where the stationary droplet size R_stat increases linearly with system size R_sys implying that R_stat is set by R_sys.  


To explore these two behaviors, we numerically solved the dynamic equations with a continuous interface (Eqs. (<ref>)) in three dimensions. 
system_size(b) shows the field ϕ_A(x) in four different stationary states, corresponding to the two different cases (i,ii) and two different system sizes (L=600 and L=150 in inset). 
We consider the same values of the conserved quantity as in  system_size(a).
In the case of ψ̅=-0.61 (blue), the radius of the stationary active droplet for the larger system is almost identical to that of the smaller system in the inset. 
In the case of ψ̅=-0.425 (red), however, the stationary droplet radius in the larger system is much larger than in the smaller box. 
We conclude that the solutions of the continuous dynamic Eqs. (<ref>) confirm the trend obtained from our single droplet study in the sharp interface limit discussed in  system_size(a). 

The qualitative change in the scaling of the stationary radius R_stat with system size R_sys 
upon changing the conserved quantity  ψ̅ suggests that there is a bifurcation between the regimes of intensive (blue) and extensive (red) active droplets. 
To determine the nature of this bifurcation, we use the sharp interface model of a single droplet and consider the limit of an infinite system; details see Appendix <ref>.
We calculate the stable (solid) and unstable (dashed) stationary radii R_stat as a function of the conserved quantity ψ̅ (system_size(c)). 
We find that R_stat diverges at the threshold values (see Appendix <ref> for details on the derivation)

    ψ̅_crit^± = ±ϕ_0/2( cos(α) + sin(α) )
 
In the vicinity of these critical values, the stationary droplet radius diverges with R_stat∝ |ψ̅-ψ̅_crit^±|^-1 within the blue domain, see App. <ref>.

Furthermore, stationary active droplets can only be found within a certain range of the conserved quantity ψ̅∈ [-ψ̅_bin,ψ̅_bin] (dashed-dotted black line in system_size(c)) with

    ψ̅_bin^± = ±ϕ_0/2cos(α+β) + sin(α+β)/cos (β)

This range in the conserved quantity graphically corresponds to steady states for which the reaction nullcline lies within the binodal domain of the phase diagram (ph_d). 
Moreover, for -ψ̅_bin<ψ̅<0, we find A-rich droplets in a solvent-rich phase, while for 0<ψ̅<ψ̅_bin, we find solvent-rich droplets in a A-rich phase.

By means of Eq. (<ref>) and Eq. (<ref>), we can study how the underlying molecular interactions (parameterized by the compositional angle α) and the activity parameter β affect the existence of droplets and the critical bifurcation; see Fig. <ref>(a,b).
Neglecting the finite size effects, 
these results are independent of the specific values of kinetic parameters such as Γ and K, and otherwise only depend on thermodynamic parameters in the free energy (f). 
We find that for β>0, the active driving leads to a regime of intensive droplets for all values -π/4 ≤α< π/4, see Fig. <ref>(a). Recall that in the segregative case, α = -π/4 corresponds to an effective binary model of components A and B, which phase separate from each other and are converted into each other via chemical reactions; this case was studied in Ref. <cit.>. Only in this limiting case are droplets always intensive, and the conserved quantity does not affect the stationary radii of chemically active droplets. For all values -π/4<α< π/4, there is an additional domain with extensive droplets. The transition between these two regimes occurs at a value of the conserved quantity  expressed in Eq. (<ref>). 
In the associative case α=π/4, =, i.e., no intensive droplets can be found, independently of the value of the activity parameter β. In this case, droplets behave like droplets in passive systems despite the presence of an active chemical reaction. Chemical reactions still drive diffusive fluxes of A and B between the phases. However, these fluxes do not affect the spatial distribution of the stationary conserved field  ψ(x).
In the passive case with β=0, we find that ψ̅_bin = ψ̅_crit independent of the compositional angle α, see Fig. <ref>(b) leftmost panel (green domain). Indeed, passive droplets are always extensive since they behave as thermodynamic phases that scale with the system size. 
However, for β>0, ψ̅_crit^±≶ψ_bin^±. Therefore, intensive droplets of finite size exist inside a domain that increases for larger activity parameter β (see Fig. <ref>(b) from left to right).

In summary, our analysis of single active droplets shows that there are two different classes of stationary states in the limit of an infinitely large system: 
intensive active droplets and extensive active droplets. While an intensive active droplet adopts a finite size in a large system, an extensive active droplet scales with the system size. 
There is a critical transition between both types of non-equilibrium steady states. 
This transition is controlled by the interactions among the components, characterized in our model by the compositional angle α, the conserved quantity ψ̅, and the activity parameter β. 
In the next section, we  study emulsions composed of many active droplets and explore their behavior in the parameter regimes corresponding to intensive and extensive active droplets. 





§ INTERACTIONS OF MANY DROPLETS
 


In this section, we study the dynamics in active emulsions by numerically solving the dynamic equations (dyn_phi) for the continuous fields ϕ_A(x,t) and ϕ_B(x,t)  with a continuous interface in three dimensions.
We initialize the system in a homogeneous state at local chemical equilibrium, with the average concentration values ϕ̅_A and ϕ̅_B inside the binodal domain but outside of the spinodal domain. 
Moreover,  N_init small spherical droplets above the critical nucleation radius are randomly positioned in the system with ϕ_A=ϕ_0 inside (Eq. (<ref>)). 
To avoid the fusion of such initially placed droplets right after initialization, we only accept random configurations with inter-droplet distances above a threshold value. 

We now discuss representative results corresponding to the intermediate values of the compositional angle α=0 and the activity parameter β=π/8.
With this choice, psi_bin and psi_crit imply that the bounds on the conserved quantity ψ̅ for the existence  of droplets are given by ψ̅_bin^± = ±ϕ_0/√(2) and the critical bifurcation between intensive and extensive droplets is ψ̅_crit^± = ± 0.5  ϕ_0. 
The initial concentrations outside are chosen such that the total conserved quantity in the system is either ψ̅=-0.6  ϕ_0, or ψ̅=-0.4  ϕ_0. 
These two values lie, respectively, in the regime where chemically  active droplets are intensive (ψ̅<ψ̅_crit^-), or extensive (ψ̅_crit^-<ψ̅<ψ̅_crit^+).


The left-most panels of many_d(a) (grey) depicts the initial state of the active emulsion  corresponding to an initial droplet number of N_init = 8  and N_init = 96, respectively. 
In the remaining columns of many_d(a) (blue and red),
we show the
the stationary states of the active emulsion for these two initial conditions and the aforementioned conserved quantities  ψ̅=-0.6  ϕ_0 (blue) and ψ̅=-0.4  ϕ_0 (red). 
We find that for ψ̅=-0.6  ϕ_0 (intensive regime), all droplets grow to the same size, i.e., a mono-disperse active emulsion emerges. 
During this kinetics to the stationary state, a few droplets vanish compared to the initial condition. 
Interestingly, for  ψ̅=-0.6  ϕ_0, the stationary droplet radius changes only slightly for the two initial conditions. 
In the right column of many_d(a) (red), we show the stationary states for ψ=-0.4  ϕ_0 (above the critical transition, extensive regime) for the same initial conditions as above. For the initial state composed of eight droplets only (first row), droplets elongate and sometimes branch, thereby forming tube-like structures. However, some of these branched structures can separate during the dynamics such that, in the end, N=17 different tubes of different lengths exist. 
When initializing 96 droplets (second row), some of them elongate, while others stay almost spherical, dependent on the initial distance to other droplets. However, none of them dissolve, such that we find N=96 objects in the stationary state. When even more droplets are initialized, some of the initial droplets dissolve, while the rest stay spherical and try to maximize their distance towards each other, see  many_d(d).

Despite these complex changes in droplet morphologies, there is a common principle in the regime where single active droplets are  extensive (red): 
Despite varying the initial conditions and different numbers of droplet-like domains in the stationary state, the total phase volume is almost constant (many_d(b)). 
In other words, when there are more droplet-like domains, they are smaller such that the total phase volume is approximately conserved. 
In contrast, in the regime where active droplets are intensive (blue), the total volume increases with the  amount of initialized droplets. Most importantly, in this regime, the average droplet volume is roughly constant (many_d(c)) and approximately equal to the  volume of one single droplet up to finite size effects.


Monodisperisty in chemically active emulsions can be explained as follows:
Intensive droplets (blue) adopt a fixed size independent of the system size. When just a few droplets are initialized in a large system, they can be far apart, and whenever they are separated by a distance much longer than the reaction-diffusion length scale in the outside phase, they hardly interact with each other. As a consequence, their dynamics becomes stationary and their average radius takes a value similar to the stationary radius R_stat obtained from the single droplet analysis (see Sect. <ref>).
When more droplets are initialized, such that they are closer than this reaction-diffusion length scale, they weakly interact and get stationary at smaller sizes. This leads to the slight decrease of average droplet volume for large N, seen in many_d(c).
In all cases, droplets reach identical radii.
For extensive droplets (red), monodispersity is only found when many droplets are initialized. Indeed, when only a few droplets are initialized, they tend to grow until they feel the presence of their neighbors. The more droplets are initialized, the sooner this arrest occurs, decreasing the radii of droplets. 
If the number of initialized droplets is too small, droplets grow to sizes where shape instabilities can occur. 
Tubes form via an elongation instability, discussed in the next section, or spherical shells form via a spinodal instability at the center of large droplets <cit.>. 
Crucially, when enough droplets are present, the arrest of growth happens earlier than shape instabilities can occur. As a consequence, droplets are spherical and reach identical radii (many_d(d)).
In summary, for extensive droplets in an emulsion (red),
the history of the emulsion matters for the final stationary mono-dispersed state or the emerging shape instabilities. 


To confirm the validity of these arguments relying on a single droplet analysis,  we compare the phase volume and the average droplet size with the analytic results obtained from the single droplet case in the sharp interface limit.
To this end, we considered N identical subsystems arranged in a hexagonal close-packed lattice that fills the total system volume V_sys, and calculated phase volumes and average sizes. 
These analytically derived results (solid lines in many_d(b,c)) are in good agreement with the measured values in the numerical simulations of the active emulsions. 
In the case of extensive droplets (dashed lines in many_d(b,c)), the total volumes of the phases are well captured by considering the limit of fast diffusion <cit.>. 






§ SHAPE INSTABILITIES AND DROPLET DIVISION



In the previous section, we have seen that for extensive active droplets in an emulsion, an elongation instability can occur. 
This instability is reminiscent of the Mullins-Sekerka instability that occurs in the diffusive growth of a single droplet without reactions in a large system. 
It requires a (weak) deformation of a droplet's spherical shape and relies on the material depositing at interfacial domains of larger mean curvature.
As a result, these domains grow faster, further deforming the droplet shape. 
Deformations are counteracted by surface tension, which tends to flatten the mean curvatures of the interface. However, this  effect weakens the larger the droplet becomes. Thus, there is a critical radius above which droplets constant growth and deform (Mullins-Sekerka instability).


A similar instability can occur for both extensive and intensive active droplets.
Indeed, in a binary mixture where droplets are always intensive, such instability was shown to lead to droplet division <cit.>.
Here, we show that droplet dynamics after elongation vastly differ between intensive and extensive active droplets.
In the upper row of division (blue), we show snapshots of the dynamics of two intensive droplets (ψ̅=-0.6 ϕ_0) and compare it to the dynamics of two extensive droplets (ψ̅=-0.4 ϕ_0) in the lower row (red). For both cases, we choose identical kinetic parameters.
However, we have chosen the overall reaction rate K (overall_reaction_rate) such that the stationary droplet radius in the intensive regime exceeds the critical radius for onset of the elongation instability.
In both cases, the two droplets are initially separated by a distance smaller than the reaction-diffusion length in the outside phase. Thus, the presence of the neighboring droplets leads to asymmetries that induce initial deformations, which get amplified by the elongation instability.
However, the division predominantly happens for intensive droplets; see division (upper row). Extensive droplets elongate (division lower row) and form tubes, as seen before.
This difference can be explained by focusing on the neck region. For intensive droplets, there is almost no growth in this region, such that a Plateau-Rayleigh instability can lead to a pinch-off. However, extensive droplets constantly grow, including in the neck region, thereby inhibiting the Plateau-Rayleigh instability.

Interestingly, with such cycles of growth and division of intensive droplets, the total phase volume can trivially grow just by increasing the number of droplets. As a result, these shape instabilities allow intensive droplets to fill the space. Thus, macroscopic phase volumes are reached, even in infinite systems, similar to the regime of extensive droplets.



§ MOLECULAR INTERACTIONS AFFECT THE ARREST OF RIPENING



In section <ref>, we have shown that extensive droplets grow until the interaction with neighboring droplets arrests their growth. However, in section <ref>, we discussed that in the associative case (α = π/4) where A and B have identical interactions with the solvent, the conserved field decouples from the active driving, leading to passive phase-separating dynamics. 
In this section, we study how the emulsion dynamics is affected as molecular interactions approach the limit of the associative case (α→π/4). We find that extensive droplets can initially ripen like in the passive case until the effects of the active driving manifest and ripening arrests.  
 
In ripe(a), we show snapshots of the dynamics of concentration field ϕ_A for four different values of the compositional angle α; see Supplementary Material for videos on the ripening dynamics. All numerical calculations are conducted with the same kinetic parameters, an activity parameter β=π/4, and a conserved quantity of ψ̅ = -0.35 ϕ_0. For these settings, all systems are populated by extensive droplets.
We initialize at t=0 (not shown in the figure), the concentrations homogeneously at chemical steady state with small fluctuations (around 0.1% of ψ̅). 
For the values of the compositional angle α considered, the system is within the spinodal regime of a non-reacting phase separating mixtures. For triggering the nucleation of many droplets, seen in the first row in ripe, we did not allow for chemical reactions in the early times of the spinodal instability (t≈ 2). Insights on the change of the nucleation dynamics in chemically active systems can be found here <cit.>.
While the nucleation process is similar for all the compositional angles α considered here, the long-time ripening dynamics depend strongly on α. 
For α=0.19  π (red/first row in ripe(a)), some of the initial droplets dissolve, but the interaction between droplets causes the ripening to arrest.
As the interaction of A and B with the solvent gets more similar, i.e., α→π/4 (associative case), this arrest of ripening occurs at a later stage of the dynamics. Prior to the arrest of ripening, smaller droplets dissolve while larger ones grow. Thus, as the arrest occurs and the emulsion becomes monodisperse, stationary droplet radii are larger as α→π/4 (associative case), see second (α = 0.23 π) and third (α = 0.24 π) row of ripe(a), in red. Only in the associative case  α = π/4, the arrest is absent and chemically active droplets ripen as passive emulsions (green/last row in ripe(a)). This means that the classical power-laws of Ostwald ripening hold in the associative case (α = π/4) with ⟨ R ⟩∝ t^1/3 and N∝ t^-2/3 <cit.>,
while for decreasing values of α < π/4, the transition to such classical ripening laws happens at later stages of the dynamics (see ripe(b,c)).
 
 



§ DISCUSSION


We have introduced and studied a model of a minimal ternary mixture of the components A, B, and solvent that can phase separate while A and B convert into each other. 
This conversion is maintained away from chemical equilibrium 
by choosing reaction fluxes independent of the chemical potentials. 
This choice prevents the system from simultaneously reaching chemical and phase equilibrium and allows it to adopt a non-equilibrium steady state.
The model is minimal because it exhibits a non-trivial conserved quantity of the reaction with an associated conserved field that is generally position-dependent, also in the non-equilibrium steady state. 
Within this minimal model for active chemical reactions in phase-separating mixtures, we found a critical transition, determining whether the size of a single droplet scales with system size (extensive active droplet), like in passive systems,
or whether it reaches a finite, characteristic size (intensive active droplet). 

The transition found in the single droplet analysis has crucial consequences for the behavior of emulsions composed of many droplets. For a value of the conserved quantity in the regime of an intensive droplet, a mono-disperse emulsion forms where the average size is approximately given by the stationary size of a single droplet. The reason for this approximate agreement is that  droplets interact only weakly in the intensive droplet regime. 
In contrast, extensive droplets in an emulsion constantly grow until they come close enough to interact strongly with each other. Once they interact, there can be initial ripening, but growth arrests. 

The morphology of the resulting structures depends crucially on initial conditions, i.e., the history of the emulsion.  In other words, few droplets deform via shape instabilities, while many droplets remain spherical and monodisperse. However, the  total phase volume is independent of the number and sizes of droplets and scales with the system volume similar to passive thermodynamic phases.

Our findings highlight the importance of several components in chemically active systems.
For example, accounting for a solvent component fundamentally alters the dynamics of chemically active droplets, reflected in the  appearance of the critical transition when varying the conserved quantity. 
We propose that understanding the effects of the conserved quantities(s) in reacting system away from equilibrium is crucial to correctly interpret emulsion dynamics in biological and chemical systems <cit.>.
Specifically, our theory can unravel the conserved quantity's effects in synthetic chemical systems, such as droplet size control in multi-component mixtures leading to chemically active coacervates <cit.>. 
According to our theory, droplet division is more likely to occur in  the regime α≃π/4  corresponding to intensive active droplets and the case where A segregates from B from each other (segregative case). 
We propose  studying a single active droplet in a reaction container for different container sizes <cit.> and test whether or not
the active droplet scales with the container size (intensive active droplets, see Fig. <ref>(a)).
In the regime of intensive active droplets, division is more likely to be obtained.  This procedure could be used to test synthesized actively reacting components for their division propensity and may thus pave the way to observe successive division events in larger, mass-supplied containers containing many active droplets. 






§ APPENDIX



§ CHEMICAL REACTION RATES



In the following, we will use the overall rate K = k_AB + k_BA and the activity parameter β as parameters for the chemical reactions. These parameters are related to the forward rate k_BA and backward rate k_AB, appearing in Eqs. (<ref>) via the following relationships: k_BA = K sin(α+β)/[cos(α+β)+ sin(α+β)] and k_AB = K cos(α+β)/[cos(α+β)+ sin(α+β)], where α is the compositional angle. Using this parameterization, the contour lines of ψ̅, where the conserved quantity is constant, have a fixed slope of -1 that is
independent of the compositional angle α.  However, note that the slope of the chemical steady state changes when varying α. Moreover, changes in α do not affect the activity in the system, in the sense that they do not alter the angle between the reaction nullcline and the tie lines. Thus, with this new parameterization, we can vary the compositional angle α for studying how different interactions affect active droplets while keeping the activity parameter β fixed. The rate K sets only a time scale.




§ THE SHARP-INTERFACE LIMIT



Similar to the droplet dynamics of passive phases (see Ref. <cit.> for example), we can study the dynamics of a single active droplet in the limit of a sharp interface. 
This approach for general multicomponent mixtures is discussed in detail in Ref. <cit.>. In this section, we derive the stationary state conditions of a chemically active droplet in the sharp-interface limit for a ternary mixture described by the free energy given in Eqs. (<ref>). With this choice of free energy, analytical results for the stationary radius and the critical value for the transition between intensive and extensive active droplets can be determined.

For a single droplet with radial symmetry sitting at the origin (r being the radial coordinate), we split the field into two domains (inside/outside, abbreviated in/out in the following) and couple these domains via boundary conditions at an interface at position r=R. 
Furthermore, we linearize the dynamics of the fields ϕ_A and ϕ_B in dyn_phi around the corresponding ϕ-values at the interface. 
In this approximation, the stationary state of a chemically active droplet in our ternary mixture is determined by the ϕ-values of A/B on both sides of the interface, named Φ_A^, Φ_A^, Φ_B^, and Φ_B^, and the interface position R. In the following, we explain how these five values are determined.


Locally, at the interface, we assume an equilibrium of phase separation, i.e., identical chemical potentials and osmotic pressure across the interface.
Furthermore, global conservation laws and vanishing flux differences across the interface fix a unique phase equilibrium locally at the interface and, therefore, a unique non-equilibrium steady state of the droplet, see Ref. <cit.>.

For our simple free-energy model (<ref>), we can obtain the four values of Φ_i^in/out generally as a function of R and α.
The determination of the stationary position of the interface R_stat  additionally needs flux relations (see below). Assuming local phase equilibrium at the interface constrains its concentrations to the two straight binodal lines shown in ph_d(b) for an infinitely large system. Laplace pressure effects of an interface with surface tension γ of finite-sized droplets are corrected up to linear order via the so-called Gibbs-Thomsom coefficients c <cit.>. We parameterize this phase equilibrium by the intersection of its tie line with the  ϕ_2-axis in its ϕ_1 and ϕ_2 representation (as defined in f). This intersection we call ϕ_2^inter, and write

    [ Φ^in/out_A; Φ^in/out_B ]   =
    [  cosα  sinα; -sinα  cosα ][ ± + γ c H; ϕ_2^inter ]

where we assumed that the conserved quantity is enriched inside and diluted outside. For such a spherical symmetric droplet, the mean curvature H is given by H = 2/R. If the dilute phase builds up the spherical droplet, the two labels in/out in interface_cond_app have to be swapped and H=-2/R.
Finally, we have to determine the intersection ϕ_2^inter via a global conservation law of the conserved quantity. Therefore, we note that while ϕ_A and ϕ_B have gradients in space, the conserved quantity is constant in space but jumps at the interface.
Thus, for a finite radial-symmetric system (system size R_sys) with an average amount of the conserved quantity ψ̅, we know

    1/2(Φ^_A+Φ^_B) R^3 + 1/2(Φ^_A+Φ^_B) (R_sys^3 - R^3) = ψ̅ R_sys^3

Using interface_cond_app, we find

    ϕ_2^inter = 
       2 R R_sys^3 ψ̅ + (2  R^4 - (2 c γ+  R) R_sys^3) (cos(α) + sin(α))/R R_sys^3 (cos(α) - sin(α)).

In an infinite system, however, the finite-sized droplet does not contribute to the average. Thus, the outside concentrations must fulfill
 
    ψ̅ = 1/2(Φ^_A+Φ^_B)

and  therefore, again by using interface_cond_app, we find

    ϕ_2^inter = 2 R ψ̅  - (2 c γ +  R) (cos(α) + sin(α))/R (cos(α) - sin(α))

We note by passing that it is only possible to derive the interface concentrations independently of the diffusivities and kinetic reaction rates due to our very symmetric form of the free energy density. Thus, these kinetic coefficients have to be taken into account only for the derivation of the stationary radius. In general, all five values have to be determined in parallel, a step that typically requires numerical solving schemes, see Ref. <cit.>.
    
An interface can only be stationary when the diffusive fluxes j_i^in/out are balanced across the interface, i.e. j_A^(R)= j_A^(R), and j_B^(R)= j_B^(R). 
These diffusive fluxes in the stationary state arise from the stationary concentration profiles ϕ_i^in/out(r) via j_i^in/out= - D_i^in/out∂_r ϕ_i^in/out(r), where D_i^in/out is the diffusion coefficient of component i in the corresponding domain obtained by linearizing. By construction, it is guaranteed that j_A^ = - j_B^ and j_A^ = - j_B^, such that once one of the flux balances is fulfilled, the second follows trivially.
Therefore, we derive only j_A^in/out in the following.

After linearizing, we can determine these fluxes through the analytical stationary ϕ_A-profiles in the two domains (solving the corresponding inhomogeneous Laplace equation in radial symmetry). From here, we can derive the total flux profiles, including the fluxes at the interface. We find

    j_A^(R)    =  (Φ_B^(R) k_AB^ - Φ_A^(R) k_BA^)
        
       ×( (λ^ R)/λ^ - 1/(λ^)^2 R)

where the rates k_ij^ are the linrearized reaction rates from dyn_phi and λ^ = √((D_A^ k_AB^ + D_B^ k_BA^)/(D_A^ D_B^)). Please note, that we use the dependency of the interface concentration Φ_i^(R) on the position of the interface R derived above.
However, in the outside domain (r>R), the solutions read


    j_A(R)^(R)    = (Φ_B^(R) k_AB^ - Φ_A^(R) k_BA^)
       ×λ^ (R - R_sys) cosh(λ^ (R - R_sys)) + ((λ^)^2 R R_sys-1) sinh(λ^ (R - R_sys))
        /(λ^)^2 R
        (λ^ R_syscosh(λ^ (R - R_sys)) + sinh(λ^ (R - R_sys)))


for a finite radial symmetric system with a no-flux boundary condition at the system size R_sys, or 

    j_A(R)^(R)     = -(Φ_B^(R) k_AB^ - Φ_A^(R) k_BA^) 
       × k_BA^ + D_A^/D_B^ k_AB^ +  (k_AB^ + k_BA^) λ^ R
        /(k_AB^ + k_BA^)(λ^)^2 R

for infinite large systems, where again k_ij^ are the linrearized reaction rates from dyn_phi and λ^ = √((D_A^ k_AB^ + D_B^ k_BA^)/(D_A^ D_B^)).
The stationary droplet radii can now be found by using the interface concentrations stated in interface_cond_app and Bintersec_finit_app or Bintersec_infinit_app, and numerically searching for the positions R at which jAin_app equals jAout_finite_app or jAout_infinite_app, depending on the system size.



§ THE CRITICAL VALUE Ψ_CRIT AND THE SCALING OF THE STATIONARY RADIUS IN ITS VICINITY



In App. <ref>, we explained how the the stationary droplet radii can be obtained in the sharp interface limit. Here, we argue first how we can obtain the value of the conserved quantity at the transition in an infinite system and, second, show how the stationary radius scales in the vicinity of this transition.

To find the critical value of the conserved quantity, we balance the fluxes across the interface for an infinitely large droplet,

    lim_R→∞ j_A^in(R) = lim_R→∞ j_A^out(R)

from jAin_app and jAout_infinite_app. When interface_cond_app and Bintersec_infinit_app are applied in these expressions, a dependency on the average conserved quantity ψ̅ follows. 
We find that there is only one value of this conserved quantity for which crit_trans_con_app is true. This value reads


    ψ̅_crit =
    -ϕ_0/2 (D_A^ k_AB^ + D_B^ k_BA^) (k_AB^ + k_BA^) + D_A^ D_B^λ^λ^ 
    ((k_AB^ - k_BA^) cos(2α) + (k_AB^ + k_BA^) sin(2α)
                )/sin(α)(
    (D_A^ k_AB^  + D_B^ k_BA^) (k_BA^ + k_AB^(α))
    + D_A^ D_B^λ^λ^  
    (k_BA^ 
    + k_AB^(α))
                )


When we apply our  parameterization of the reaction rates, i.e.,   
k_AB^in/out = K^in/outcos(α+β)/[cos(α+β)+ sin(α+β)] and k_BA^in/out= K^in/outsin(α+β)/[cos(α+β)+ sin(α+β)], see App. <ref>, we find 
ψ̅_crit as a function of the kinetic rates K^in/out, diffusivities D_i^in/out, and the angles α and β. With the assumptions D_A^ = D_A^, D_B^ = D_B^ and K^= K^, this expression simplifies and is solution of  psi_crit. 
Furthermore, we considered a droplet of the dense phase in a dilute environment, which is the case for ψ̅<0 (- branch). For ψ̅>0, the interface concentrations  swapped resulting in + branch of the solution of psi_crit.

Furthermore, we can check the scaling of the stationary radius of intensive active droplets in the vicinity of the critical transition. For this, we expand j_A^ (R)= j_A^ (R) for large R. We find the stationary radius R_stat as a function of ψ̅ that scales like R_stat∝ |ψ̅ - ψ̅_crit|^-1. The general solution can be obtained from the equations above straightforwardly. Due to its length, however, in this appendix we restrict ourselves to the special case of D_A^ = D_A^ = D_B^ = D_B^ = D and K^ = K^ =K. We find for the case of a droplet of the dense phase in a dilute environment (ψ̅<ψ̅_crit<0)

    R_stat = -  Λϕ_0 tan(β) cos(α) - sin(α)
        /ψ̅- ψ̅_crit

where we defined the length scale Λ = √(D/K) and used the fact that for large droplets, the Laplace pressure becomes negligible. 




§ PARAMETER CHOICES AND METHODS USED FOR FIGURES


To numerically solve the dynamics equations in this work, we rescale time t · K with K denoting the  overall rate (overall_reaction_rate), 
and position x/ℓ where ℓ=√(ϕ_0^2 κ_A/b_1), thus, half of the interface width in the continuous model. 
For simplicity, we consider the overall rate equal inside and outside for all studies (K^in=K^out=K). 
These rescalings yield the following non-dimensional parameters in our numerical studies:
D_i^in/out /(ℓ^2K), γ c/ ℓ, and Γ/(ℓ^2K b_1 ). 




 §.§ Figures showing results obtained in the sharp interface model


In the previous sections, we used the sharp interface limit to calculate the results shown in 
Fig. <ref>(a,b).
Using a standard root-finding scheme, we can numerically solve the specific values of the stationary interfaces shown in Fig. <ref> in the main text. The parameter values used to produce these figures are D_A^/(ℓ^2K)=0.0014^-1, D_B^/(ℓ^2K)=0.0014^-1, α = 0, β = 0.12 π, γ c/ℓ = 1/6,  and ϕ_0 = 1.
 


 §.§ Figures showing results obtained in the continuous model

All the results obtained in the continuous model were obtained by using an implicit-explicit Runge-Kutta solver of the second-third order. We considered periodic boundary conditions and approximated the higher-order derivatives using pseudo-spectral methods on a regular lattice.
For simplicity, we used b_1/ϕ_0^2=b_2=ϕ_0=1, α = 0, κ_A = κ_B and β=0.12π were not otherwise explicitly mentioned in the figure caption.
The mobilities were set to Γ/(ℓ^2K b_1)=0.0014^-1 in Fig. <ref>, Γ/(ℓ^2K b_1 )=0.0028^-1 in Fig. <ref>, Γ/(ℓ^2K b_1 )=0.00112^-1 in Fig. <ref>, and  Γ/(ℓ^2K b_1 ϕ_0^2)=0.004^-1 in Fig. <ref>.
Furthermore, we used N_1 × N_2 or N_1 × N_2 × N_3 grid points for two-dimensional or three-dimensional systems of sizes of L_1 × L_2 or L_1 × L_2  × L_3, respectively.
The specific values for the corresponding figures are listed in Table <ref>.


