pubs.acs.

org/Langmuir Article

Tween-80 on Water/Oil Interface: Structure and Interfacial Tension

by Molecular Dynamics Simulations
Arthur Mussi Luz, Gabriel Barbosa, Carla Manske, and Frederico Wanderley Tavares*
Cite This: Langmuir 2023, 39, 3255−3265 Read Online

ACCESS Metrics & More Article Recommendations *

sı Supporting Information
See https://pubs.acs.org/sharingguidelines for options on how to legitimately share published articles.

ABSTRACT: Surfactants are used in many fields of the chemical

industry in a wide range of applications. Generally, surfactants on
water/oil interfaces reduce interfacial tension, enabling the
Downloaded via OXITENO S.A. - MATRIZ on July 14, 2023 at 14:52:31 (UTC).

formation of emulsions or providing greater stability to the

emulsion formed. Although molecular dynamics has been
extensively used and has achieved remarkable success in describing
thermodynamic and molecular properties of systems with surface-
active compounds, the traditional molecular simulation force fields
considerably constrain the system size and the time scale of
simulations. Here, we propose a coarse-grained model of
polysorbate 80, a nonionic surfactant commercially known as
Tween-80. Based on the proposed coarse-grained model, we
evaluate the influence of the more internal ethylene oxide chain on the properties of the water/Tween-80/decane system. We verify
with the simulation results that the model can reproduce the expected decrease in interfacial tension as the surfactant quantity
increases in the simulations. Furthermore, we observe changes in the surfactant orientation as their quantity in the interface
increases, indicating a preferential orientation for these molecules in the adsorption layer. We also assessed the partition coefficient
of the surfactant between the two bulk phases by performing free energy calculations, which showed a higher affinity of Tween-80
surfactants with the water phase.

■ INTRODUCTION
Exploration of oil reservoirs is the main factor in the
dynamics simulation can provide insights into surfactant
behavior, such as surfactant self-assembly and aggregation at
development of the chemical industry. Due to this, optimizing interfaces.18 Even though an atomistic representation can
the processes of exploration and refining of petroleum-based provide structural details, the computational cost associated
products is a constant goal of researchers in this field. In this with this description is much higher than those of coarse-
context, surfactant injection is a well-established enhanced oil grained models. Thus, coarse-graining allows for longer
recovery technique;1−5 however, the formation of stable crude- simulation times and bigger box sizes, making it possible to
oil emulsions during oil production and processing is investigate micellization behavior, emulsion formation, geom-
associated with operational problems mainly related to flow etry transitions, and other processes whose dynamics are
assurance. Understanding the surfactant’s action in emulsion usually in the microsecond scale.19−24 All-atoms and coarse-
stabilization is essential to optimizing the process throughout grained approaches can even be combined to provide faster
the oil production chain.6,7 simulations and more accurate results, such as done by
On the other hand, surfactants are widely used to stabilize Peroukidis et al.25 who obtained transport properties for
emulsions due to their capacity to reduce the water/oil micellar solutions by combining coarse-graining and reverse
interfacial tension.8 Specifically, polysorbate 80 (Tween-80) is all-atoms mapping.
a nonionic surfactant with several applications: bioremediation Regarding molecular modeling of Tween-80, Amani et al.26
of soils,9,10 water desalinization,11 nanoparticles stabilization modeled a micelle of this surfactant in water using the Martini
with antibacterial purposes,12 and hydrate formation.13 Even Force Field.27 Results for the radius of gyration and diffusion
so, the molecular structure of this surfactant is not univocal
since ethylene oxide (EO) groups are randomly distributed in
the multiple sites available for ethoxylation during its Received: November 2, 2022
synthesis.14 Revised: February 16, 2023
Molecular simulation is a powerful tool to investigate the Published: February 24, 2023
correlation between different structures and their effect on the
properties of the system since we can isolate one single
structure for the simulations.15−17 In particular, molecular

© 2023 American Chemical Society https://doi.org/10.1021/acs.langmuir.2c03001

3255 Langmuir 2023, 39, 3255−3265
Langmuir pubs.acs.org/Langmuir Article

Figure 1. (a) Structure of Tween-80 molecules with x, y, and z indicating the sites of ethoxylation; (b) Martini representation of the 3-H structure
(no internal ethylene oxide chain) with the selected beads; (c) Martini representation of the 4-H structure (equally distributed ethoxylated units)
with the selected beads.

coefficients agreed with the experimental data. Tang et al.28 beads were parametrized to represent different functional groups, with
used a united-atoms force field to improve efficiency and a total of 16 beads separated into four main categories: apolar (C),
represent different surfactant structures, analyzing micellar nonpolar (N), polar (P), and charged (Q). The nonbonded
solutions and air/water and oil/water interfaces. They verified parameters of each bead are obtained by fitting the hydration and
transfer free energies of several base compounds, allowing for
they could reproduce radial distribution functions for micelles transferability and accuracy for obtaining thermodynamic proper-
preassembled in a fixed aggregation number and the effect of ties.27 Papavasileiou et al.32 compared the efficiency of this approach
the number of tails in both surface and interfacial tension. for n-alkane mixtures to several all-atoms force fields, showing good
Posocco et al.20 used DPD techniques to verify the synergistic predictions for viscosity, liquid density, and surface tension with
interface stabilizing effects of Tween-80 and Span 20 due to coarse-grained models. For polarizable water and n-decane bulk
their complementary structures. For the same surfactants phases, we used the standard parameters provided by Marrink et al.27
(Tween-80 and Span 20), Liang et al.29 evaluated the stability and Yesylevskyy et al.31 The polarizable water model is better suited
of oil-in-water emulsions using DPD models. for this type of investigation due to the inclusion of partial charges,
In this work, we develop a model for Tween-80 molecules causing a better reproduction of the dielectric constant difference
based on the Martini 2 Force Field, evaluating the effect of between bulk phases. We chose not to include partial charges for n-
decane and Tween-80 in conformation with the parametrization
increasing surfactant quantities in water/Tween-80/decane scheme for the beads in the Martini Force Field.
interfacial systems. We also assess structural changes in the Here, two different distributions of ethylene oxide (EO) units in
surfactant conformation and provide insights into the Tween-80 (Figure 1a) structures were chosen: one with equally
interfacial organization of these structures. We evaluate two distributed units (x = y = z = w = 5, 4-H structure) and the other with
different molecular configurations for these surfactants, one no internal chains (w = 0, x = 6, y = z = 7, 3-H structure). We assign
with the ethylene oxide units equally split between the possible beads through the guidelines provided by the Martini Force Field. For
sites and the other with no internal ethylene oxide chain. This nonethylene oxide beads, apolar C1 and C4 beads were chosen for the
way, we can analyze the differences between the two hydrophobic tail, while nonpolar Na and SNa beads were chosen to
distributions to conclude if these variations are significant to model the oxygenated ether and ester groups in the core of the
sorbitan ring. We employed smaller beads to preserve the ring
the full description of these systems. Furthermore, the
structure of the sorbitan portion of the molecule. Ethylene oxide units
Weighted Histogram Analysis Method (WHAM)30 was used were modeled with customized beads using bonded and nonbonded
to obtain the potential of mean force of a single surfactant parameters from Grunewald et al.33 and Rossi et al.34 In a
model molecule as it moves from the bulk phases to the transferrable way, these parameters were obtained for ethoxylated
interface. This property provides information regarding CnEm surfactants by reproducing the phase behavior of these
partition between bulk phases and adsorption-free energy to surfactants in water, in terms of both self-assembly and lamellar
the interface. phases. Notably, the model proposed by Rossi et al.34 successfully
described the thermodynamic and structural properties of ethoxylated

