Janus Droplet Formation via Thermally Induced Phase Separation: A Numerical Model with Diffusion and Convection

Microscale Janus particles have versatile potential applications in many physical and biomedical fields, such as microsensor, micromotor, and drug delivery. Here, we present a phase-field approach of multicomponent and multiphase to investigate the Janus droplet formation via thermally induced phase separation. The crucial kinetics for the formation of Janus droplets consisting of two polymer species and a solvent component via an interplay of both diffusion and convection is considered in the Cahn–Hilliard–Navier–Stokes equation. The simulation results of the phase-field model show that unequal interfacial tensions between the two polymer species and the solvent result in asymmetric phase separation in the formation process of Janus droplets. This asymmetric phase separation plays a vital role in the establishment of the so-called core–shell structure that has been observed in previous experiments. By varying the droplet size, the surface tension, and the molecular interaction between the polymer species, several novel droplet morphologies are predicted in the development process of Janus droplets. Moreover, we stress that the hydrodynamics should be reckoned as a non-negligible mechanism that not only accelerates the Janus droplet evolution but also has great impacts on the coarsening and coalescence of the Janus droplets.


INTRODUCTION
The Janus droplet is an anisotropic material that is composed of two hemispheres possessing different physical, chemical, biological, and electric properties. 1−7 Due to its anisotropic physicochemical interaction with the surrounding fluid, the Janus particle has captured extensive scientific attention and exhibits an extraordinary potential in many fields, such as display, 8,9 microsensing, 6,10,11 drug delivery, 12−14 and micromotor. 15−20 Pertaining to the prospect for the broad application of Janus droplets, several manufacturing methods, including droplet microfluidics, 1,21 flow lithography, 22 electrohydrodynamic co-jetting, 23,24 and micro-molding, 25,26 have recently been explored. Besides these above-mentioned effective fabrication techniques, Janus droplets can also be produced via thermally induced phase separation (TIPS), which has been widely investigated and utilized due to its convenient setup, high productivity, and low equipment cost. 27−29 However, there are several disadvantages of the TIPS method, which limit its utilization in some highly sophisticated fields, for example, high droplet polydispersity, sustainable reproducibility, and diverse morphologies. 27 In order to improve the stability for the production of the Janus particles via TIPS, it deserves more extensive studies to understand the kinetics for the formation of the Janus particles.
In previous research studies, many complex systems are experimentally investigated for both inorganic 30−33 and organic 12,34−39 materials. Considering the complicated recipes adopted, as well as the huge amount of adjustable material parameters, such as surface tension, viscosity, critical temperature, and diffusivity, it is almost impossible to avoid repetitive labor works for achieving a comprehensive understanding of the final morphologies of the Janus droplets. Meanwhile, the crucial but transient intermediate states during the Janus droplet evolution are highly associated with the kinetics of the system, which can be strictly restricted by the experimental methods and equipment. With the advance in computational technology, the material simulation methods provide us another prospect to carefully heed the behaviors of Janus particles. For instance, the hydrodynamics of the Janus droplet and its stability have been simulated via fluid dynamic methods, including the Lattice−Boltzmann method 40−44 and the arbitrary Lagrangian−Eulerian method. 45 The selfassembly properties of the Janus droplets are mainly discussed with the help of molecular dynamics, 2,46−52 dissipative particle dynamics, 53−57 and Brownian dynamics. 58−62 In this work, we adopt the phase-field method to elucidate the Janus formation via TIPS, which has been widely introduced to simulate the phase separation of several multicomponent systems. 63−67 Most crucially, the microscopic material evolution can be investigated in detail by solving a series of rigorous governing equations obeying thermodynamic principles. Meanwhile, the microscopic mass transfer during the phase separation inside the Janus droplet, such as diffusion and convection, can be extensively analyzed, which may be difficult to capture under microscopes. Moreover, by adjusting key parameters, including the droplet size, initial polymer composition, interfacial tension, and molecular interaction of the Janus system, we can scrutinize the underlying mechanisms that dominate the kinetics of the Janus droplet formation. Taking the advantages of computational simulation, several interesting morphologies are observed, which may help to design and manipulate the micro-droplet with complex structures for applications in biomedical, pharmaceutical, and food industries. In addition, we explore the crucial impact of hydrodynamics on the coarsening and coalescence mechanisms of the Janus droplets, which will shed light on the control of droplet polydispersity and help to design optimum conditions for the production of Janus particles.

