Confinement-Induced Fractionation and Liquid–Liquid Phase Separation of Polymer Mixtures

The formation of (bio)molecular condensates via liquid–liquid phase separation in cells has received increasing attention, as these aggregates play important functional and regulatory roles within biological systems. However, the majority of studies focused on the behavior of pure systems in bulk solutions, thus neglecting confinement effects and the interplay between the numerous molecules present in cells. To better understand the physical mechanisms driving condensation in cellular environments, we perform molecular simulations of binary polymer mixtures in spherical droplets, considering both monodisperse and polydisperse molecular weight distributions for the longer polymer species. We find that confinement induces a spatial separation of the polymers by length, with the longer ones moving to the droplet center. This partitioning causes a distinct increase in the local polymer concentration near the droplet center, which is more pronounced in polydisperse systems. Consequently, the confined systems exhibit liquid–liquid phase separation at average polymer concentrations where bulk systems are still in the one-phase regime.


Introduction
The liquid-liquid phase separation (LLPS) of (bio)molecules in living cells has attracted much attention as a mechanism for intracellular organization via the formation of biomolecular condensates [1,2]. To elucidate the underlying mechanisms of LLPS, bulk solutions of purified biomolecules from cells have been analyzed extensively over the past decade. Although such in vitro studies facilitate the analysis of LLPS and comparisons with theoretical models [3], they systematically ignore confinement effects and the molecular diversity encountered in cellular environments [4]. For example, broad molecular weight distributions can have a profound impact on the phase behavior of polymers, e.g., leading to the self-assembly of monodisperse micelles from polydisperse surfactants [5] or the fractionation of polymer chains by molar mass in solutions near criticality [6,7]. Typically, the demixing of polymers becomes more pronounced in confinemenet due to the associated spatial inhomogeneity [8]. For example, simulations of binary polymer mixtures in spherical droplets have revealed an entropy-driven spatial segregation of the confined polymers based on their molecular weight, stiffness, and/or topology [9][10][11]. It is thus important to understand the phase behavior of confined systems for designing experiments and understanding the physical mechanisms driving LLPS in cellular environments.
Experimentally, confinement effects have been studied using DNA-based protocells [12,13] or water-in-oil droplets containing polymer mixtures such as polyethylene glycol (PEG) and bovine serum albumin, PEG and DNA, or PEG and Dextran, serving as a synthetic cytoplasm [14][15][16][17]. At sufficiently high concentrations, these polymer mixtures separate into two coexisting aqueous compartments, which can be used to localize additional components such as proteins [14,15]. Recently, Watanabe et al. systematically studied the phase coexistence of PEG-Dextran mixtures in cell-sized droplets, discovering that the two-phase coexistence region in small droplets extends to much lower PEG and Dextran concentrations compared to bulk systems [18]. The bulk behavior was only recovered for rather large droplets with radii R > 20 µm. They speculated that this R-dependent phase separation stemmed from a confinement-induced partitioning of the polymers. To elucidate the origin of this confinement-induced LLPS, we simulate PEG-Dextran mixtures in spherical confinement at two different droplet radii.