■ METHODOLOGY
System Setup. The Martini Force Field27,31 with polarizable
polymers. A restricted angle approach was also used to improve
numerical stability.35 The remaining bonded parameters between
other bead types were obtained through optimization in the SWARM-
water was used to model both the water and oil phase (n-decane). CG software, comparing bond, angle, and dihedral distributions with
The coarse-graining process is based on a building block approach, those obtained through all-atoms simulations.36 Following the
mapping preassigned beads to groups of up to 4 heavy atoms. These SWARM-CG approach, the proposed 3-H coarse-grained structure

3256 https://doi.org/10.1021/acs.langmuir.2c03001
Langmuir 2023, 39, 3255−3265
Langmuir pubs.acs.org/Langmuir Article

Figure 2. Parameterization scheme for the Tween-80 model. Bonded distance and angle distributions were first obtained by all-atoms simulations,
and coarse-grained bonded parameters were optimized to match the distribution of the mapped configuration.

Figure 3. Visual representation of the end-to-end distances and the defined vectors for the orientation profile for (a) 3-H structure and (b) 4-H
structure. The represented distances are the end-to-end distances within each group for the fully stretched chain.

was mapped into the atomistic trajectory generating distribution 10125 decane molecules as bulk phases, with two interfaces between
profiles for bonds, angles, and dihedrals. A direct Boltzmann inversion them. It is worth highlighting that, in the Martini Force Field
generated the initial guess and is further optimized through particle framework, each CG unit of water involves four molecules. Tween-80
swarm optimization. The model was then validated by comparing the molecules for both structures (4-H and 3-H) were added in the
bonds, angles, and dihedral distributions of the mapped all-atoms following quantities: 20, 50, 80, 110, 140, and 170 molecules per
simulations to the coarse-grained simulations and comparing interface. Initially, the volume of the simulation boxes was 15 × 15 ×
characteristic distances such as the radius of gyration and end-to- 45 nm3. All simulations were conducted using the GROMACS 2021.1
end distances within the tail and head groups. The obtained package.43,44 All runs were performed at 298 K using a velocity rescale
parameters were then used for the 4-H structure, only changing the thermostat.45 The equilibration scheme is as follows: first, all systems
disposition of the beads. had their energies minimized through the steepest descent
Figure 2 illustrates the methodology used to obtain the bonded minimization algorithm. Then, simulations of 100 ps with a 1 fs
parameters for the coarse-grained. First, all-atoms simulations of a time step were carried out in an isotropic Berendsen barostat46 set to
single surfactant structure in water were conducted and generated the 1 bar. After this, the time step was increased to 10 fs. Simulations ran
distance and angle distributions within the mapped structure. These through 5 ns with a semi-isotropic Berendsen barostat keeping the
profiles were then used to fit the bonded parameters through a series interfacial area constant and the perpendicular pressure set to 1 bar.
of coarse-grained simulations by the particle swarm optimization The time step was again increased, reaching 20 fs, and the system was
method. further equilibrated for 40 ns with the same settings as before. Finally,
In our all-atoms simulations, the system with a single Tween-80 production runs employed a Parrinello−Rahman42 semi-isotropic
molecule in water was modeled using the General Amber Force barostat for 400 ns. Lennard-Jones interactions were calculated with a
Field37 and the TIP3P water model.38 The partial charges of the Potential-shift scheme with radii of 0.9 and 1.2 nm. A Particle-Mesh-
surfactant structure were calculated using the AM1-BCC scheme Ewald method47 was used to treat electrostatic interactions beyond
provided in Antechamber.39 The simulations were carried out with a the 1.2 nm cutoff.
single Tween-80 molecule surrounded by 1668 water molecules at For the free energy calculations, 9000 water CG units and 2500
298 K and 1 bar, simulated using the Nosé−Hoover40,41 thermostat decane molecules were placed in a 10 × 10 × 20 nm3 simulation box.
and Parrinello−Rahman barostat.42 The cutoff for short-range The Umbrella Sampling method48 was used to sample several
interaction was set to 1.2 nm, and the Coulombic interactions configurations of the surfactant across the simulation box, with one
beyond this radius were obtained with the Particle-Mesh-Ewald single Tween-80 molecule in various initial positions separated by 0.2
method. The time step was 1 fs for all simulations, and the nm each. The center of mass of the decane phase was used as the
equilibration and production parts lasted 2 and 50 ns, respectively. point of reference as described by Engin et al.49 This way, we
The topology for the proposed model of both Tween-80 structures is obtained the relative free energy of solubility in the form of the
available in the Supporting Information. potential of mean force (PMF) profile through the WHAM method.
Based on the obtained coarse-grained (CG) model, we then The used spacing coupled with a spring constant of 1000 kJ/mol
analyzed the interfacial behavior of the studied surfactants. Simulation proved to be sufficient to obtain good overlap in the configurations
boxes were built with 60000 coarse-grained (CG) units of water31 and sampled in each simulation, which were conducted for 100 ns.