NUMERICAL MODEL
In this work, we apply a phase-field model of Cahn−Hilliardtype coupling with the Navier−Stokes equations to simulate the formation process of Janus droplets. Janus droplets have been experimentally observed in polymer solutions consisting of polymer species A and B, and a solvent species (like water) with molecular weight much less than those of A and B. 33 A Janus particle usually has two distinct faces that are immiscible with each other as well as with the solvent-rich matrix. Each face of the Janus particle is composed of distinct polymer species. 28 The space-and time-dependent composition of polymer A, polymer B, and the solvent is defined as c 1 (x, t), c 2 (x, t), and c 3 (x, t), respectively, and is subjected to the constraint of the incompressible condition ∑ i=1 3 c i = 1. The starting point of our model is the free energy functional of the Cahn−Hilliard type, 68 where Ω indicates the spatial domain occupied by the system. The first term inside the bracket denotes the bulk free energy density, with σ* being the characteristic surface tension. The modeling parameter ϵ controls the interface width between different phases. The second term depicts the gradient energy density, where the coefficient σ i is related to the interfacial tension between the immiscible phases in the ternary A−B− solvent phase diagram.
To obtain the phase diagram for the ternary system, we adopt the Flory−Huggins model 69 with the following free energy density formulation r  1   3   , 1,2,3   123  1 2 3 (2) Here, T denotes the temperature, T r is a reference temperature that is chosen to be the half of the critical temperature T c for the binary A−B miscibility gap, and N i indicates the degree of polymerization for component i. The Flory parameters χ ij and χ ijk are supposed to be temperature dependent as follows 70,71 The degree of polymerization N i and the temperature coefficients for the Flory parameters a ij , a ijk , b ij , and b ijk are chosen to be consistent with the following experimental observations: (i) when the temperature is greater than the critical temperature T c , that is, T > T c , polymers A and B are well miscible with each other and form a homogeneous droplet. The homogeneous droplet consisting of polymers A and B is immiscible with the solvent matrix, which is attributed to a relatively strong repulsive interaction between the solvent and the two polymer species. (ii) As T < T c , the repulsive force between polymers A and B increases with a reduction in temperature. The enhanced repulsive force at low temperatures leads to a miscibility gap involving polymer A-rich and polymer B-rich phases, which results in the phase separation inside the homogeneous droplets formed above the critical temperature. Such kind of features [(i) and (ii)] for the phase diagram can be found in the ternary system of hexane−perfluorohexane− water 33 and many other ternary polymer solutions. The parameters N i , a ij , b ij , a ijk , and b ijk for modeling the ternary phase diagram are tabulated in Table 1.
The phase diagram is constructed by solving the equation system for the binodal composition pair c eq = (c 1 eq ,c 2 eq ,c 3 eq ) and c eq * = (c 1 eq *,c 2 eq *,c 3 eq *) as : The spinodal composition c s = (c 1s , c 2s , c 3s ) is defined by the locus of  Two typical phase diagrams at temperatures greater and less than T c are shown in Figure 1I. The binodal and spinodal lines are depicted by the solid and dashed lines, respectively, in the contour plot for the bulk free energy density. The temporal evolution for the compositions follows the Here, the mobility M ij is assigned with the Onsager's relationship 72,73 as M ij = M 0 c i (δ ij − c j ) (see the derivation in ref 74), in which δ ij is the Kronecker delta. The parameter M 0 scales the mobility and is formulated as ∑ i=1 3 D i c i , where D i stands for the diffusivity of pure component i. In this work, we assume D 1 = D 2 = D 3 = 1.0 × 10 −9 m 2 /s. Moreover, the Gaussian white noise term ξ i for the diffusion flux is essential to trigger the spinodal decomposition, which follows where A ξ is the amplitude of the perturbation. Concerning the hydrodynamic effect in the system, we use incompressible Navier−Stokes equations (eqs 8 and 9) to calculate the flow velocity u in the ternary system. The parameters ρ, η, and P denote the fluid density, viscosity, and pressure, respectively. The surface tension effect is coupled into eq 8 by the Korteweg stress tensor as Ä 2.1. Non-dimensionalization. All the physical parameters are non-dimensionalized by the characteristic length x* = 2 × 10 −10 m, reference surface tension σ* = 1.0 × 10 −2 N/m, and diffusivity D* = 1 × 10 −9 m 2 /s. Thereafter, we have the following scaling factors for the physical parameters, as shown in Table 2 (more details on non-dimensionalization are available in ref 75).
Substituting the scaling factors into the CHNS model (eqs 7−9), we obtain ∇̃* · *̃= x u u ( ) 0 (12) After simplification, the non-dimensionalized form of the Cahn−Hilliard−Navier−Stokes equation reads Ä The dimensionless quantities Re, We, and Péare calculated as As we are interested in the phase separation coupled with fluid flow, the convection and diffusion processes become comparable. In other words, the relaxation time for the thermodynamic process (x*) 2 /D* is in the same magnitude with hydrodynamics x*/u*, which gives rise to the Pećlet x* characteristic length 2 × 10 −10 m σ* characteristic interfacial tension Langmuir pubs.acs.org/Langmuir Article number Pé= u*x*/D* = 1.0, while Re and We vary with the magnitude of the hydrodynamic effect. In our work, the finite difference method and the explicit Euler scheme are implemented to solve the evolution equations eq 13 with the equidistant Cartesian mesh. Periodic boundary conditions are applied for the composition and fluid velocity.

