Phase Behavior of Alkyl Ethoxylate Surfactants in a Dissipative Particle Dynamics Model

We present a dissipative particle dynamics (DPD) model capable of capturing the liquid state phase behavior of nonionic surfactants from the alkyl ethoxylate (CnEm) family. The model is based upon our recent work [Anderson et al. J. Chem. Phys.2017, 147, 094503] but adopts tighter control of the molecular structure by setting the bond angles with guidance from molecular dynamics simulations. Changes to the geometry of the surfactants were shown to have little effect on the predicted micelle properties of sampled surfactants, or the water–octanol partition coefficients of small molecules, when compared to the original work. With these modifications the model is capable of reproducing the binary water–surfactant phase behavior of nine surfactants (C8E4, C8E5, C8E6, C10E4, C10E6, C10E8, C12E6, C12E8, and C12E12) with a good degree of accuracy.


Bead Volumes Corresponding to the DPD Model
presents the bead volumes for the chemical species studied in the main article. Here we report on work undertaken to understand the limitations of the model from the main text in this respect, focussing on C 12 E 12 . We explore multiple options relating to the modification of our model to begin to explore what factors lead to surfactant freezing.
Analysis is carried out using protocols developed in our previous work for freezing (waxing) of alkanes. 2 These protocols involve monitoring a nematic order parameter and mean squared displacement (MSD) of the molecules over a time interval ∆t = 500 DPD time units, for identifying the crystalline phase (see Bray et al. 2 for details). We additionally assessed the extent of segregation which occurred between the alkyl and ethoxylate groups which appears to be a necessity for locking of molecules into the solid. We explore four options for the bonded and non-bonded interactions of C 12 E 12 : The parameter set detailed in the main text as the revised model.

(Rigid)
Bond stiffness increased to k b = 5000 k B T and corresponding r 0 increased to maintain constant equilibrium length (as prescribed by Bray et al. 2 ).

(Straight)
All equilibrium angles set as θ 0 = 180°to straighten out the ethoxylate chain.    Table S2.  Figures S2(e-i)). Straightening the molecules results in high alignment with neighbouring molecules (Figures S2(c-d,g-h)). Both these two requirements seem to be necessary to cause the molecular diffusion to drop and the crystal to form (which occurred within 50,000 DPD time units, see Figures S2(g,h)). Bond rigidity appears to be of secondary importance and leads to more stable conformations with lower molecular diffusion than the more flexible case (Figures S2(b,d,f,h)).

S5
In this section we present a catalog of the typical phase morphologies seen in simulation for the different phases encountered in the main article. Identifying phases from static images can often be challenging so we provide rough descriptions and comments where appropriate.
We hope the typical images presented here will aid future work.

S6
Image Phase Description Comments L 1 (low concentration) Typical presentation is multiple aggregates of low aspect ratio.
L 1 also corresponds to high concentration large disordered aggregates (see below).
Many aggregates of the same size aligned along a particular direction.
May not be completely aligned, look for order in the simulation box.
H 1 A series of tubes/rods aligned along a particular direction.
There can often be bridging between rods and some defects in the structure.
V 1 / T A space filling phase which presents similarly to the high concentration L 1 phase, sometimes displaying order.
Usually presents with increased segregation in the eigenvalues of the isosurface normals than the L 1 phase.
L α Easy to identify as regular ordered slabs in the simulation cell.
May present as bridged or perforated lamellar.
Isotropic space filling phase.
Isosurface normal eigenvalues will all be broadly equal.

