Abstract
Liquid–liquid phase separation leads to the formation of condensed phases that coexist with a fluid. Here we investigate how the positions of a condensed phase can be controlled by using concentration gradients of a regulator that influences phase separation. We consider a mean field model of a ternary mixture where a concentration gradient of a regulator is imposed by an external potential. A novel first order phase transition occurs at which the position of the condensed phase switches in a discontinuous manner. This mechanism could have implications for the spatial organisation of biological cells and provides a control mechanism for droplets in microfluidic systems.
Export citation and abstract BibTeX RIS
Original content from this work may be used under the terms of the Creative Commons Attribution 3.0 licence. Any further distribution of this work must maintain attribution to the author(s) and the title of the work, journal citation and DOI.
1. Introduction: positioning of condensed phases
Phase separation of a mixture refers to the formation of a liquid condensed phase that coexists with a dilute phase of lower concentration [1, 2]. Such demixing is the result of a first order thermodynamic phase transition where the concentration difference between the phases changes discontinuously. It can be observed in many forms in everyday life, for example when oil is added to water. The occurrence of a transition from the homogeneous mixture to a system with coexisting phases can be controlled by temperature or by changing the composition of the mixture. Condensed phases are influenced by surfaces possibly causing wetting transitions [3–5]. Furthermore, phase separation can be affected by external forces such as gravity causing sedimentation.
A key question is how liquid condensed phases such as droplets are positioned in systems with external cues like concentration gradients or external fields. The study of positioning of phases provides general insights in the physics of phase separation of spatially inhomogeneous systems. Understanding the underlying mechanism of the positioning may open the possibility of applications in microfuidic devices. Positioned liquid condensed phases could be used to seal and open junctions at specific locations in the microfluidic device, or simply position chemicals that partition into the condensed phase. The positioning of liquid condensed phases in a complex mixture also plays a role in cell biology. In particular, positioned droplets are used to segregate molecules during asymmetric cell division [6–9].
Here we study the equilibrium physics of the positioning of two liquid condensed phases in inhomogeneous systems. We present a simplified model that provides the basic mechanism for the positioning at thermal equilibrium which can be further extended to non-equilibrium processes such as the kinetics of droplet formation and ripening [10]. In our model, phase separation of two components is subject to a concentration gradient of a regulator component where the gradient is generated by an external field. The regulator component affects demixing of the two components but does not phase separate itself. The system then relaxes to a spatially inhomogeneous thermodynamic equilibrium state with two coexisting phases positioned by the regulator gradient. The spatial distributions of the three concentration profiles at thermal equilibrium are determined by minimising a mean field free energy functional. We find that as a function of an interaction parameter the position of the condensed phase switches discontinuously from a position in the region of large regulator concentration (correlated state) to the region of low regulator concentration (anti-correlated); see figures 1(a), (b). This switching of position corresponds to a novel, equilibrium first order phase transition at which an order parameter jumps discontinuously.
2. Equilibrium model for spatial regulation of phase separation
In our equilibrium model for spatial regulation of phase separation we consider three components [11]: two components which can demix from each other, A and B, and a regulator R that interacts with these components. The regulator affects phase separation but does not demix from A and B. Demixing and interactions with the regulator are described by the Flory–Huggins free energy density for three components ([12, 13] and appendix A for a derivation):
where ϕi denotes the volume fraction field of component i = A, B, R. For simplicity, we consider an incompressible system where the molecular volumes are equal to ν and . The logarithmic contributions in equation (1) correspond to the mixing entropy, while the quadratic terms describe the molecular interactions between the components; χij is the interaction parameter between component i and j. The gradient terms represent contributions to the free energy density associated with spatial inhomogeneities. They introduce two length scales, and . The regulator R is subject to an external field described by a position-dependent potential U(x). In the following we consider a potential that varies solely along the x-coordinate, , where s > 0 characterises the slope of the potential and its inverse corresponds to a third length scale in our model. In the absence of A molecules and for diluted regulator, , ϕR(x) attains a concentration profile that is linear in space with a slope proportional to s. For simplicity, we also consider a one dimensional system of size L and two type of boundary conditions: (i) Neumann boundary conditions, and (), where the primes denote spatial derivatives, and (ii) periodic boundaries with ϕi (0) = ϕi (L) and . The conditions (i) imply that there is no explicit energetic bias to wet or dewet the boundary. However, such a boundary mediates a coupling between the slopes of the volume fraction fields for A and R at the boundary. In contrast, the periodic conditions (ii) allow to study the system in the absence of boundaries. For the considered case of an external potential U(x) varying only along the x-coordinate, the restriction to a one dimensional, phase separating system represents a valid approximation for large system sizes, where the interface between the condensed phases becomes flat.
3. Equilibrium concentration profiles
To calculate the equilibrium profiles and , we minimise the free energy
Due to particle number conservation, two constraints are imposed for the minimisation: each field (i = A, R) obeys , where are the average volume fractions and . Variation of the free energy equation (2) with the constraints of particle number conservation implies (i = A, R):
where λR and λA are Lagrange multipliers, and the primes are derivatives with respect to the position x. The boundary terms vanish for both, Neumann and periodic boundary conditions. Using the explicit form of the free energy density (equation (1)), the Euler–Lagrange equations can be derived:
where χ = χAR − χAB − χBR. We solve these equations using a finite difference solver (e.g., bvp4c [14]). As control parameters we consider the three interaction parameters χAR, χAB and χBR, the slope of the external potential s and the mean volume fraction of A-material, . The mean regulator material is fixed to a rather dilute value of in all presented studies. Moreover, we focus on the limit of strong phase segregation where the interfacial width is small compared to the system size, i.e. . In this limit, we verified that our results depend only weakly on the specific values of κi.
4. Discontinuous phase transition in the positions of the condensed phases
Solving the Euler–Lagrange equations (4) and (5) with Neumann boundary conditions (i), we find two spatially inhomogeneous solutions for component A, which we denote and , and the two corresponding solutions for the regulator component R, are denoted and (the profile of B follows from volume conservation). The phase separating material A is either accumulated close to the right boundary of the system ( and ) and correlated with the concentration of the regulator material (figure 1(a)), or it is accumulated at the left ( and ) and anti-correlated with the regulator (figure 1(b)). Upon varying the interaction parameter χBR in figures 2(a), (b), the free energies of the correlated and the anti-correlated states, and , are different. They intersect at a certain value of the interaction parameter, χBR = , which defines the transition point of the system (figure 2(a)). At this point the lowest free energy exhibits a kink, which means that the system undergoes a discontinuous phase transition when switching from the spatially anti-correlated ('left') to the spatially correlated ('right') solution with respect to the regulator. This transition point does not depend on the slope of the regulator s, while increases linearly with s (figures 3(a), (b)). Similar results can be found when fixing the value of χBR and considering χAR as control parameter.
Download figure:
Standard image High-resolution imageDownload figure:
Standard image High-resolution imageThe emergence of the correlated or the anti-correlated state can be qualitatively understood when considering the energetic interactions with the regulator. In the case of a negative value of the interaction parameter χBR, the B-particles prefer the neighbourhood of the regulator R. This preference leads to an accumulation of A-particles at low regulator concentration and thus to an anti-correlated equilibrium state. If, however, the value of the interaction parameter χBR has a positive value, the B-particles have a tendency to avoid the vicinity of the R-particles causing an accumulation of A-particles at high regulator concentration. As a result, the equilibrium profile of A is correlated with the regulator.
The transition between the correlated and anti-correlated state can be described by the following set of order parameters:
where the squared normalisation with , denoting the variance and , where is the Heaviside step function. This normalisation ensures that −1 < ρij < 1 and ρij = ±1 if . The derivative of the free energy with respect to the interaction parameter χij generates the covariance between the spatially dependent fields ϕi(x) and ϕj(x). If the fields are spatially correlated, ρij > 0, and if they are anti-correlated, ρij < 0. For homogeneous fields with , ρij = 0. Varying the interaction parameter χBR (figure 2(b)), the order parameters ρBR and ρAR jump at the transition point , while in the absence of a regulator gradient (s = 0), they change smoothly (figure 2(b), grey lines). The jump of both order parameters in the presence of a regulator gradient indicates that the spatial correlation of A and B with respect to R changes abruptly, which is expected in the case of a first order phase transition.
By means of the order parameter ρBR (equation (6)) we can now discuss the phase diagrams as a function of the interaction parameters for different volume fractions of the demixing material, . We find three regions (figures 4(a)–(c)): a mixed region (M), where volume fraction profiles are only weakly inhomogeneous and no phase separation occurs. In addition, there are two regions, (C) and (AC), where components A and B phase separate and A is spatially correlated or anti-correlated with the regulator R, respectively. There exists a triple point where all three states have the same free energy. For , the shape of the transition line between correlated and anti-correlated states is straight and the transition point is independent of χAB (figure 4(b)). If is decreased, the correlated state is favoured while for increasing , the anti-correlated state is preferred. The transition line to the mixed states is horizontal for (figure 4(b)). For both, larger and smaller -values, it becomes curved and moves towards larger χAB interaction parameters. This behaviour can be qualitatively understood by the upshift of the demixing threshold χAB once deviates from 1/2, as known for binary systems. Since the regulator R is considered to be dilute, this analogy to binary systems is a good approximation ( in equation (1)). Both trends explain the parabolic shape of the positions of the triple point in the phase diagrams when the mean volume fraction of the demixing material is varied (figure 4(d)).
Download figure:
Standard image High-resolution imageThe transition line in the phase diagrams between the correlated and anti-correlated solution as a function of the interaction parameters can be estimated analytically. In the absence of a regulator gradient (s = 0), the free energies of both solutions are the same for all interaction parameters for which phase separation occurs. In the presence of a regulator gradient, however, the free energies corresponding to the correlated and the anti-correlated solutions are unequal for most points in the phase diagram. The reason is that the external potential U(x) forces the regulator to form a gradient, and thus the interactions with the regulator lead to different free energies of the correlated and anti-correlated states. Only along the transition line between both states, the free energies are equal:
This condition can be used to estimate the transition line for varying interaction parameters and the slope of the external potential s. To estimate ΔF we parametrise the profiles of the stationary solutions and using physical assumptions that are in agreement with our numerical results. First, we idealise the already narrow interface of the demixed component as sharp, which can be realised by a strong phase separation far away from the critical point. Second, we use only one profile, denoted as , for both regulator states because the regulator is dilute and maintained by the external potential. By means of the numerical solution, we actually confirm that close to the transition line. In addition, we approximate the regulator profile as linear function of slope m, neglecting spatial nonlinearities that can be seen in figures 1(a), (b). In appendix B we show that for small enough κR, the amount of regulator inside this peak become negligible. Finally, the low volume fractions outside the condensed phase of the demixed binary A-B system are approximated as constant values . The value of is determined by the Flory–Huggins parameter χAB and can be calculated from the binary phase diagram by solving for . The larger volume fraction (inside) shows a weakly linear profile (figures 1(a), (b)). For diluted regulator, the volume fraction inside the condensed phase can be well described as , where is the constant volume fraction inside the condensed phase of the binary A-B mixture (figure 5). In summary, in the case of a diluted regulator, a linear regulator profile and strong phase separation the approximated profiles are:
The conservation of A determines the domain sizes of the phase separated region,
which depends on the slope of the regulator m. For the special case of zero slope, the domain sizes left and right are equal, . To calculate ΔF (equation (7)), the free energy density (equation (1)) is integrated in the domain [0, L]. Using the approximated profiles (equations (8)–(10)), we find
where the function depends only on the parameters of the simplified solutions (see equations (8)–(10)), and reads
In the limit , the function as the domain sizes of the phase separated regions become equal, l(0) = r(0), for vanishing regulator slope m. To leading order in the regulator slope m,
The expression above indicates that for small regulator slopes m, the asymmetry of the domain sizes of the phase separated region l and r is not essential for the free energy difference ΔF. Consistently, according to equations (12) and (14), the free energy difference between the correlated and anti-correlated state vanishes (ΔF = 0), if there is no regulator gradient (m = 0). In the presence of a regulator gradient (), the free energy difference ΔF is zero only if χBR = χAR, which corresponds to the transition line between the correlated and anti-correlated state according to the approximate profiles (equations (8)–(10)). This prediction is in very good agreement with our numerical results for see figures 4(b) and 6(a).
Download figure:
Standard image High-resolution imageThe condition for the transition line, χAR = χBR, for the case (see figure 4(b)) can also be understood by symmetry arguments. For and dilute regulator, switching the identity of A and B leads to the same free energy density. Thus, in the presence of an external potential acting on the regulator, the difference in free energy between the correlated and anti-correlated state ΔF vanishes at equal interaction parameters with respect to the regulator, χAR = χBR.
By means of the approximated profiles (8)–(10) and the definition of the order parameter (6), we can estimate how the jump of the order parameter (definition see figure 2(b)) at the transition point depends on the model parameters:
We find that the estimated as a function of the slope of the regulator and the interaction parameter χAB almost perfectly describes the data obtained from the numerical minimisation of the free energy (figure 6(b)). This agreement shows that the proposed stationary profiles (equations (8)–(10)) are a consistent approximation to describe the discontinuous phase transition in the case of strong phase separation and a linear and diluted regulator profile. We could also show that the asymmetry of the phase separated domains, l and r, is not essential for the jump of the order parameter (equation (14)). Instead the jump is determined by the slope of the regulator profile where the jump height is affected by the mean amount of regulator material and the degree of phase separation characterised by χAB (figure 6(b)).
Download figure:
Standard image High-resolution image5. Discontinuous phase transition in a periodic domain and the presence of fluctuations
The phase diagrams (figures 4(a)–(c)) depend on the boundary conditions raising the question whether the boundary may play a role for the existence of the phase transition. To this end we considered a periodic system without boundaries. As for the non-periodic system, we are minimising the free energy (equation (1)), now using periodic boundaries with ϕi (0) = ϕi (L) and . In the periodic domain, we also use a periodic external potential:
where ω denotes a phase shift. The value of the phase is chosen such that the region of segregated A-material is placed at x = 0. The logarithmic form of the potential ensures that a sinus distribution of the regulator is obtained in the dilute limit.
We find the same main results as for the non-periodic system with Neumann boundary conditions, namely the existence of a discontinuous phase transition. In particular, we find two stationary solutions of different spatial correlations with respect to the regulator. They switch at by a discontinuous phase transition (figures 7(a)–(c)). Therefore, a boundary of the system is not a necessary requirement for the emergence of the reported discontinuous phase transition since it also exists in the absence of boundaries. Thus the transition is not induced by boundaries as for example in the case of wetting transitions [3–5].
Download figure:
Standard image High-resolution imageWe have also scrutinised the robustness of the phase transition in position considering Monte Carlo studies corresponding to the mean field free energy density equation (1). We could confirm that the positioning mechanisms is robust against the fluctuations arising from the probabilistic Monte Carlo update and that the phase diagrams of correlated and anti-correlated states coincide qualitatively.
6. Experimental verification and outlook
The discontinuous switching of phase separation could be tested experimentally. We suggest to use a soluble salt of high magnetic susceptibility in order to create and maintain the regulator concentration gradient by applying an inhomogeneous magnetic field [15]. Phase separation in a regulator gradient could then be observed by introducing components that phase separate in a salt dependent manner. In particular, a pre-formed droplet could be added to an existing regulator gradient or the regulator gradient is created after coarsening via Ostwald-ripening and coalescence is completed [10, 16–18]. The phase transition could be triggered by changing the concentrations of the phase separating material, by changing the temperature or by adding additional components that influence the interaction parameters. The validity of our theory could be probed by comparing the order parameter with corresponding experimental measurements. In particular, the order parameter jump at the transition point could be determined for different amounts of regulator material (see figure 6).
Our finding of a phase transition in the position of coexisting phases could also be relevant for applications. As the composition of condensed phases provide a distinct chemical environment (e.g. for chemical reactions), our work suggests a novel mechanism to control and switch chemical environments in microfluidic devices by the use of a phase transition in position.
Acknowledgments
We would like to thank Martin Elstner and Omar Adame Arana for fruitful and stimulating discussions. This project was supported by the Center for Advancing Electronics Dresden (cfAED). Christoph A Weber thanks the German Research Foundation (DFG) for financial support. Samuel Krüger and Christoph A Weber contributed equally to this work.
Appendix A.: Gradient contributions in the ternary Flory–Huggins free energy density
A.1. Derivation using a mean field approximation
To derive the gradient contribution in the ternary Flory–Huggins free energy density, we start from the local mean-field free energy on the lattice and calculate the continuum limit of this free energy as shown in [19] for a binary system. The local free energy density of the three component system is derived in [20, 21] using a mean-field approximation:
where ν is the molecular volume. The greek indices α and β indicate the positions on the lattice. The first two lines describe the entropy of the mixture. The remaining contributions stem from the interactions between the components at neighbouring lattice sites.
In the next steps we will perform the continuum limit. In case of the entropic contribution, we can simply replace ϕi(α) → ϕi(x). In order to perform the continuum limit for the energetic contributions, we rearrange the corresponding terms as follows:
Each contribution can be rewritten as
We can identify the Flory–Huggins interaction parameter as . In the continuum limit we can introduce the gradient of the volume fractions as , where a denotes the lattice size. We finally obtain the free energy with the free energy density given as
where
The parameters characterising the 'penalty' corresponding to spatial inhomogeneities are κi = a2χi B, i∈{A, R}, and .
A.2. Phenomenological derivation
In the Ginzburg–Landau free energy the penalties corresponding to spatial inhomogeneities are phenomenologically introduced based on symmetry considerations:
where since spatial inhomogeneities are unfavored. Moreover, f0 is the free energy density that only depends on the volume fractions ϕi, i ∈ A, B, R. However, only two volume fraction fields are independent due to particle conservation and incompressibility, 1 = ϕA + ϕB + ϕR. Thus we can write ∇ϕB = −∇ϕA − ∇ϕR, leading to
Here, , and .
A.3. Choice of the parameters κi
In the presented studies, we have chosen κA = κ for simplicity. Please note that the derivation presented in appendix A.1 is based on a mean field approximation and therefore it should only serve as an estimate for the values κi. We chose the values for the parameters κA and κR that are consistent with these estimates (see figure captions).
Appendix B.: Regulator peak at the interface
The numerically obtained regulator profiles show a significant peak at the interface between the A-rich and the B-rich phase (figures 1(a), (b)). The emergence of the regulator peak can be understood by entropic and energetic considerations of the free energy. For large and positive χAR and χBR (corresponding to a repulsive tendency with respect to the regulator), the energy of the system decreases as regulator accumulates at the interface. Moreover, the entropy decreases as the composition of the interfacial region of all three components is closer to a well-mixed state.
The amount of regulator material that is accumulated at the interface is strongly influenced by the κi-parameters; see figure B1. In figure B2(a), the peak area is shown for varying κi-parameters. For simplicity, we chose κA = κR = κ. The peak area vanishes as the κi-parameters approach zero. This behaviour is expected since these parameters determine the size of the interface between the phase separated phases. In this limit, the approximated profiles (equations (8) and (10)) accurately describe the numerical solutions and thus the corresponding predictions for the phase boundaries coincide well with the phase boundaries obtained from the numerical calculations.
Download figure:
Standard image High-resolution imageDownload figure:
Standard image High-resolution imageHowever, the peak height and thereby the existence of the peak is approximately independent of κi (figure B2(b)). This indicates that the existence of the peak may depend on the interaction parameters for example. Since we also observed that the peak is more pronounced at the transition line between anti-correlated state and correlated state, we investigated the energetic influence on the peak height along the transition line. As derived in section 4, the transition line is governed by the condition χAR = χBR for . We find that the peak height increases as a function of the energetic parameters χAR = χBR (figure B2(c)). Large and positive values of χAR and χBR correspond to a repulsive tendency with respect to the regulator. This indicates that the energetic contribution to the free energy decreases as regulator accumulates at the interface.