3257 https://doi.org/10.1021/acs.langmuir.2c03001
Langmuir 2023, 39, 3255−3265
Langmuir pubs.acs.org/Langmuir Article

Thermodynamic and Structural Properties. Density profiles distances to verify how well the coarse-grained model matches
along the direction perpendicular to the interface for all species in the the behavior of the all-atoms representation. These parameters
systems were obtained to investigate the relative position of bulk can be found in Table 1. We verify that both radius of gyration
phases and surfactant molecules. Besides, distance histograms
between the ends of the head and tail were calculated to obtain a
more detailed description of their conformational behavior. These Table 1. Average Values of Radius of Gyration and End-to-
distances are represented in Figure 3. Similar histograms were End Distances for the Mapped All-Atoms and Coarse-
obtained for the distance between the tail and each headgroup to Grained Modelsa
evaluate their relative position. This was done by obtaining the
3-H structure 4-H structure
distances between the centers of mass of each group.
Interfacial tension (γ) was obtained through the diagonal stress All-atoms Coarse- All-atoms Coarse-
tensor components using the Irving-Kirkwood50 formalism: reference grained reference grained
Rg 0.84 ± 0.10 0.70 ± 0.08 0.75 ± 0.03 0.69 ± 0.10
Lz Pxx + Pyy
= Pzz Tail 1.28 ± 0.35 1.25 ± 0.28 1.16 ± 0.42 1.12 ± 0.31
2 2 (1) Head 1 1.10 ± 0.24 1.05 ± 0.24 0.96 ± 0.21 0.90 ± 0.28
Head 2 1.17 ± 0.32 1.13 ± 0.29 1.03 ± 0.27 1.00 ± 0.30
in which Pii is the ii component of the pressure tensor and Lz is the
Head 3 1.15 ± 0.30 1.12 ± 0.30 1.07 ± 0.24 1.02 ± 0.33
box length in the z-dimension; notably, we divide the original
equation of Irving and Kirkwood50 by two since there are two Head 4 0.99 ± 0.30 1.02 ± 0.35
a
interfaces in each system. We further investigate the orientation All values are in nanometers.
profile along the z-axis by calculating the order parameter Sz to
evaluate the orientation of the head and tail groups of the surfactants
regarding the interface. This parameter takes into account the angle θ and end-to-end distances are in agreement for both models
between the vector defined by the first and last bead of each group with overlapping intervals. This result shows that our
and that perpendicular to the interface. Namely, parametrization, combined with the parameters for ethoxylated
units of surfactants obtained from Grunewald et al.,33 is well
1
Sz = 3 cos2 1 suited to reproduce the behavior of Tween-80. This is further
2 (2) verified by the analysis of the tensoactive behavior of the
From eq 2, one can notice that Sz approaches 1.0 when the vector is model surfactants in the next section.
perpendicular to the interface, −0.5 when it is parallel to the Density Profiles and Interfacial Tension. In the
interface, and 0 when there is no preferential orientation. The defined previous section, we found out explicitly that the CG model
vectors are represented in Figure 3. obtained here can accurately describe the structural properties
Free energy-related quantities were also obtained through the of Tween-80. Thus, based on this model, we now investigate
differences between regions of the PMF profile. The partition the structural and thermodynamic properties of the Tween-80
coefficient for the surfactant between water and decane can be
obtained by a free energy difference between the surfactant in water
surfactant in water/decane interfacial systems. Figure 4 shows
and decane bulk phases: the density profile for the interfacial systems with 20 and 170
molecules of 3-H and 4-H structures. These profiles have not
ÅÄÅ É
Gwater oil ÑÑÑÑ
co Å been manually symmetrized over the two interfaces; the
P o/w = w = expÅÅÅÅ ÑÑÑ
c ÅÅÇ RT ÑÑÖ (3) symmetry was naturally achieved over the course of simulation
averages. Results for all simulated quantities of surfactant can
in which Po/w is the oil/water partition coefficient of the surfactant, co be found in Figures S4 and S5. As observed by inspecting
and cw are its solubilities in n-heptane and water, and Gwater oil is the Figure 4, the simulation approach leads to reasonably
relative molar free energy of solubility of the surfactant between oil equilibrated systems, in which the interfacial regions and
and water phases. When the coefficient is large, the surfactant shows a bulk phases are well-described. Notably, the tensoactive
higher affinity for the oil phase; the inverse is valid for the water
phase. We also obtained standard adsorption free energies for the behavior of both surfactant structures is verified in simulation
surfactant going from an infinite dilute bulk phase to the interface by results once these molecules tend to stay in the interfacial
taking the free energy differences of the surfactant in bulk phases and region between water and decane bulk phases. As seen in
at the interface. This adsorption free energy is directly related to the Figure 4, an analysis of the density profiles confirms this
drop in interfacial tension at infinite dilution through Gibbs behavior, showing the relative position of surfactants and bulk
adsorption isotherm. A more detailed deduction of partition phases. For low interfacial occupations, we verify the formation
coefficients and adsorption free energy can be found in the of a surfactant monolayer on the two interfaces of the
Supporting Information. simulation box.

■ RESULTS AND DISCUSSION

Validation of the Coarse-Grained (CG) Model. The
Although the density for both bulk phases is close to the
experimental values, with 989 kg/m3 for water and 755 kg/m3
for n-decane, density distributions are affected by the presence
obtained distributions for bond length, angles, and dihedrals of surfactants. In lower surface occupations, the monolayer
are represented in Figures S1, S2, and S3. The coarse-grained causes a density peak in the water phase close to the interface,
model reproduced the bond lengths with good accuracy, suggesting that these molecules are more attracted to this
although some distributions were less sharp than those of the region due to preferential interactions with surfactant
all-atoms mapped representation. The angles and dihedral molecules. A less significant effect is observed for the n-decane
distributions are also in good agreement with all-atoms results. phase, whose density decreases near the interfacial region.
Some limitations are inherent to the coarse-graining process, As the number of surfactants in the box increases, we verify a
such as the bimodal distributions that cannot be reproduced higher occupation of the interfacial region, leading to an
due to the loss of all-atoms details. increase in the surfactant density in the region between bulk
In addition to these distributions, we obtained geometric phases. A secondary surfactant density peak appears at higher
parameters such as the radius of gyration and end-to-end coverages, suggesting the formation of a second surfactant
3258 https://doi.org/10.1021/acs.langmuir.2c03001
Langmuir 2023, 39, 3255−3265
Langmuir pubs.acs.org/Langmuir Article

