Modeling droplet aggregation and percolation clustering in emulsions

2.1 Langevin impulse integrator

The time evolution of single droplets and aggregates is calculated by the integration of the Langevin equation. The motion of droplets and aggregates is described using the Langevin impulse integrator (Skeel and Izaguirre, 2002) of the form:

Drag acts on the droplet according to Stokes’ law. The friction coefficient ζ is given as

The attraction energy (U_m) with partially electromagnetic retardation can be written as (Israelachvili, 2011) follows:

The electrostatic repulsion energy between droplets (U_el) is expressed as (Derjaguin et al., 1987)

The attraction and repulsion forces are calculated by differentiating the equations of potential energy with respect to the coordinate. The force of gravity acting on single drops and clusters is also taken into account.

The Brownian motion is represented by a random force with components R₁ and R₂ that satisfy the following statistical characteristics (Skeel and Izaguirre, 2002):

Starting values of particle velocities are obtained according to the Maxwell velocity distribution at a given temperature:

The formation of aggregate is set to occur when droplets approach each other at a distance of less than the distance of coagulation (taken to be equal to 1 nm), with the droplet aggregation being assumed to be reversible. The sums of the attractive, repulsive, and gravity forces for each pair of droplets in each aggregate are taken to be equal to the sum of these forces at the distance of coagulation. This sum is kept unchanged even when droplets approach each other and the distance between droplets in an aggregate becomes smaller. If the distance between particles in an aggregate becomes larger than the coagulation distance, then the pair is disaggregated.

2.2 Rattle algorithm adapted for use with Langevin impulse integrator

The Rattle algorithm (Andersen, 1983) can be used to optimize the distances between neighboring particles in an aggregate to prevent them from overlapping. In this work the algorithm is adapted for use in combination with Langevin impulse integrator. The Rattle algorithm consists of the following steps:

• Step 1: Calculate half velocity values $v_{n+\frac{1}{2}}$ for all particles with no bonds taken into account. For each i-particle:

  (12)

  $$v_{i\text{,}n+\frac{1}{2}}=v_{i\text{,}n}+\frac{1}{2}{a}_{i\text{,}n}\mathrm{\Delta }t$$ where a_i,n is particle acceleration at step n, Δt is a time gap between steps n and n + 1, n + ½ is an intermediate step between steps n and n + 1.

• Step 2: Execute the correction of half velocities $v_{n+\frac{1}{2}}$ of all collided particles in clusters.

• Step 2.1: Select an unchecked pair of particles i and j. Calculate the desired distance d_ij between the centers of droplets:

  (13)

  $$d_{\mathit{ij}}=r_i+r_j$$

• Step 2.2: Calculate the relative shift between particles in pair:

  (14)

  $$s_n=(x_{i\text{,}n}+\mathrm{\Delta }{x}_{i\text{,}n}^{\prime })-(x_{j\text{,}n}+\mathrm{\Delta }{x}_{j\text{,}n}^{\prime })$$ where x_n is a particle position, Δx′ is a particle shift that is calculated by the eq. (2) using the adjusted value of $v_{n+\frac{1}{2}}$ to simulate the expected displacement of the particle at the next step.

• Step 2.3: If $\left|s_n^2-d_{\mathit{ij}}^2\right|⩽{\delta }_x$ then go to Step 2.5, otherwise calculate the half velocity correction coefficient:

  (15)

  $${\beta }_n=\frac{s_n^2-d_{\mathit{ij}\text{,}n}^2}{2s_n(x_{i\text{,}n}-x_{j\text{,}n})\left(\frac{1-e^{-\frac{{\zeta }_i}{m_i}\mathrm{\Delta }t}}{{\zeta }_i}+\frac{1-e^{-\frac{{\zeta }_j}{m_j}\mathrm{\Delta }t}}{{\zeta }_j}\right)}$$