Simulation Model
We perform dissipative particle dynamics (DPD) simulations [19][20][21] using a coarsegrained polymer model in an explicit solvent. The simulations contain three different particle types, i.e., P (PEG), D (Dextran), and W (water), which are equal in their diameter a and mass m. To approach experimental length-and time-scales, we model each PEG chain as a single spherical bead of diameter a, as shown schematically in Figure 1. To facilitate a direct comparison with experiments, we have chosen polymer sizes similar to previously conducted experimental studies [18,22]. To establish a mapping between experiments and simulations, we first estimated the characteristic size of a PEG chain in solution: the mass and length of a Kuhn segment of PEG are m P K = 0.137 kg/mol and b P K = 1.1 nm, respectively, resulting in N K ≡ M w /m K ≈ 44 Kuhn segments for a PEG chain with molecular weight M P w = 6 kg/mol. At θ-conditions, the root-mean-square radius of gyration can be estimated as R g = b K (N K /6) 1/2 , resulting in a value of R P g ≈ 3.0 nm. The DPD particles representing water have the same diameter a, and thus contain roughly 3700 water molecules. The Dextran chains are modeled as linear chains consisting of DPD particles with diameter a because the branches of Dextran are, on average, typically shorter than three glucose units [23] and therefore cannot be resolved at this level of coarse-graining (see Figure 1). We determined the number of DPD beads per Dextran chain, N D , by matching R D g from single-chain simulations to experimental R g measurements; in ref. [24], R g was derived from self-diffusion coefficient measurements of Dextran chains in water at T = 293 K, with molecular weights ranging between 4 kg/mol and 464 kg/mol. By extrapolating their data, we estimated R g = 23.8 nm for a Dextran chain with 500 kg/mol, which leads to N D = 80 for our simulation model. Hence, roughly 40 Dextran monomers are mapped to one DPD bead. Thus, each DPD particle effectively represents a coil-like polymer segment, which typically interact with each other via soft and bounded pair potentials [25][26][27]. Therefore, we use the standard soft repulsion for the conservative forces acting between bonded and non-bonded DPD beads where r is the distance between the two particles, andr is the unit vector connecting the two. The parameter A ij controls the repulsion strength between particles of type i and j and has been set according to [21] with A ii = 25 k B T/a and Flory-Huggins interaction parameter χ ij (the values for χ ij are discussed in Section 3.1 below). This specific value for A ii was originally determined by Groot and Warren to match the compressibility of water [21], but since then, this choice has been widely used to model other (incompressible) liquids and polymer melts [28][29][30][31].
Neighboring monomers within a Dextran chain are bonded through harmonic springs with force Note that the individual beads should not be regarded as single monomers, but rather as fluid elements containing several chain segments. Thus, the harmonic springs between DPD particles do not represent covalent bonds between monomers, but instead ensure the connectivity of the Dextran chains. We use a soft spring constant k = 4 k B T/a 2 [28], which is consistent with the typical free energy associated with deforming a polymer in the blob model [32].
In addition to these two conservative forces, all particles are subjected to pairwise dissipative and random forces with drag coefficient γ ij , and velocity difference between two particles ∆v. The parameter ξ is a uniformly distributed random number drawn for each particle pair, with zero mean ξ(t) = 0 and variance ξ(t)ξ(t ) = 2k B Tδ(t − t ) to satisfy the fluctuation-dissipation theorem. For simplicity, we used the same drag coefficient γ ij = γ = 4.5 √ mk B T/a for all particles, and the standard DPD weight function [21] ω(r) = (1 − r/a) 2 , r ≤ a 0, r > a.
For the droplet simulations, all beads are confined to a spherical container with radius R by applying a purely repulsive Weeks-Chandler-Andersen (WCA) potential [33] U WCA (r ) = 4k B T a r 12 − a r 6 + 1 4 , r ≤ 2 1/6 a 0, where r is the distance between the droplet surface and the center of a bead. In all simulations, the particle number density was set to ρ = 3 a −3 . The equations of motion were integrated using a time step of ∆t = 0.02 τ, with the unit of time being τ. Each simulation was run for at least 10 7 time steps, and three independent simulations were performed for each parameter set to improve the statistics and determine measurement uncertainties.