Figure 4. Density profiles along the direction perpendicular to the interface for both structures at the minimum and maximum surface coverages
analyzed. All simulations were conducted at 298 K and 1 bar.

Figure 5. Snapshot of the trajectory of the 3-H structure with 170 surfactants positioned in each side of the box. Blue, water; green, head groups;
black, tail groups; orange, decane.

adsorption layer. Upon further visual inspection of the the interfacial region, with faster saturation for the 3-H
trajectory (Figure 5), we verify a significant penetration of structure.
surfactant molecules into the water bulk phase due to the We also obtained specific density profiles separating head
saturation of the interface. In these simulations, surfactant and tail groups within each model structure, represented in
molecules self-organize in the water bulk by forming micelle- Figures S6 and S7. We verify a single sharp density peak for
both groups in simulations with low surface coverage. As
like aggregates, with a tail core in the center and the head
expected, the tail groups are closer to the decane phase, while
groups facing outward. Interestingly, this is the case only for the head groups are facing toward the water phase. Since there
the 3-H structure, as the density profile shows a negligible are more head groups in the structure of Tween-80, we verify
concentration of 4-H molecules at the core of the water phase, that the density peak is considerably higher when compared to
while for 3-H structures the density is nonzero. This shows that of the tail group. The head peak splits into two for higher
that structural changes affect the packing of surfactants within interfacial occupations. This behavior is caused by the
3259 https://doi.org/10.1021/acs.langmuir.2c03001
Langmuir 2023, 39, 3255−3265
Langmuir pubs.acs.org/Langmuir Article

saturation of the interface, causing surfactant molecules to Table 2. Comparison between Total Number of Surfactants
adjust their conformations to maximize exposure. Out of the per Interface in the System and Average Number of
four possible head groups, the ones near the tail are closer to Surfactants at Each Interface for Both Analyzed Structures
the n-decane phase, while head groups in the sorbitan ring are
Average number of surfactants at the
often found closer to water. interface
In terms of interfacial tension, the simulations can reproduce
Total surfactants in system 3-H structure 4-H structure
the decrease in this property when the surfactant coverage is
increased. Figure 6 shows this behavior for both analyzed 20 19.976 20.000
50 49.677 49.951
80 78.272 80.00
110 100.92 109.46
140 105.83 115.57
170 111.31 129.89

interfacial tension with the same number of surfactants in the

system as in the 3-H.
Structural Analysis. Given that the structure of Tween-80
can vary significantly due to the different possible ethoxylate
distributions of ethylene oxide units in the head groups, we can
obtain a greater understanding of the structure of the
interfacial region through analysis of distances and angles
within the molecules or in regard to the interface. In this work,
we verify the different conformations adopted by these
molecules when considering the presence of an internal
chain of ethylene oxide units.
Figure 7 shows the histograms of end-to-end distance within
the tail and head groups. We verify two unimodal profiles for
Figure 6. Interfacial tension behavior as a function of the surfactant
concentration in the water bulk phase for both structures. Error bars
the structure with equal distribution of ethylene oxide units (4-
are the minimum and maximum values of the calculated interfacial H structure). For both tail and head groups, the most probable
tensions for three independent runs in each concentration. distance is lower than when the structure is stretched, as
presented in Figure 3, showing a folding behavior around the
center bead of the chains, represented as a bending of the
chain. This observation is more subtle in the case of tail
structures. Interestingly, the calculated interfacial tensions are groups, whose distance peak is around 1.62 nm compared to
consistent with those from the DPD simulations of Posocco et 2.0 nm for the stretched structure. However, the distribution is
al.20 and the GROMOS53A6 force field of Tang et al.,28 sharper for head groups, showing there are fewer possible
showing that the coarse-grained model based on the Martini configurations for these groups. There is no significant
approach agrees with other models in the literature. difference in the folding behavior between head groups,
Still in Figure 6, the different surfactant molecules have little demonstrating that their positioning within the molecular
effect on the interfacial tension values, as both structures (4-H structure is not a relevant factor for their internal organization.
and 3-H) have similar tendencies for this property. We observe For the structure without an internal headgroup (3-H
slightly higher interfacial tension values for the structure with structure), we verify the same behavior for the tail group seen
the internal headgroup (4-H structure) at lower concen- for the other surfactant. The headgroup closer to the tail,
trations, but the opposite behavior is verified at higher labeled head 1, presented a bimodal distribution. This profile
concentrations, in which the 4-H structure shows lower indicates two equally possible configurations that result from
interfacial tension values compared to the 3-H structure. bending the chain around the center of mass of the group. In
Notice that, although simulations with both structures were the case of heads 2 and 3, although a bimodal distribution is
conducted with the same quantity of surfactant molecules and also present, the higher ethoxylation number slightly shifted
initial conditions, there is a slight difference in terms of the second peak position to a higher value than head 1.
equilibrium interfacial occupation due to the effect of the Notably, in the case of the heads 2 and 3, the distributions
ethoxylate distribution. Saturation is reached faster for the 3-H points to the configuration in which the groups are more
structure, with little difference in number of molecules at the stretched. We also verify a shorter average end-to-end distance
interface between the last three simulations, while the 4-H for head 1 when compared to heads 2 and 3. However, this is
structure shows higher occupancy of the interface, suggesting expected, as there is one less ethylene oxide unit in this
that this structure allows for more dense packing in the headgroup.
monolayer formed. This is illustrated by the comparison There is no significant difference in the end-to-end distances
between surfactants in the system and surfactants at the when the surface coverage is increased, as represented in
interfacial region in Table 2. These values were obtained by Figures S8 and S9. In other words, interfacial packing does not
integrating the density profile for both primary and secondary influence the internal disposal of tail and head groups
peaks of surfactant density. This is a crucial difference between regardless of the structure analyzed. This result could be
the two structures, as the higher interfacial occupation is caused by the coarse-graining process, which causes several
translated to lower interfacial tension values, as observed by different atoms to be grouped, reducing their flexibility and the
the last simulation of the 4-H structure, which achieved lower number of possible configurations.
3260 https://doi.org/10.1021/acs.langmuir.2c03001
Langmuir 2023, 39, 3255−3265
Langmuir pubs.acs.org/Langmuir Article