• Step 2.4: Calculate a new value of $v_{n+\frac{1}{2}}$ :

  (16)

  $$v_{i\text{,}n+\frac{1}{2}}=v_{i\text{,}n+\frac{1}{2}}-{\beta }_n·\frac{x_{i\text{,}n}-x_{j\text{,}n}}{m_i}$$

  (17)

  $$v_{j\text{,}n+\frac{1}{2}}=v_{j\text{,}n+\frac{1}{2}}+{\beta }_n·\frac{x_{i\text{,}n}-x_{j\text{,}n}}{m_j}$$

• Step 2.5: If relative particle displacements in all particle pairs are within some tolerance or the upper limit of iterations has been reached, then go to Step 3. Else go to Step 2.1.

• Step 3: Calculate next particle coordinates for each i-th particle:

  (18)

  $$x_{i\text{,}n+1}=x_{i\text{,}n}+\mathrm{\Delta }{x}_{i\text{,}n}$$ using Eqs. (2), (8), (9) and (10).

• Step 4: Calculate accelerations a_n+1 of all particles for the positions x_n+1 using the sum over all interdroplet forces F.

• Step 5: Calculate particle velocities at the n + 1 step for each i-th particle using Eqs. (3) and (7).

• Step 6: Execute the velocity v_n+1 correction loop for all bonds between particles.

• Step 6.1: Select an unchecked pair of particles i and j. Calculate the desired distance d_ij between the centers of particles by the Eq. (13).

• Step 6.2: If $(v_{i\text{,}n+1}-v_{j\text{,}n+1})(x_{i\text{,}n+1}-x_{j\text{,}n+1})⩽{\delta }_v$ then go to Step 6.4, otherwise find the correction coefficient:

  (19)

  $${\sigma }_{n+1}=\frac{(v_{i\text{,}n+1}-v_{j\text{,}n+1})(x_{i\text{,}n+1}-x_{j\text{,}n+1})}{d_{\mathit{ij}}^2\left(\frac{1}{m_i}+\frac{1}{m_j}\right)}$$

• Step 6.3: Calculate new velocity values v_n+1:

  (20)

  $$v_{i\text{,}n+1}=v_{i\text{,}n+1}-{\sigma }_{n+1}·\frac{x_{i\text{,}n+1}-x_{j\text{,}n+1}}{m_i}$$

  (21)

  $$v_{j\text{,}n+1}=v_{j\text{,}n+1}+{\sigma }_{n+1}·\frac{x_{i\text{,}n+1}-x_{j\text{,}n+1}}{m_j}$$

• Step 6.4: If relative particle displacements in all particle pairs are within some tolerance or the upper limit of iterations has been reached, then go to Step 7. Else go to Step 6.1.

• Step 7: End of the algorithm. For all particles, the values of coordinate x_n+1, velocity v_n+1 and acceleration a_n+1 are known.

2.3 Implementation details and emulsion data

A finite number of droplets (N) are randomly distributed in a cubic box of 10 μm in edge length L. The volume fraction of aqueous droplets (ϕ) ranges from 0.05 to 0.50. The distance between just generated droplets is ⩾20 nm at ϕ < 0.50 and ⩾10 nm at ϕ = 0.50.

The selection of the time step is dictated by the numerical stability of the simulated system during the calculation, that is why a variable step Δt is used to enhance the stability. In every step each droplet can displace on a distance no larger than 25 nm that being the radius of the smallest droplet r_min in simulating emulsion. If any displacement exceeds this distance, Δt is halved.

In this study water-in-oil emulsions, i.e., emulsions wherein aqueous droplets are dispersed in the oil phase, have been modeled. The sizes of aqueous droplets in a typical W/O emulsion measured in the laboratory (Koroleva et al., 2011) are taken as the droplet size distribution in the simulations. The droplet radii range from 25 to 400 nm, and the average droplet radius is 160 nm.

It is assumed that the percolation cluster has sizes of no fewer than L - r_min along all three axes.