Parameterization and Bulk Phase Behavior
To faithfully reproduce experimental conditions, we first needed to determine the interaction parameters for the PEG (P), Dextran (D), and water (W) particles (see Section 2). Following previous simulation studies [34,35], we used a Flory-Huggins interaction parameter of χ P−W = 0.3 for the PEG-water interactions, which was extracted from experimental phase coexistence measurements by Saeki et al. [36]. We took χ D−W = 0.50 for the Dextranwater interactions, derived from experimental vapor-pressure measurements at T = 298 K conducted by Bercea et al. [37]. This value is in excellent agreement with previous findings by Clark [38], who applied a Flory-Huggins theory-based analysis to experimental tie line data of PEG-Dextran mixtures in water at T = 298 − 300 K [39,40]; Clark also extracted the PEG-Dextran interaction parameter from his analysis, i.e., χ P−D = 0.031 ± 0.007. Note that the physically relevant quantity for a pair of polymers is the combined Flory-Huggins parameter χN [32], and thus a scaled interaction parameter χ eff P−D ≈ 1.23 needs to be used in the simulations to reach the same (χN) P−D ≈ 100 as in the experiments, since we mapped about 40 Dextran monomers onto one DPD particle. Given the uncertainties in extracting χ ij from experiments, and the high degree of coarse-graining of our model, it was, however, unclear whether this initial parameterization would faithfully reproduce the interactions between the PEG and Dextran chains in water.
To test (and tune, if necessary) the interaction parameters of our coarse-grained model, we first attempted to reproduce the experimentally known phase behavior [18,41] of aqueous PEG-Dextran mixtures in the bulk. The groups of Dimova [41] and Yanagisawa [18] determined the binodals of mixtures containing short PEG chains (M P w = 6 or 8 kg/mol) and long Dextran chains (M D w ≈ 500 kg/mol), finding a critical concentration of slightly below 4 wt% for PEG and 4 wt% for Dextran. Initially, we simulated bulk systems at two concentrations, i.e., c P = 3 wt% for PEG and c D = 3 wt% for Dextran, where we expected a single mixed phase, and at c P = 4 wt% and c D = 4 wt%, where the polymers should phase separate. We performed simulations in a cubic box with an edge length of L ≈ 360 nm and applied periodic boundary conditions to all three Cartesian directions. The systems were initialized by placing all PEG particles and Dextran chains in opposite halves of the simulation box and were then run until the density profiles did not change anymore. For χ eff P−D = 1.23, we found a single phase at both concentrations, which indicates that the initially chosen χ eff P−D value was too small. Therefore, we iteratively increased χ eff P−D until we observed aggregation of the Dextran chains at the higher concentration and observed a fully mixed system at the lower concentration. This was achieved for χ eff P−D = 18.6 ( Figure 2), which is about 15 times larger than our initial estimate for the PEG-Dextran interactions. To quantify the size and concentration of the Dextran condensates, we performed a cluster analysis using the density-based spatial clustering of applications with noise (DBSCAN) algorithm [43]; here, Dextran monomers are assigned to the same aggregate if their distance is smaller than 7 nm, which roughly corresponds to the position of the first minimum of the radial distribution function g(r) between Dextran particles and all other particle types. To establish a baseline, we performed additional simulations of ideal mixtures by setting χ ij = 0. Figure 2c shows the radial concentration profile c D (r) of Dextran monomers belonging to the largest cluster identified in the system. For the system at c P = 4 wt% and c D = 4 wt%, this analysis revealed one large droplet with a Dextran concentration of about 35 wt% in its core. In contrast, the cluster identified at c P = 3 wt% and c D = 3 wt% was much smaller and more diluted, and the resulting concentration profile resembled that of an ideal Dextran chain in solution. To characterize the size of the aggregates in more detail, we also computed the average number of Dextran monomers in a cluster, M . For the the ideal mixtures, this analysis yielded M 0 ≈ 100, which is comparable to the number of monomers per Dextran chain (N D = 80). For the less concentrated system with χ ij = 0, we found M ≈ 400, which indicates the formation of small intermittent clusters. In contrast, at c P = 4 wt% and c D = 4 wt%, we found a much larger value of M ≈ 4800.
The discrepancy between the initial estimate and the final value of χ eff P−D is rather large, and we can only speculate about its origin: Clark extracted χ P−D from experimental coexistence curves using Flory-Huggins solution theory [38], which ignores the polymer architecture and thus the branching of the Dextran chains. Further, a monodisperse molecular weight distribution is assumed in his treatment, although generally available Dextran polymers typically have a broad molecular weight distribution [7,17,18,22,44]. Finally, we used a rather coarse-grained description, which maps about 40 monomers onto a single bead, resulting in a larger entropy of mixing in the simulations compared to the experiments [32]. To compensate for these effects, we will use χ eff P−D = 18.6 in the following. Figure 3 shows the resulting phase coexistence results from our simulations compared to experiments [18], which are in excellent agreement with each other, thus further corroborating the appropriate parameterization of our simulation model.  [18]. The dotted lines indicate the estimated binodal. Simulation data are generated from monodisperse systems with M P n = 6 kg/mol and M D n = 500 kg/mol, while experimental data were gathered for blends with M P w = 6 kg/mol and M D w = 500 kg/mol.