Figure 7. End-to-end distance histograms within the same groups for each surfactant structure at the minimum and maximum surface coverages
analyzed. All simulations were conducted at 298 K and 1 bar. Dashed lines indicate fully stretched values, which can be found in Figure 3.

In addition to the information within the same group, we We also analyzed the orientation of surfactant molecules in
also obtained histograms for distances between the same regards to the direction perpendicular to the interface through
molecule’s tail and head groups. This way, we can evaluate order parameter profiles, as seen in Figure 9. For the 4-H
both the internal structure of chains within the same molecule structure, all groups show positive values of the order
and the relative position of these groups. Figure 8 shows the parameter, indicating that these molecules show preferential
tail−head distances for both structures. Regarding the 4-H orientations in the interfacial region. For lower occupancies,
structure, we verify a sharp peak at low distances, with an tail groups are the most tilted followed by the head groups with
average value of approximately 1.3 nm in the tail−head 4 less preferential behavior. At higher occupancies, all groups
distribution. This group is a direct continuation of the tail achieve higher values of order parameter: interface packing
chain, which explains the low variance observed. For heads 1, leads all groups to adopt tilted orientations to accommodate
2, and 3, we verify a broader distribution ranging from 0.5 to the high local concentration of surfactants. Specifically, the
3.5 nm, showing several possible configurations for the internal head chain adopts the tilted configuration, reaching
molecules at the interfacial region. These profiles become angles of about 58.9° with the interfacial plane. These profiles
narrower toward higher distances as we increase the surface also show two different orientations depending on the
coverage due to the packing effect at the interface. molecule’s position within the surfactant monolayer. Groups
The profiles for the 3-H structure show lower distances closer to the decane phase are more tilted in regards to the
overall compared to the 4-H structure due to the lack of the interface than those next to the water phase. We also verify
internal headgroup separating the tail and head groups. At low small order parameters for surfactants located farther from the
surface coverages, for head 1, the profile is sharp with a mean interface, corresponding to those present in aggregates in the
distance of 1.5 nm, while the distributions are broad for the water bulk phase. Since these aggregates suffer lower effects
head groups bonded to the sorbitan ring around the same from the interface, there is no preferential orientation.
value. These profiles also shift when surface coverage increases, We can make similar remarks regarding the orientation of
with a more significant change than that observed for the 4-H groups within the 3-H structure. At low surface coverages, all
structure. The head groups connected to the sorbitan ring groups show a single preferential orientation as their average
locate themselves farther from the tail group as an effect of the position regarding the interface stays the same throughout the
surfactant packing at the interface, shifting their average simulation, with less tilted configuration compared to high
configuration. surface coverages. In this case, we observe two possible
3261 https://doi.org/10.1021/acs.langmuir.2c03001
Langmuir 2023, 39, 3255−3265
Langmuir pubs.acs.org/Langmuir Article

Figure 8. Head−tail distance histograms for both structures at the minimum and maximum surface coverages analyzed. All simulations were
conducted at 298 K and 1 bar.

positions for the tail groups with different orientations within ethylene oxide units and a relatively small tail group in
the interfacial region. Their orientations are slightly more tilted comparison.
to maximize interactions with the decane phase when located We also verify that the interfacial region is the most
farther away from this bulk phase. However, there is a single favorable regarding surfactant molecule positioning, as
preferential orientation for head groups. This behavior is due illustrated by the deeper well of free energy when changing
to the lack of an internal headgroup, causing these groups to be
more affected by the interfacial region. This effect is also from one bulk phase to the other. This result also confirms the
verified through the order parameters at the bulk water phase: tensoactive behavior of these coarse-grained models of Tween-
aggregates formed with the 3-H structure are more oriented 80.
toward the interface than those for the 4-H structure. Finally, analysis of this profile shows no significant difference
Potential of Mean Force (PMF) and Partition. As in terms of free energy between the 4-H and 3-H structures. A
described previously, the PMF of a single surfactant molecule slight difference in the relative solubility is attributed to the
as a function of the distance between the center-of-mass of the inherent uncertainties of the coarse-grained force field; i.e.,
surfactant and decane phase was calculated to investigate the without population effects, both structures are equivalent and
relative affinity of the different structures throughout the
can be used to describe the solubility of Tween-80 in different
simulation box without population effects. Notably, these PMF
profiles are inherently related to the solubility of the surfactant phases. This result corroborates the trend in the literature to
model in both decane and water bulk phases and also in the use only a single representative structure when modeling the
interfacial region. surfactant’s properties.26,52 However, the structure of the
As shown in the PMF profile of Figure 10, there is a surfactant is relevant to the packing of these molecules in the
significant free energy difference between decane and aqueous interfacial region, as previously evidenced by the presence of 3-
phases. Both surfactant structures analyzed show a much H structures in the water bulk phase. Therefore, although both
higher solubility in water. The results in Table 3 reveal the
structures show the same free energy contributions when
magnitude of this solubility difference. Although this difference
is very high, the Martini Force Field tends to overestimate the isolated, the molecular structure, specifically the distribution of
free energy of solubility of this type of surfactant.51 Even so, a ethoxylated units, contributes to the efficiency in reducing
high difference in solubility between bulk phases was already interfacial tension by allowing a higher number of possible
expected, as the structure of Tween-80 contains many polar configurations at the interfacial region.
3262 https://doi.org/10.1021/acs.langmuir.2c03001
Langmuir 2023, 39, 3255−3265
Langmuir pubs.acs.org/Langmuir Article