RESULTS AND DISCUSSION
3.1. Formation of Janus Droplets. In the following discussions, we consider the scenario where the phase separation via the diffusion process dominates the Janus droplet formation and the weak convection can be neglected as reported. 76 Initially, at the temperature T = 1.05T c , a homogeneous solution of polymers A and B is dissolved in the abundant solvent (c 1 /c 2 /c 3 = 15/15/70). The initial composition of the ternary polymer solution is marked by the triangle symbol in Figure 1I. As shown in the phase diagram, this initial composition is located between the black dashed lines, that is, inside the spinodal region. A Gaussian white noise with amplitude A ξ = 0.1 for the diffusion flux of all components, as depicted by eq 13, is used to model the composition fluctuation. Triggered by the composition fluctuation, phase separation begins and polymer-rich droplets (light green) start to form and grow from the solvent-rich matrix (dark purple), as illustrated in Figure 1IIi. With time, these droplets coarsen with each other dominated by the Ostwald ripening effect due to their radius difference. Meanwhile, the coalescence of relatively large droplets is observed as highlighted by the dashed square. Consequently, dispersed droplets with various sizes are obtained. The equilibrium composition in the droplets is marked by the green circle on the binodal line in the left phase diagram. Because at T > T c , polymers A and B are miscible with each other, the Janus particle is not observed at this temperature. As demonstrated in the following section, the droplet with miscible polymers A and B formed at T > T c is the precursor for the formation of the Janus droplet.
Next, we reduce the temperature till T = 0.95T c , where the miscibility gap between the two polymer species begins to appear, which results in the spinodal decomposition inside the precursor polymer droplet. To simulate the Janus droplet formation and to ensure a stable three-phase region, we assign the following interfacial tension parameters σ 1 = σ 2 = 0.64 and σ 3 = 3.6. As presented in the snapshots of Figure 1IIii, at t = 5 s, the two polymer species get separated to form a polymer Arich region (blue) and a polymer B-rich region (red). These two regions are divided by the light green interfaces. Accompanied by the proceeding spinodal decomposition, the interfaces between the blue and the red regions are stabilized (see the snapshot at t = 10 s). Decided by the wetting effect at the triple junction involving the polymer A-rich phase, polymer B-rich phase, and the solvent-rich phase, Janus droplets with r < ≈ 15 μm are generated at t = 10 s. Afterward, the phase separation subsides, and the minimization of free energy functional is mainly manifested by the reduction in the interfacial energy. Then, the Ostwald ripening effect and coalescence come on to the stage. For instance, the tiny droplet in the white circle vanishes, because of its huge curvature difference from its surrounding large droplets. As highlighted by the orange-dashed squares in Figure 1IIii, three tiny Janus droplets coalesce into a large one due to their relatively narrow apart distances. Moreover, a large Janus droplet with a transient hamburger structure 77,78 is captured in the white-dashed circle at t = 50 s. This observation indicates that the size of the precursor polymer droplet plays an important role in the morphological evolution of the Janus droplet, which will be discussed in the following sections.
3.2. Radius and Initial Composition. In the manufacturing process of Janus particles, the final product with different sizes can be achieved by controlling the holding time at T = 1.05T c above the critical temperature. As shown in Section 3.1, the Ostwald ripening effect together with the droplet coalescence results in precursor droplets with various sizes. However, with the descending temperature, the phase separation of polymers A and B cannot simultaneously engender the perfect Janus droplet for differently sized precursor droplets. In order to have a heedful look at this size effect, we place three precursor droplets with the initial radius r 0 = 20, 40, and 80 μm in the solvent-rich matrix and reduce the temperature to T = 0.95T c . As shown in Figure 2, the smallest droplet with r 0 = 20 μm forms a perfect Janus at 50 s, while the largest precursor droplet with r 0 = 80 μm spends 20 times longer for the transformation into a perfect Janus particle. The reason for this size effect can be explained by viewing the intermediate stage of the Janus droplet development.
As depicted in Figure 2, more complex bicontinuous structures are produced inside the relatively large droplets by phase separation. The final perfect Janus droplet is formed via the coalescence of the tortuous bicontinuous structure. As demonstrated in the third snapshot of Figure 2, for all the three phases, namely, the blue-highlighted polymer A-rich phase, the red-colored polymer B dense phase, and the matrix, the equilibrium compositions are reached. Hence, at this stage, the phase separation plays hardly any role now and the Ostwald ripening is the main mechanism for the Janus droplet evolution. According to the classic theory of the Ostwald ripening, 79−81 the mean average radius of the phase separation structure follows the quintessential LSW (Lifshitz−Slyozov− Wagner) power law r ̅ ∼ t 1/3 . Therefore, the transformation of the droplet with r 0 = 80 μm into a Janus shape takes about 4 3 times longer than the tiny one with r 0 = 20 μm.
In some experiments, 82−85 non-equal compositions of polymers A and B are adopted to synthesize different droplet morphologies, for example, vesicle, uneven Janus, and so forth. In Figure 3I, we perform simulations with unequal polymer compositions, namely, c 1 /c 2 = 3:7 and c 1 /c 2 = 7:3, for a precursor droplet with r 0 = 80 μm. Because the composition ratio between the two polymers is off 5:5, in lieu of bicontinuous structures, sub-droplets are engendered inside the precursor droplet via the spinodal decomposition. To minimize the interfacial energy, those sub-droplets merge into a joint phase. This routine for the Janus droplet formation is in good agreement with the confocal microscopy images for the dextran−PEG−water system, as illustrated in Figure 3II. It is also noteworthy that the interface between polymers A and B in our simulation (7:3 in Figure 3I) has a different convexity from the experiment (PEG-riched case in Figure 3II). As observed in the hexane−perfluorohexane Janus system, 33 , the Janus structure is extremely sensitive to the compositiondependent interfacial tension parameter σ i . Thus, the slight  Langmuir pubs.acs.org/Langmuir Article polymer composition changes can result in the huge difference in the interface convexity.