Phase Behavior in Confinement
In recent experiments by Watanabe et al., the confinement-induced phase separation of PEG-Dextran mixtures was observed for droplets with radii R < 20 µm [18]. Simulating such large droplets is computationally infeasible, even at the employed level of coarse-graining (see Methods section), as roughly 7 × 10 9 particles would already be needed to represent a droplet with R = 5 µm. Therefore, we performed simulations at two smaller radii, i.e., R ≈ 260 nm and R ≈ 380 nm, which should still allow us to capture the effect of confinement on the phase behavior. Further, we considered mixtures containing either monodisperse or polydisperse Dextran chains, since high-molecular-weight Dextran usually has a broad molecular weight distribution; for example, the Dimova group used Dextran chains with M w ≈ 380 − 490 kg/mol and dispersities in the range of Ð ≡ M w /M n ≈ 1.8 − 2.2 [7,41,44], while the Yanagisawa group used Dextran with M w ≈ 500 kg/mol and Ð ≈ 3.1 [17,18]. In our simulations with polydisperse Dextran chains, we drew the molecular weight of each polymer from a Gaussian distribution, targeting Ð ≈ 1.5 and M n ≈ 500 kg/mol. The PEG chains were kept monodisperse throughout, which is consistent with the rather small polydispersity of Ð ≈ 1.1 reported in the experimental literature for low-molecular-weight PEG [7,17,44]. In all simulations, we selected the number of PEG and Dextran chains so that c P = 3 wt% and c D = 3 wt%, averaged over the entire droplet volume, which lies in the mixed one-phase regime of the bulk phase diagram (see Section 3.1).
To study the spatial distribution of the PEG and Dextran polymers in the droplet, we first calculated the radial monomer concentration profiles c(r) of the two species (Figure 4), which reveal several important features: (i) in all cases, there is a local surplus and layering of PEG near the droplet surface, which is typical for short molecules close to hard walls [45]. In contrast, the long Dextran polymers are depleted from the droplet surface because of the associated loss in conformational entropy in that region [11]. In the monodisperse case, the width of this depletion zone is roughly 2 R D g ≈ 50 nm, with R D g ≈ 23.8 nm being the radius of gyration of a Dextran chain at infinite dilution. In contrast, the excluded region is much wider for the polydisperse case due to the broader R D g spectrum; (ii) as a result, the Dextran concentration in the droplet center becomes distinctly larger than the average value (3 wt%), reaching almost 8 wt% for the polydisperse case in the small droplets (see Figure 4). By comparison, the concentration of PEG chains near the droplet center is only slightly below the average. (iii) As expected, the effect of confinement is significantly more pronounced in the smaller droplet, since the region close to the droplet surface occupies a larger volume fraction, i.e., 1 − (R − 2R D g ) 3 /R 3 ≈ 0.45 for R = 260 nm vs. ≈ 0.33 for R = 380 nm.  We analyzed the shape of the confined Dextran chains by computing their radius of gyration tensor with ∆r i being the position of monomer i relative to the polymer's center of mass. The root-mean-square radius of gyration is then R g ≡ R 2 g 1/2 = G n + 2G t 1/2 , where G n and G t are the components of G normal and tangential relative to the droplet surface, respectively. Figure 5 shows these components for monodisperse Dextran vs. the distance between the polymer's center of mass and the droplet surface. Polymers at distances larger than ≈ 2 R g have isotropic shapes, G n 1/2 = G t 1/2 ≈ 12 nm, whereas they become increasingly compressed along the normal direction as they approach the droplet surface (the tangential component is nearly constant). Further, chain conformations in the small and large droplets are almost indistinguishable.  Figure 5. Radial profiles of the radius of gyration of monodisperse Dextran, R g , and its normal ( G n 1/2 ) and tangential ( G t 1/2 ) components relative to the droplet surface. Data plotted against the distance between the droplet surface and the polymer's center of mass, r . Solid and dashed lines show results for R ≈ 260 nm and R ≈ 380 nm, respectively.
The confinement-induced increase in the polymer concentration near the droplet center could induce a phase separation of the PEG and Dextran chains there, as the local polymer concentration might cross the binodal of the mixtures (cf. Figure 3). Figure 6 shows simulation snapshots for the monodisperse and polydisperse systems confined in large droplets: while the monodisperse system appears to be still in the mixed one-phase regime, we can clearly see large Dextran aggregates in the polydisperse case, which is in excellent agreement with recent experiments by Watanabe et al. [18]. Looking more closely at the simulation snapshot, we can see that the aggregates in the polydisperse systems primarily consist of longer Dextran chains, whereas the shorter ones still remain well dispersed. Zhao et al. observed a similar molar mass fractionation in aqueous two-phase polymer solutions of PEG and Dextran [7], where the longer Dextran chains accumulated in the Dextran-rich phase, while the shorter Dextran chains were contained in the PEG-rich phase.
To better understand the distinct differences between the behavior of the monodisperse and polydisperse mixtures, we determined the probability P(r) of finding Dextran chains at center-of-mass position r, itemized by their molecular weight. These results are shown in Figure 6b,c for both droplet sizes, revealing a distinct spatial fractionation of the Dextran chains: short polymers with M n < 250 kg/mol are distributed almost homogeneously throughout the droplet and also come much closer to the droplet surface compared to the longer chains. In contrast, longer Dextran chains are moving to the droplet center to maximize their conformational entropy. These findings are consistent with recent experiments [18], where Watanabe et al. inferred from surface tension measurements that short Dextran chains (M w 500 kg/mol) preferentially adsorbed to the droplet surface. This radial partitioning of short and long Dextran chains promotes their phase separation, as longer chains have a distinctly smaller entropy of mixing compared to their shorter counterparts [32]. M n < 250 kg/mol 250 kg/mol < M n < 500 kg/mol 500 kg/mol < M n < 750 kg/mol M n > 750 kg/mol Monodisperse Finally, we quantified the size of the aggregates through the radial concentration profile c D (r) of Dextran monomers belonging to the largest cluster, and through the average number of Dextran monomers in a cluster, M , as explained in Section 3.1 above. To establish a reference, we performed additional simulations of ideal mixtures by setting χ ij = 0 and determined the corresponding mean aggregation number M 0 . The concentration profiles are shown in Figure 7, revealing that the Dextran clusters in the polydisperse confined systems are much larger and more concentrated compared to the monodisperse systems in bulk as well as in confinement.  Table 1 summarizes the results of M for all droplet simulations, which will be discussed in the following: For ideal mixtures with monodisperse Dextran chains, we find M 0 ≈ 100, which is comparable to the number of monomers per Dextran chain (N D = 80). This value is sensible given that χ ij = 0 and that the search radius of the clustering algorithm is slightly larger than the average segment length of our Dextran model (b ≈ 5.5 nm). For polydisperse Dextran in the smaller droplets (R = 260 nm), the mean aggregation number increases to M 0 ≈ 380 due to an accumulation of longer Dextran chains in the droplet center, even at ideal conditions. This effect is considerably less pronounced in the larger droplets (cf. Figure 6), where we find M 0 ≈ 110 instead. For the non-ideal mixtures with monodisperse Dextran, we find M ≈ 1300 − 1600, which indicates the existence of small aggregates consisting of 15 − 20 Dextran chains (there are, in total, 348 and 1175 Dextran chains in the small and large droplets, respectively). The mean aggregation number becomes significantly larger in polydisperse PEG-Dextran mixtures, reaching values up to M ≈ 10, 000. Interestingly, M is about two times smaller in the smaller droplets (see Table 1), which is likely a finite-size effect, as the smaller droplets contain about three times fewer Dextran chains compared to the large droplets. Nevertheless, both systems show clearly that the confinement-induced fractionation of short and long Dextran chains drives phase separation, as observed in recent experiments [18].