Figure 9. Order parameter profiles for tail and heads groups for both structures at the minimum and maximum surface coverages analyzed. All
simulations were conducted at 298 K and 1 bar.

Table 3. Free Energy Differences (kcal/mol) of Transferring

the Surfactant between Bulk Phases and Interface and Oil/
Water Partition Coefficient of the Surfactant at 293 K and 1
Bar
Surfactant ΔG̅ oil−int ΔG̅ water−int ΔG̅ water−oil ln Po/w
3-H −45.16 −14.02 31.13 −52.54
4-H −45.01 −14.73 30.28 −51.11

analyzed two different structures to account for the possible

distributions of ethylene oxide units at the headgroup. We
verify that the developed models are consistent with the
interfacial tension behavior expected for these molecules,
showing an accumulation at the interface and reduction of the
interfacial tension of the system as the surface coverage of
surfactants increases. We also observed an aggregation
Figure 10. Potential of mean force of a single surfactant structure as a behavior of the 3-H structure when the interfacial monolayer
function of its position along the z-axis of the simulation box. The is saturated, causing molecules to migrate to the water bulk
distance is measured between the center of mass of the decane phase phase and form micelle-like aggregates, suggesting that the
and the surfactant. Thus, the decane bulk is at the left side while the
critical micelle concentration is within the range of surface
water bulk is at the right side. All simulations were conducted at 298
K and 1 bar. coverages analyzed. Both structures showed the same behavior
in terms of interfacial tension. However, the 4-H structure,
which contains the internal headgroup, achieved higher surface

■ CONCLUSION
Here, we obtain a coarse-grained model of Tween-80
coverages due to more efficient packing at the monolayer,
causing lower values of interfacial tension when compared to
the same amount of molecules as the 3-H structure.
surfactant based on the Martini 2 Force Field. Using this Regarding structural behavior at the interface, both
model, we investigate the surfactant behavior at the water/n- structures exhibit a rigid behavior for the end-to-end distances
decane interface through molecular dynamics simulations. We within head and tail groups. This result shows that the coarse-
3263 https://doi.org/10.1021/acs.langmuir.2c03001
Langmuir 2023, 39, 3255−3265
Langmuir pubs.acs.org/Langmuir Article

graining process reduces the degrees of freedom within chains

and causes similar folding behavior regardless of the number of
■ ACKNOWLEDGMENTS
The authors gratefully acknowledge Conselho Nacional de
surfactants. The presence of the internal head chain causes the
Desenvolvimento Científico e Tecnológico (CNPq), CAPES,
overall tail−head distances to increase, so the interface effect at
FAPERJ, and Oxiteno for financial support. We also acknowl-
the head groups bonded to the sorbitan ring is modified. This
edge NACAD-UFRJ and LNCC for the computational
way, more conformations can be accessed by the 4-H structure,
resources.
resulting in broader distance distributions. In terms of
orientation, the saturation of the interfacial monolayer induces
a tilted orientation to the interface of both tail and head groups
to accommodate the high density of surfactant molecules
■ REFERENCES
(1) Raffa, P.; Broekhuis, A. A.; Picchioni, F. Polymeric surfactants for
present. enhanced oil recovery: A review. J. Pet. Sci. Eng. 2016, 145, 723−733.
Finally, analysis of the solvation-free energy confirms the (2) Mandal, A. Chemical flood enhanced oil recovery: A review.
dependence of the interfacial tension with surfactant International Journal of Oil, Gas and Coal Technology 2015, 9, 241−
concentration. It also shows that both structures are much 264.
more soluble in water than in n-decane, as evidenced by the (3) Alvarado, V.; Manrique, E. Enhanced oil recovery: An update
review. Energies 2010, 3, 1529−1575.
high order of magnitude of the partition coefficient.
(4) Sheng, J. J. Review of Surfactant Enhanced Oil Recovery in
Furthermore, the PMF and the interfacial tension profiles Carbonate Reservoirs. Advances in Petroleum Exploration and Develop-
show no significant difference between both analyzed ment 2013, 6, 1−10.
structures in terms of interfacial free energy. These results (5) Zhou, Y.; Yin, D.; Chen, W.; Liu, B.; Zhang, X. A comprehensive
indicate that a single representative structure can be used to review of emulsion and its field application for enhanced oil recovery.
model this surfactant at infinite dilute concentrations. As the Energy Science and Engineering 2019, 7, 1046−1058.
number of surfactant molecules increases, the molecular (6) Sousa, A. M.; Pereira, M. J.; Matos, H. A. Oil-in-water and water-
structure affects packing at the surfactant monolayer, causing in-oil emulsions formation and demulsification. J. Pet. Sci. Eng. 2022,
the 4-H structure to be more efficient in reducing the 210, 110041.
interfacial tension of the system. (7) de Oliveira, M. C. K.; Miranda, L. R. O.; de Carvalho, A. B. M.;
Miranda, D. F. S. Viscosity of Water-in-Oil Emulsions from Different

■
*
ASSOCIATED CONTENT
sı Supporting Information
American Petroleum Institute Gravity Brazilian Crude Oils. Energy
Fuels 2018, 32, 2749−2759.
(8) Bhardwaj, A.; Hartland, S. Applications of Surfactants in
The Supporting Information is available free of charge at Petroleum Industry. J. Dispersion Sci. Technol. 1993, 14, 87−116.
https://pubs.acs.org/doi/10.1021/acs.langmuir.2c03001. (9) Cheng, M.; Zeng, G.; Huang, D.; Yang, C.; Lai, C.; Zhang, C.;
Liu, Y. Tween 80 surfactant-enhanced bioremediation: toward a
Theoretical background regarding the calculation of
solution to the soil contamination by hydrophobic organic
partition coefficients and relating adsorption-free en- compounds. Critical Reviews in Biotechnology 2018, 38, 17−30.
ergies to the change in interfacial tension; bond length, (10) Cheng, M.; Zeng, G.; Huang, D.; Yang, C.; Lai, C.; Zhang, C.;
angle, and dihedral distributions; density profiles; Liu, Y. Advantages and challenges of Tween 80 surfactant-enhanced
distance histograms; orientation profile; GROMACS technologies for the remediation of soils contaminated with
topology files for each surfactant structure (PDF) hydrophobic organic compounds. Chem. Eng. J. 2017, 314, 98−113.
(11) Bide, Y.; Fashapoyeh, M. A.; Shokrollahzadeh, S. Structural