Interfacial Tension.
In the previous part, we elucidated the formation of Janus droplets via liquid−liquid phase separation, which transforms the bulk free energy of the mixture into the interfacial energy. The latter one determines the pattern of the microstructure around the triple junction. In this section, we address the microstructure pattern at the region of the triple junction by varying the interfacial tension. In experiments, the adjustment of interfacial tensions is achieved by adding various surfactants into the system. 33 As a rule of thumb, the surfactant quantity less than 0.1% v/v has a minor effect on altering the ternary phase diagram. 76 In our simulations, the effect of specific surfactants can be mimicked by changing the interfacial tension parameter σ i . To be consistent with the sharp interface model, the interfacial tension between the component A-rich phase and the component B-rich phase in the Cahn−Hilliard approach (see the right column of Figure 4I) is calculated as where s delineates the integral routine that follows the energy minimum principle. As displayed by the solid arc with a white open arrow in Figure 4I, the integration starts from the white open dot with one bulk equilibrium composition (c 1 eq ,c 2 eq ,c 3 eq ) inside phase α to the other bulk equilibrium composition (c 1 eq *,c 2 eq *,c 3 eq *) inside phase β. The term Δf is defined as  (16) which measures the excess free energy density referring to the equilibrium state f(c 1 eq ,c 2 eq ,c 3 eq ) + ∑ i=1 3 μ i eq c i eq . From this expression, the interfacial tension of the system not only is decided by the interfacial tension parameters σ i but also relies on the bulk free energy density of the mixture. Owing to the addition of the third component (C in this example), the integral routine at the equilibrium state does not take the dotdashed straight line from α to β in the phase diagram but follows the deterministic path that corresponds to the fundamental energy minimization principle as addressed by the Cahn−Hilliard equation (eq 7). In this way, the interfacial tensions in the ternary system can deviate largely from the one in the binary system, crucially depending on the Flory parameters χ ij and χ 123 . It appears that the addition of the third component dramatically increases the complexity for measuring and calculating the interfacial tension. In most cases, there is a paucity of experimental data for the interfacial tensions to validate the model. Here, we systematically discuss the effect of the interfacial tension on the pattern formation of the Janus particles.
First, we analyze the interfacial tension parameter σ i affecting the formation of the Janus droplet. The interface tension parameters σ 1 and σ 2 both vary from 0.2 to 4.2, while all other parameters are the same as that in Section 3.2. This setup can be comprehended as introducing a tiny amount of a surfactant that can drastically reduce the interfacial tension of the system, while the phase diagram is not altered. As displayed in Figure  4II, an increase in σ 1 and σ 2 results in different morphological evolutions during the phase separation. For instance, the threefold symmetric droplet is observed in the setup with σ 1 = σ 2 = 1.0, which shows high similarities with the experimental observations for the patchy droplet structures of self-assembly materials. 86 With an increase in σ 1 and σ 2 , the time expense for the formation of the Janus droplet via the Ostwald ripening is sharply reduced and the final contour of the Janus droplet transforms from a joint sphere to a quasi-ellipse. In Figure 4III, the contact angle θ between two Janus hemispheres is measured with two different methods. (1) Orange open dots: the interfaces between the Janus hemispheres and the solvent matrix are extracted from the simulation by the criterion: c 1 = 0.5 for the polymer A hemisphere and c 2 = 0.5 for polymer B. As sketched by the inset in Figure 4III, we fit the interfaces between the polymer A-rich hemisphere and the solvent matrix and between the polymer B-rich hemisphere and the solvent matrix with circles. The interfaces intersect at the triple junction where the contact angle θ is computed by the included angle of the two tangent lines of the fitted circles.
(2) Black crosses: with eq 15, the three interfacial tensions form the so-called Neumann triangle relation, which constrains the θ at equilibrium as Showing good agreement in these two methods, the contact angle θ decreases by enlarging the interfacial tension parameter σ 1 (=σ 2 ). The comparison not only indicates the droplet Langmuir pubs.acs.org/Langmuir Article contour changing from Janus (θ ∼ 180°) to quasi-dumbbell structure (θ ∼ 155°) controlled by the interfacial tension but also validates the Cahn−Hilliard model for the Janus droplet system. Second, we assign a different Flory parameter χ 12 with the same interfacial tension parameters σ 1 = σ 2 = 1.0 in the simulation. This setup may be achieved by applying different polymer species in the system. As demonstrated in the free energy landscapes in Figure 5I, an increase in χ 12 not only modifies the interfacial tension between polymers A and B but also broadens the miscibility gap, as depicted by the light red shadow regions. Hence, a more pronounced phase separation accelerating the production of polymer A/B interfaces is expected for a larger value of χ 12 , as demonstrated in Figure  5II. Most interestingly, in the blue highlighted simulation snapshots of Figure 5II, intermediate droplets with polygon and triangle shapes are observed. This observation implies that the thermally induced liquid−liquid phase separation may have the potential to synthesize droplets with some special morphologies, which shall be investigated in the future.
3.4. Asymmetric Phase Separation. In the previous sections, the Flory parameters χ 13 and χ 23 are assumed to be the same, so do the interfacial tension parameters σ 1 = σ 2 . These setups usually give rise to symmetrical kinetics for the formation of the Janus droplet, that is, different faces of the Janus droplet and the diffusion path of distinct polymer species following identical kinetics. Actually, in most real systems, the Flory parameters of polymer A and polymer B with the solvent are not equal, for example, χ 13 ≠ χ 23 , denoting different attractive/repulsive potential of the polymer species with the solvent molecules. In addition, the interfacial tensions between polymer species and the solvent should be distinct from each other, which may be caused by unequal interfacial tension parameters, σ 1 ≠ σ 2 . The Janus droplet formation is attributed to the interplay of both liquid−liquid phase separation and the interfacial energy minimization, the kinetics of which is essentially affected by both the Flory parameters and the interfacial tension parameters. Thus, the distinguishing properties of polymers A and B may lead to asymmetrical kinetics for the microstructure evolution of the Janus droplet formation as well as a dissymmetric diffusion path, which will be discussed in this section.
To elucidate the symmetric kinetics, we place one homogeneous droplet with r 0 = 30 μm and initial composition c 1 = c 2 = 0.47 inside the polymer lean matrix (c 1 = c 2 = 0.00025). In this setup, we adopt the same interfacial tension and Flory parameters for polymers A and B, that is, σ 1 = σ 2 = 0.6 and χ 13 = χ 23 = 4.2. From eq 15, we have σ 13 = σ 23 . Triggered by the composition fluctuation, the phase separation generates a polymer A-rich region (blue colored) and a polymer B-rich region (red colored), which comprise the bicontinuous structure resulting from the spinodal decomposition, as displayed in Figure 6Iii. The kinetics of the Janus droplet formation is characterized by tracing the maximal polymer A composition c 1 (turquoise line) and the maximal polymer B composition c 2 (red line) starting from the green circle in the phase diagram. As can be noticed in Figure 6Ii, once the polymer A-rich blue region forms, the turquoise trajectory of the maximal c 1 falls exactly on the black solid binodal line. This overlap implies that there exist a series of transient pseudo-binary equilibrium states between the polymer A-dense blue region and the polymer lean solvent matrix (black open square). One typical exemplary pseudobinary equilibrium is shown by the blue dot-dashed tie line. This kind of pseudo-binary equilibrium is unstable over time because of the presence of the adjacent polymer B-rich red region, which is not in equilibrium with the polymer A-dense region. Driven by the free energy minimization, the immiscible polymer A in the adjoining red region is rejected. For the low solubility of polymers in the solvent matrix, phase separationinduced mass transformation via the surrounding matrix is in vain. The precipitated polymer A from the red region can diffuse only across the A−B interface, resulting in an enrichment of polymer A in the blue-colored regions. Consequently, the maximal composition c 1 continuously moves leftward along the binodal line and finally reaches the ternary equilibrium composition labeled by the blue Janus symbol in the phase diagram in Figure 6Ii. At this right stage, a perfect Janus droplet is completely formed. The same kinetics happens inside the polymer B region, because of the exactly equal interfacial tension and Flory parameters of polymers B and A. As depicted in Figure 6Iiii, we observe the same tendencies of the maximal polymer A composition (blue solid line) and maximal polymer B composition (red solid line) changing with time. The small ups and downs are in line with the coarsening and coalescence of the phases that change the curvature and consequently the composition. In this way, we denominate this kind of kinetics as "symmetric", which denotes the identical kinetic pathways for polymers A and B.
In our second simulation, we set the interfacial tension parameters to be asymmetric for polymers A and B as σ 1 = 1.0 and σ 2 = 0.6. The Flory parameters are χ 13 = χ 23 = 4.2, as the one in Figure 6Ii. With eq 15, it can be shown that the interfacial tension between the polymer B-rich phase and solvent-rich phase σ 23 is smaller than that between the polymer A-rich phase and solvent-rich phase σ 13 . Hence, at the very beginning of the phase separation, in order to reduce the interfacial energy, polymer B prefers to contact with the solvent matrix and encircles polymer A to prevent the formation of the interface between polymer A and the solvent Langmuir pubs.acs.org/Langmuir Article matrix. Such a kinetic process leads to centrosymmetric diffusion fluxes of polymers A and B. Accordingly, an intermediate core−shell morphology is established, which has been observed in many micro-capsule systems. 37,87−90 As shown in Figure 6IIi for the kinetic pathway in the phase diagram, the core−shell morphology is accompanied by a deviation in the red trajectory (the maximal c 2 ) from the black solid binodal line. This deviation is explained as follows. From the intermediate state M 1 to M 2 , due to spinodal decomposition, the blue polymer A-rich core forms and quickly reaches the equilibrium composition on the binodal line. However, in the shell region, the kinetics is totally different. Because both the polymer A-rich phase and the solvent matrix contact with the shell, there exist two interfaces, namely, the polymer A−B interface and the polymer B−solvent interface. The composition in the shell region evolves via not only the phase separation with respect to the polymer-A rich core but also the binodal decomposition pertaining to the solvent matrix. The non-equal diffusion fluxes at the core−shell interface via spinodal decomposition and at the shell−solvent interface via binodal decomposition give rise to the deviation in the maximal c 2 trajectory from the binodal line in the polymer B-rich region.
In Figure 6II, between the time points M 2 and M 3 , the maximal c 2 hardly increases with time and the core−shell structure stagnates. The reason for this utterly slow morphological transformation can be explained by the kinetic pathway in the phase diagram in Figure 6IIi. At M 2 , the red trajectory passes through the black dotted spinodal line and the composition in the shell falls outside the spinodal region. Thereafter, the microstructure evolution is dominated by the binodal decomposition instead of spinodal decomposition. Because the former is much slower than the latter, the composition c 2 in the shell region increases with time at a relatively slow rate, leading to a stagnation state from M 2 and M 3 . Noteworthily, the core−shell structure is not the energy minimal state. In the shell region, the everlasting Gaussian composition fluctuation gives rise to small composition gradients and produces a circular diffusion for the polymer B component, leading to the breakup of the core−shell structure at M 4 . When the core−shell structure collapses, the interfacial tension rapidly reshapes the polymer B-rich layer into a hemisphere (at E), which is the final equilibrium state of the system. We name the observation in this setup "asymmetric kinetics" to make the distinction from the previous "symmetric" case, which has also been studied in many experiments, 83−85 as the microscopic images shown in Figure  6III.
3.5. Fluid Dynamics. In the previous Sections 3.1−3.4, we considered the Janus droplet formation via the diffusioncontrolled phase separation process and observed a stagnation stage for the intermediate core−shell microstructure in the formation process of the Janus droplet. In this section, we will shed light on the effect of hydrodynamics on the morphological evolution of the Janus particle as well as the core−shell microstructure. Considering the kinetics of the polymeric phase separation, the polymer diffusivity can vary from 10 −9 to 10 −12 m 2 /s, which is decided by the degree of polymerization, temperature, solvent property, and so on. Thus, when the polymer system becomes semidilute, or even dense, the diffusivity of polymer chains decreases drastically. When D* drops to the order of 10 −12 m 2 /s, the hydrodynamic effect plays a non-negligible role compared with diffusion. As written in eq 11, the composition gradient ∇c i in the inhomogeneous fluid system can propel the fluid flow via the capillary tensor term Θ in the Navier−Stokes (NS) eq 8. This mechanism is called the Marangoni effect. By adjusting the dimensionless numbers, (1) Weber number: We and (2) Reynolds number: Re, in the non-dimensionalized NS eq 13, we will show distinct hydrodynamic behaviors for the microstructure evolution of the Janus droplet.
In a first setup, we perform a single droplet simulation with the asymmetric interfacial tension σ 13 > σ 23 as discussed in Section 3.4. To couple the hydrodynamic effect with the phase separation, the Weber number We and Reynolds number Re are both set as 1. As can be noticed in both the morphological evolution and the kinetic pathway in Figure 7, the transient core−shell structure survives 10 times shorter than the simulation without the hydrodynamics shown in Figure 6. By viewing the pressure distribution around the droplet in Figure  7Iii, the pressure has its maximum (dark red colored) in the polymer B-dense shell region, the minimum inside the solvent matrix, and the intermediate value (yellowish) in the polymer A-rich core. Here, the pressure is induced by the surface tension force subjected to the incompressible condition. In order to reduce the surface energy, the convection resulting from the pressure difference on both sides of the shell structure together with the diffusion leads to the breakup of the shell, which is faster than the morphological transformation solely via diffusion. The stagnation period of the core−shell structure around the spinodal composition lasts only several seconds in this case.
The Marangoni flow has two origins: (i) phase separation driving force. During the phase separation process, the composition has not yet reached the equilibrium value, where excess free energy density of the mixture results in the surface tension force ∇·Θ in the N−S equation. As demonstrated in eq 10, the term ∇·Θ is proportional to the chemical potential gradients, which cannot be balanced by the pressure p around the droplet. Consequently, the convection takes place. (ii) Minimization of the interfacial energy. The phase separation leads to the creation of new interfaces. The Langmuir pubs.acs.org/Langmuir Article non-uniform curvature along the interface as well as the surface tension force enforcing the Young's contact angle at the triple junction also gives rise to a convection. The former occurs at the early stage of the phase separation. The latter appears once new interfaces are established and dominates the evolution when the bulk composition reaches the equilibrium value. It is noteworthy that both mechanisms also exist when the Janus droplets are produced by the phase separation solely via diffusion (see Figure 1). The main difference is the way that the energy is minimized or, more suffice it to say, the energy dissipation associated with different kinetics. When the hydrodynamics is coupled, the excess free energy and surface energy can be transformed into the kinetic energy of the fluid flow. As displayed in Figure 7Iiii, the fluid velocity around the droplet increases simultaneously with the proceeding spinodal decomposition. After reaching the equilibrium composition and the breakup, the non-uniform curvature and the free energy minimization at the triple junction lead to a further increase in the fluid velocity. When the curvature becomes uniform around the Janus droplet and the contact angle at the triple junction reaches the equilibrium value, the fluid flow dissipates and tends to vanish, as shown in Figure 7I at t = 16 s. Next, we simulate the Janus droplet formation via a two-step phase separation. At step (i), homogeneous polymer droplets are produced at T = 1.05T c . At step (ii), the temperature is reduced to T = 0.95T c , where Janus droplets are formed. Here, two magnitudes of the hydrodynamic effect are considered, namely, We = 1 and 0.01, to elucidate its influence on the multi-droplet system. By the definition of the Weber number, a large We indicates a strong capillary effect. By comparing the morphological evolution in Figure 8 with the previous case only with diffusion (Figure 1), the phase separation process is magnificently accelerated. For the phase separations at a high temperature T = 1.05T c , it always takes less time t e for the droplet to reach its equilibrium concentration when the hydrodynamic effect becomes more pronounced. For We = ∞ (diffusion only), t e = 8.5 s; We = 1.0, t e = 2.0 s; We = 0.01, t e = 0.5 s. As T = 0.95T c , concerning the size effect discussed in Section 3.2, droplets with similar sizes are compared. The Janus droplet production also becomes quicker with the stronger convection. For We = ∞, t e = 30.0 s; We = 1.0, t e = 18.0 s; We = 0.01, t e = 2.25 s.
Also, the size distribution of the Janus droplet becomes more dispersed with the enhancing hydrodynamic effect. The underlying mechanism can be elucidated with the Tanaka− Golovin theory 91 that has two important contributions to the polydispersity of the Janus droplet. On the one hand, the hydrodynamic force is produced by the overlap of the diffusion potential and propels the droplet coalescence. Hence, the droplets with relatively small apart distances are prone to impinge with each other, forming a large droplet, as illustrated in the gold dashed squares of Figure 8. On the other hand, the hydrodynamics can act as a supporter for the tiny droplet and curtails the consumption of small droplets by its large peers  Figures 1 and 8 confirms the supporting effect of hydrodynamics. The pronounced coalescence effect for large droplets and the shielding effect for the small droplets result in the polydispersity of Janus droplets. This observation is also consistent with the theoretical derivation. 94 In addition, because both the phase separation and the fluid mechanics are considered in our model, the Janus production via the microfluidic process can also be simulated and will be discussed in future works.
Finally, we simulate the Janus droplet formation where the surface tensions of polymers A and B are identical. As shown in Figure 9, the perfect Janus droplets with symmetric hemispheres of A and B are produced. Similar to the previous simulations with unequal surface tensions (Figure 8), the Janus sizes increase with the increase in the hydrodynamic effect. The size distribution also turns to be more dispersed by enhancing the convection inside the system, which is attributed to the two reasons discussed in the previous paragraph: (i) the imbalanced surface tension force ∇·Θ induces convection that propels the droplet coalescence; (ii) the strong Marangoni flow also acts as a counter-force that curtails the diffusion process. Under this circumstance, the LSW mechanism plays only a subordinate role and the mean radius of the Janus droplet does not follow the 1/3 scaling law anymore. 91