S7
Here we present the phase diagrams of C 12 E 6 and C 8 E 4 as determined by the DPD model for non-ionic surfactants developed by Johnston and co-workers (Figures S4 and S5). 4 In this model, the C n E m molecules are governed by only two bead types (CH 2 CH 2 & OCH 2 CH 2 ) versus the four adopted by the model presented in the main article. Specifically, the Johnston model does not incorporate the terminal polyethylene oxide (PEO) chain functionality of CH 2 OH and their specification of the PEO chain beads are subtly different, preferring to adopt beads of OCH 2 CH 2 versus the CH 2 OCH 2 we adopt in our work. This difference facilitates our model incorporating the CH 2 OH and CH 3 terminal beads as presented in the main article. Parameters defining the Johnston model can be found in their original article. Figure S4: Eigenvalues of the second moment of the isosurface normal distribution, as a function of concentration, for C 12 E 6 , determined via simulations using the Johnston et al. 4 DPD model. Lines are µ 3 (green -top), µ 2 (red -middle) and µ 1 (black -bottom). Figure S5: Eigenvalues of the second moment of the isosurface normal distribution, as a function of concentration, for C 8 E 4 , determined via simulations using the Johnston et al. 4 DPD model. Lines are µ 3 (green -top), µ 2 (red -middle) and µ 1 (black -bottom).

S8
Simulations were set up at 10 wt% spacing and, as in the main text, we sample over 20 million steps whilst applying the same phase identification criteria as we have throughout.
The Johnston model results in phase behavior well aligned to the experimental phase diagrams of the sampled surfactants. Here the Johnston model yields a single L 1 phase across the concentration ranged scanned for C 8 E 4 and manages to reproduce the L 1 , H 1 , and L α phases of C 12 E 6 . For C 12 E 6 the Johnston model presents a transitional region at 60 wt% that appears as a perforated lamellar phase when visualising the outputs of simulation. The H 1 phase appears to be slightly narrower than the experimental phase behavior although it would be prudent to sample the phase behavior with more simulated points to be certain.
The L α phase extends to a too high concentration where experimentally this phase would transition to L 1 at approx 85 wt%.

S9
Here we present the phase diagrams of C 12 E 6 and C 8 E 4 as determined by the DPD model for non-ionic surfactants developed by Lavagnini and co-workers ( Figures S6 and S7). 5  Here the Lavagnini model yields a single L 1 phase across the concentration ranged scanned for C 8 E 4 and manages to reproduce the L 1 , H 1 , and L α phases of C 12 E 6 . Figure S6: Eigenvalues of the second moment of the isosurface normal distribution, as a function of concentration, for C 12 E 6 , determined via simulations using the Lavagnini et al. 5 DPD model. Lines are µ 3 (green -top), µ 2 (red -middle) and µ 1 (black -bottom). S10 Figure S7: Eigenvalues of the second moment of the isosurface normal distribution, as a function of concentration, for C 8 E 4 , determined via simulations using the Lavagnini et al. 5 DPD model. Lines are µ 3 (green -top), µ 2 (red -middle) and µ 1 (black -bottom).
The onset of the H 1 phase of C 12 E 6 appears to occur at a higher concentration than reported experimentally and is narrower in it's stability in terms of concentration range.
The onset of the L α phase is in line with that expected from experimental studies and transitions to a re-entrant L 1 phase between 80 and 90 wt% as expected. For C 12 E 6 the Lavagnini model may present a transitional region at approx 65 wt%, however, we did not sample at this point in our brief screen of phase behavior of this surfactant. S11 6 Phase Diagrams of C 8 E 4 and C 12 E 6 According to the Original Model of Anderson et al. 3 Here we present the plots of concentration versus the eigenvalues of the orientational order parameter M for the C 12 E 6 and C 8 E 4 surfactants. These images support the data presented in Table 5 of the main article and the discussion therein. Figure S8: Eigenvalues of the second moment of the isosurface normal distribution, as a function of concentration, for C 12 E 6 , determined via simulations using the Anderson et al. 3 DPD model. Lines are µ 3 (green -top), µ 2 (red -middle) and µ 1 (black -bottom). S12 Figure S9: Eigenvalues of the second moment of the isosurface normal distribution, as a function of concentration, for C 8 E 4 , determined via simulations using the Anderson et al. 3 DPD model. Lines are µ 3 (green -top), µ 2 (red -middle) and µ 1 (black -bottom).