■ AUTHOR INFORMATION
Corresponding Author
investigation and application of Tween 80-choline chloride self-
assemblies as osmotic agent for water desalination. Sci. Rep. 2021, 11,
1−11.
Frederico Wanderley Tavares − Escola de Química, (12) Bekhit, M.; Abu el naga, M. N.; Sokary, R.; Fahim, R. A.; El-
Universidade Federal do Rio de Janeiro, Rio de Janeiro, Rio Sawy, N. M. Radiation-induced synthesis of tween 80 stabilized silver
de Janeiro 21941-909, Brasil; Programa de Engenharia nanoparticles for antibacterial applications. Journal of Environmental
Química - PEQ/COPPE, Universidade Federal do Rio de Science and Health - Part A Toxic/Hazardous Substances and
Janeiro, Rio de Janeiro, Rio de Janeiro 21941-909, Brasil; Environmental Engineering 2020, 55, 1210−1217.
orcid.org/0000-0001-8108-1719; Email: tavares@ (13) Abolfazl Mohammadi. Rate of Methane Hydrate Formation in
Presence of Tetra-n-Butylammonium Chloride and Tween 80:
eq.ufrj.br
Application in Air Conditioning Systems. Theoretical Foundations of
Authors Chemical Engineering 2021, 55, 451−456.
(14) Prouzet, E.; Cot, F.; Nabias, G.; Larbot, A.; Kooyman, P.;
Arthur Mussi Luz − Escola de Química, Universidade Federal
Pinnavaia, T. J. Assembly of Mesoporous Silica Molecular Sieves
do Rio de Janeiro, Rio de Janeiro, Rio de Janeiro 21941-909, Based on Nonionic Ethoxylated Sorbitan Esters as Structure
Brasil Directors. Chem. Mater. 1999, 11, 1498−1503.
Gabriel Barbosa − Department of Chemical and Biological (15) Dong, R.; Bao, L.; Yu, Q.; Wu, Y.; Ma, Z.; Zhang, J.; Cai, M.;
Engineering, University of Alabama, Tuscaloosa, Alabama Zhou, F.; Liu, W. Effect of Electric Potential and Chain Length on
35487, United States Tribological Performances of Ionic Liquids as Additives for Aqueous
Carla Manske − Escola de Química, Universidade Federal do Systems and Molecular Dynamics Simulations. ACS Appl. Mater.
Rio de Janeiro, Rio de Janeiro, Rio de Janeiro 21941-909, Interfaces 2020, 12, 39910−39919.
Brasil (16) Zhao, F.; Wang, S.; Guo, J. Experimental and Molecular
Dynamics Simulation Study on the Effects of the Carbon Chain
Complete contact information is available at: Length of Gemini Surfactants on the Inhibition of the Acid-Rock
https://pubs.acs.org/10.1021/acs.langmuir.2c03001 Reaction Rate. LangmuiEmpirical molecular modelling of suspension
stabilisationr 2021, 37, 5232−5241.
Notes (17) Maillet, J. B.; Lachet, V.; Coveney, P. V. Large scale molecular
The authors declare no competing financial interest. dynamics simulation of self-assembly processes in short and long

3264 https://doi.org/10.1021/acs.langmuir.2c03001
Langmuir 2023, 39, 3255−3265
Langmuir pubs.acs.org/Langmuir Article