Conclusions
To better understand the phase behavior and conformations of (bio)polymers in droplets, we have performed coarse-grained molecular simulations of binary polymer mixtures in spherical confinement. We have parameterized our model to closely mimic the behavior of aqueous PEG-Dextran mixtures and have considered both monodisperse and polydisperse molecular weight distributions for the Dextran chains. Simulations have been conducted for two droplet sizes at polymer concentrations lying in the mixed one-phase region in bulk. In spherical confinement, we have found a distinct spatial separation of the polymers by length, with the longer ones accumulating at the droplet center to maximize their conformational entropy. Furthermore, chains near the droplet surface became increasingly compressed along their normal direction. This confinement-induced partitioning was much more pronounced in the polydisperse systems and caused the phase separation of the two polymer species at average polymer concentrations where the bulk system was still in the one-phase regime.
Although we have chosen the model parameters to replicate PEG-Dextran mixtures, the rather generic nature of our coarse-grained model makes our results applicable to a wide range of different polymer mixtures. Our simulations demonstrate how the distribution of polymers is affected by confinement effects, even at good solvent conditions, with longer chains moving to the droplet center to maximize entropy. The resulting spatial inhomogeneity can drastically alter the phase behavior of the confined polymers, which is important for understanding, e.g., the liquid-liquid phase separation of biopolymers in cellular environments. Furthermore, our simulations provide useful guidelines for the fabrication of polymer-loaded droplets. For example, by tuning the interactions between the droplet surface and the different polymer species, one can either enhance or suppress their spatial separation and thus control the resulting phase behavior and surface tension.

Institutional Review Board Statement: Not applicable.
Informed Consent Statement: Not applicable.

Data Availability Statement:
The data that support the findings of this study are available from the corresponding author upon reasonable request.