The viscosity and relative permittivity of the continuous phase (mineral oil) are equal to 0.1 Pa s and 2.43, respectively. The Debye length is considered equal to 400 nm similar to the 1/κ value for W/O emulsions stabilized with sorbitan monooleate. In W/O emulsions the surface potential slowly decays in the continuous media, so ζ-potential is used instead of ψ₀. The ζ-potential of aqueous droplets stabilized with sorbitan monooleate was experimentally determined by microelectrophoresis and is given in calculations as −10 mV. The temperature was assumed to be equal to 298 K.

3.1 Aggregation before the percolation clustering

The fraction of aqueous phase has a great effect on the structure of aggregates in W/O emulsions. The progress of the aggregation process with respect to change in cluster number N, relative to the initial number of droplets N₀, is illustrated in Fig. 1. Induction period for aggregation decreases sharply as aqueous phase concentration increases, and at ϕ = 0.5 it is almost absent.

Close

Figure 1

Time variations in cluster numbers N relative to the initial number of droplets N₀ at different aqueous phase fractions. Bold dots on the dashed line correspond to the formation of percolation clusters.

Export to PPT

As expected, droplets aggregate faster in more concentrated emulsions than in more dilute ones. Moreover, the aggregation rate increases with time before the percolation threshold is reached. The aggregation rate constant K enhances linearly on the log-log plot (Fig. 2).

Close

Figure 2

Log-log plot for the time evolution of rate constants upon aggregation in simulated emulsions with different aqueous phase fractions.

Export to PPT

At low fractions of aqueous phase, aggregates are located at a large distance from each other. Just formed aggregates transform into subrounded ones because attractive force exceeds repulsive force in apolar media. In more concentrated emulsions chain clusters aggregates and highly branched structures grow. As a result the aggregation rate constant increases over time.

Analysis of changes of cluster sizes over time was carried out. Distributions of i-fold clusters are illustrated in Fig. 3. For an emulsion with ϕ = 0.3, the time of percolation clustering is ∼0.2, and these histograms correspond to the percolating time and half-time – t_per and t_per/2. Number of singlets gradually decreases, but they do not disappear completely and remain in the system in the moment of percolation cluster formation and even later. Over the course of time, the number of dimers, trimers and other small clusters also diminishes and the number of n-fold clusters with large n increases. There is no maximum in the histograms that indicates aggregation mainly as a result of Brownian diffusion.

Close

Figure 3

Cluster number distribution for different times of aggregation. In all cases ϕ = 0.3 and N₀ = 9450.

Export to PPT

Evolutions in the concentrations of singlets, dimers, and trimers over time relative to coagulation half-time θ are presented in Fig. 4. Curves for simulated number of single droplets and dimers differ only slightly from data calculated by the Smoluchowski equation. Small deviations are attributed to randomness in the simulation process.

Close

Figure 4

Time evolution of singlets, dimers, and trimers in modeling emulsions with different droplet concentrations and theoretical distribution obtained from Smoluchowski equation.

Export to PPT

Distinct differences are observed in theoretical and simulated results for trimers and larger aggregates. The positions of simulated maxima for trimers (Fig. 4) and n-fold aggregates with n > 3 are near to theoretical maxima calculated using the Smoluchowski equation, but the fractions of the aggregates are smaller. Droplet concentrations in the simulated emulsions are high enough. Thus the probability of simultaneous collision of three droplets or aggregates is increased. Collisions between aggregates with different sizes take place, and all these lead to a more rapid diminishing of the number trimers and larger aggregates.

3.2 Percolation clustering

Droplet aggregation in emulsions leads to formation of globular or branched agglomerates. If droplet fraction emulsion exceeds a certain value branched agglomerates can combine and create a percolation cluster. The time required to reach the percolation threshold increases with decrease in aqueous phase fraction, as shown in Fig. 5.

Close

Figure 5

The time of percolation cluster formation as a function of aqueous phase fraction.