chain cationic surfactants. Phys. Chem. Chem. Phys. 1999, 1, 5277− for coarse-grained molecular dynamics simulations. J. Chem. Theory
5290. Comput. 2013, 9, 3282−3292.
(18) Taddese, T.; Anderson, R. L.; Bray, D. J.; Warren, P. B. Recent (36) Empereur-Mot, C.; Pesce, L.; Doni, G.; Bochicchio, D.; Capelli,
advances in particle-based simulation of surfactants. Curr. Opin. R.; Perego, C.; Pavan, G. M. Swarm-CG: Automatic Parametrization
Colloid Interface Sci. 2020, 48, 137−148. of Bonded Terms in MARTINI-Based Coarse-Grained Models of
(19) Sanders, S. A.; Panagiotopoulos, A. Z. Micellization behavior of Simple to Complex Molecules via Fuzzy Self-Tuning Particle Swarm
coarse grained surfactant models. J. Chem. Phys. 2010, 132, 114902. Optimization. ACS Omega 2020, 5, 32823−32843.
(20) Posocco, P.; Perazzo, A.; Preziosi, V.; Laurini, E.; Pricl, S.; (37) Wang, J.; Wolf, R. M.; Caldwell, J. W.; Kollman, P. A.; Case, D.
Guido, S. Interfacial tension of oil/water emulsions with mixed non- A. Development and testing of a general amber force field. J. Comput.
ionic surfactants: Comparison between experiments and molecular Chem. 2004, 25, 1157−1174.
simulations. RSC Adv. 2016, 6, 4723−4729. (38) Berendsen, H. J.; Grigera, J. R.; Straatsma, T. P. The missing
(21) Illa-Tuset, S.; Malaspina, D. C.; Faraudo, J. Coarse-grained term in effective pair potentials. J. Phys. Chem. 1987, 91, 6269−6271.
molecular dynamics simulation of the interface behaviour and self- (39) Case, D.; Darden, T.; Cheatham, T., III; Simmerling, C.; Wang,
assembly of CTAB cationic surfactants. Phys. Chem. Chem. Phys. 2018, J.; Duke, R.; Luo, R.; Walker, R.; Zhang, W.; Merz, K.; et al.
20, 26422−26430. AmberTools 16; University of California: San Francisco, 2016.
(22) Li, C.; Pu, H.; Zhang, S.; Zhao, J. Effect of nanoparticles and (40) Nosé, S. A molecular dynamics method for simulations in the
surfactants on oil/Water interfacial tension: A coarse-grained canonical ensemble. Molecular physics 1984, 52, 255−268.
molecular dynamics simulation study. In SPE/AAPG/SEG Unconven- (41) Hoover, W. G. Canonical dynamics: Equilibrium phase-space
tional Resources Technology Conference 2020, URTeC 2020 2020; pp distributions. Phys. Rev. A 1985, 31, 1695.
(42) Parrinello, M.; Rahman, A. Polymorphic transitions in single
1454−1464.
crystals: A new molecular dynamics method. J. Appl. Phys. 1981, 52,
(23) Sangwai, A. V.; Sureshkumar, R. Coarse-grained molecular
7182−7190.
dynamics simulations of the sphere to rod transition in surfactant
(43) Berendsen, H. J.; van der Spoel, D.; van Drunen, R.
micelles. Langmuir 2011, 27, 6628−6638.
GROMACS: A message-passing parallel molecular dynamics
(24) Chen, M.; Lu, X.; Liu, X.; Hou, Q.; Zhu, Y.; Zhou, H.
implementation. Comput. Phys. Commun. 1995, 91, 43−56.
Temperature-dependent phase transition and desorption free energy (44) Abraham, M. J.; Murtola, T.; Schulz, R.; Páll, S.; Smith, J. C.;
of sodium dodecyl sulfate at the water/vapor interface: Approaches Hess, B.; Lindah, E. Gromacs: High performance molecular
from Molecular dynamics simulations. Langmuir 2014, 30, 10600− simulations through multi-level parallelism from laptops to super-
10607. computers. SoftwareX 2015, 1−2, 19−25.
(25) Peroukidis, S. D.; Tsalikis, D. G.; Noro, M. G.; Stott, I. P.; (45) Bussi, G.; Donadio, D.; Parrinello, M. Canonical sampling
Mavrantzas, V. G. Quantitative Prediction of the Structure and through velocity rescaling. J. Chem. Phys. 2007, 126, 014101.
Viscosity of Aqueous Micellar Solutions of Ionic Surfactants: A (46) Berendsen, H. J.; Postma, J. P.; Van Gunsteren, W. F.; Dinola,
Combined Approach Based on Coarse-Grained MARTINI Simu- A.; Haak, J. R. Molecular dynamics with coupling to an external bath.
lations Followed by Reverse-Mapped All-Atom Molecular Dynamics J. Chem. Phys. 1984, 81, 3684−3690.
Simulations. J. Chem. Theory Comput. 2020, 16, 3363−3372. (47) Darden, T.; York, D.; Pedersen, L. Particle mesh Ewald: An N·
(26) Amani, A.; York, P.; De Waard, H.; Anwar, J. Molecular log(N) method for Ewald sums in large systems. J. Chem. Phys. 1993,
dynamics simulation of a polysorbate 80 micelle in water. Soft Matter 98, 10089−10092.
2011, 7, 2900−2908. (48) Torrie, G.M.; Valleau, J.P. Nonphysical sampling distributions
(27) Marrink, S. J.; Risselada, H. J.; Yefimov, S.; Tieleman, D. P.; De in Monte Carlo free-energy estimation: Umbrella sampling. J. Comput.
Vries, A. H. The MARTINI force field: Coarse grained model for Phys. 1977, 23, 187−199.
biomolecular simulations. J. Phys. Chem. B 2007, 111, 7812−7824. (49) Engin, O.; Villa, A.; Sayar, M.; Hess, B. Driving forces for
(28) Tang, X.; Huston, K. J.; Larson, R. G. Molecular Dynamics adsorption of amphiphilic peptides to the air-water interface. J. Phys.
Simulations of Structure-Property Relationships of Tween 80 Chem. B 2010, 114, 11093−11101.
Surfactants in Water and at Interfaces. J. Phys. Chem. B 2014, 118, (50) Irving, J. H.; Kirkwood, J. G. The statistical mechanical theory
12907−12918. of transport processes. IV. The equations of hydrodynamics. J. Chem.
(29) Liang, X.; Wu, J.; Yang, X.; Tu, Z.; Wang, Y. Investigation of Phys. 1950, 18, 817−829.
oil-in-water emulsion stability with relevant interfacial characteristics (51) Huston, K. J.; Larson, R. G. Reversible and Irreversible
simulated by dissipative particle dynamics. Colloids Surf., A 2018, 546, Adsorption Energetics of Poly(ethylene glycol) and Sorbitan
107−114. Poly(ethoxylate) at a Water/Alkane Interface. Langmuir 2015, 31,
(30) Kumar, S.; Bouzida, D.; Swendsen, R. H.; Kollman, P. A.; 7503−7511.
Rosenberg, J. M. The weighted histogram analysis method for free- (52) Konkel, J. T.; Myerson, A. S. Empirical molecular modelling of
energy calculations on biomolecules. Journal of computational suspension stabilisation with Polysorbate 80. Mol. Simul. 2008, 34,
chemistry 1992, 13, 1011−1021. 1353−1357.
(31) Yesylevskyy, S. O.; Schäfer, L. V.; Sengupta, D.; Marrink, S. J.
Polarizable water model for the coarse-grained MARTINI force field.
PLoS Computational Biology 2010, 6, e1000810.
(32) Papavasileiou, K. D.; Peristeras, L. D.; Bick, A.; Economou, I.
G. Molecular Dynamics Simulation of Pure n-Alkanes and Their
Mixtures at Elevated Temperatures Using Atomistic and Coarse-
Grained Force Fields. J. Phys. Chem. B 2019, 123, 6229−6243.
(33) Grunewald, F.; Rossi, G.; De Vries, A. H.; Marrink, S. J.;
Monticelli, L. Transferable MARTINI Model of Poly(ethylene
Oxide). J. Phys. Chem. B 2018, 122, 7436−7449.
(34) Rossi, G.; Fuchs, P. F.; Barnoud, J.; Monticelli, L. A coarse-
grained MARTINI model of polyethylene glycol and of polyoxy-
ethylene alkyl ether surfactants. J. Phys. Chem. B 2012, 116, 14353−
14362.
(35) Bulacu, M.; Goga, N.; Zhao, W.; Rossi, G.; Monticelli, L.;
Periole, X.; Tieleman, D. P.; Marrink, S. J. Improved angle potentials

3265 https://doi.org/10.1021/acs.langmuir.2c03001
Langmuir 2023, 39, 3255−3265