CONCLUSIONS
In summary, we have presented a multi-component and multiphase Cahn−Hilliard−Navier−Stokes model to simulate the formation of the Janus droplet via the thermally induced liquid−liquid phase separation. By considering the initial polymer composition, Flory parameters, and the interfacial tension parameters, we elucidate the competing mechanisms for the morphological evolution, namely, (i) the diffusiondominated phase separation process and (ii) the minimization of the interfacial energy. Due to the interplay of these two crucial factors, various transient morphologies of the droplets have been observed, including the polygon-shaped structures and patchy droplets, which need to be investigated in detail in future work. Most importantly, we stress the significance of the hydrodynamic effect. Owing to the curtailed diffusivity of entangled long polymer chains, the composition inhomogeneity-induced Marangoni effect becomes comparable with the diffusion of polymeric species. Hence, not only the phase separation process of the droplet but also the droplet coalescence is magnificently accelerated by the Marangoni flow, which results in a broad range for the size of Janus droplets. Further systematic computational studies of material Langmuir pubs.acs.org/Langmuir Article properties and pattern correlations will be considered in the forthcoming researches. As the approach is already capable to be applied in three-dimensional structures, the 3D simulation of the Janus droplet will be conducted and analyzed. We believe that the present research involving Janus droplet formation via diffusion as well as fluid dynamics will shed new light on the multi-phase microfluidic manipulation technology for wide applications.