Export to PPT

At high aqueous phase fractions ϕ > 0.3 percolation clustering lasts for fractions of a second. The velocity of this process decreases with decreasing ϕ, and the time required for the percolation cluster formation increases up to 1.5 s at ϕ = 0.15 and rather higher at ϕ = 0.12 – up to 8.5 s.

In emulsions with aqueous phase fractions less than 0.12, the percolation threshold is not reached for several dozens of seconds. This can be ascribed to the fact that less branched aggregates exist in low-concentrated emulsions. Due to locations at large distances the aggregate structure is altered with time and transformation of small aggregates into subrounded ones takes place faster than their aggregation.

Fig. 6 illustrates the number of droplets in the percolation cluster and the fraction of the total amount of droplets in this cluster for emulsions with different water phase concentration. The quantity of droplets decreases with diminishing the water phase concentration, but the fraction of droplets combined into a percolation cluster is almost constant in the range of ϕ values from 0.35 to 0.50. Large branched aggregates are formed at these simulation conditions and percolation cluster arises by combining such aggregates.

Close

Figure 6

The effect of aqueous phase fraction in simulated emulsion on the number and the fraction of droplets in percolation cluster.

Export to PPT

Percolation clusters contained almost equal amounts of droplets in emulsions with ϕ spanning the range of 0.12–0.30. Relatively large aggregate is formed in such systems primarily. Then smaller aggregates and single particles attach to this large cluster until percolation clustering is reached. The fraction of droplets in a percolation cluster sharply increased.

The reason for the increase in the number of droplets in the percolation cluster at ϕ = 0.12 compared to ϕ = 0.15 is that less branched aggregates grow in the modeling cell. And at ϕ < 0.12 the percolation cluster is not formed at all.

The results of mathematical modeling enable us to conclude that percolation cluster is rapidly formed in simulated emulsions with ϕ ⩾ 0.3. In less concentrated emulsions percolation clustering is slower and the floc network involves a fewer droplets. In emulsions with aqueous phase fraction of 0.1 the percolation cluster is not formed for several dozens of seconds.

3.3 Comparison of simulated results and experimental W/O emulsion stability

Concentrated W/O emulsions are unstable toward flocculation. Dispersed phase droplets begin to aggregate almost immediately after dispersing is stopped in the process of emulsion preparation. Droplet aggregates are formed with different sizes and structures depending on the composition of a disperse system and the conditions.

Micrograph of W/O emulsion with a dispersed phase fraction of 0.3 is presented in Fig. 7. Droplet aggregation in investigated emulsions leads to the formation of flocs with a network structure. The structure of aggregates calculated according to Langevin dynamics model (Fig. 7) is similar to the structure of flocs composed of aqueous droplets in the laboratory prepared emulsion.

Close

Figure 7

Microphotograph of flocs in experimental W/O emulsion (left) and the snapshot of simulated aggregate (right) at ϕ = 0.3.

Export to PPT

The shapes and sizes of the formed flocs in emulsions noticeably affect the rate of sedimentation and, consequently, the rate of organic phase separation. The time before the visible start of the oil phase separation is almost equal to zero for emulsions with low aqueous phase fraction (Table 1). If ϕ exceeds 0.2 then the time of oil phase non-separation increases drastically. In the investigated emulsions the floc network is formed at the aqueous phase fractions above 0.25–0.30 (Koroleva et al., 2011). The higher the fraction of the dispersed phase is in the emulsion, the less time is necessary for the percolation cluster formation. The network of water droplets occupies almost entire volume of the emulsion and as a result sedimentation stability of emulsions increases greatly. According to the results of mathematical simulations and experiments, examined emulsions can be divided into the three following groups depending on their sedimentation stability (Table 2):

Overall, it can be concluded that simulation data correspond to experimental results on W/O emulsion stability. Due to formation of the percolation cluster emulsion stability significantly increases.