Molecular Dynamics Simulation of the Oil ... - ACS Publications

Subscriber access provided by UNIV OF CONNECTICUT

Article

Molecular Dynamics Simulation of the Oil Detachment Process within Silica Nanopores Hui Yan, and Shiling Yuan J. Phys. Chem. C, Just Accepted Manuscript • DOI: 10.1021/acs.jpcc.5b09841 • Publication Date (Web): 20 Jan 2016 Downloaded from http://pubs.acs.org on January 26, 2016

Just Accepted “Just Accepted” manuscripts have been peer-reviewed and accepted for publication. They are posted online prior to technical editing, formatting for publication and author proofing. The American Chemical Society provides “Just Accepted” as a free service to the research community to expedite the dissemination of scientific material as soon as possible after acceptance. “Just Accepted” manuscripts appear in full in PDF format accompanied by an HTML abstract. “Just Accepted” manuscripts have been fully peer reviewed, but should not be considered the official version of record. They are accessible to all readers and citable by the Digital Object Identifier (DOI®). “Just Accepted” is an optional service offered to authors. Therefore, the “Just Accepted” Web site may not include all articles that will be published in the journal. After a manuscript is technically edited and formatted, it will be removed from the “Just Accepted” Web site and published as an ASAP article. Note that technical editing may introduce minor changes to the manuscript text and/or graphics which could affect content, and all legal disclaimers and ethical guidelines that apply to the journal pertain. ACS cannot be held responsible for errors or consequences arising from the use of information contained in these “Just Accepted” manuscripts.

The Journal of Physical Chemistry C is published by the American Chemical Society. 1155 Sixteenth Street N.W., Washington, DC 20036 Published by American Chemical Society. Copyright © American Chemical Society. However, no copyright claim is made to original U.S. Government works, or works produced by employees of any Commonwealth realm Crown government in the course of their duties.

Page 1 of 30

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

Molecular Dynamics Simulation of the Oil Detachment Process within Silica Nanopores Hui Yan† and Shiling Yuan*‡ † School of Pharmacy, Liaocheng University, Liaocheng 252059, China ‡ Key Laboratory of Colloid and Interface Chemistry, Shandong University, Jinan 250100, China

ABSTRACT: In this work, we investigated the effect of surfactants on the oil displacement process inside a nanoscale silica pore using molecular dynamics simulations. Firstly, an oil cylinder was built inside a silica pore to mimic residual oil in the porosity of the reservoir rock after water flooding. In the simulations, we focused on a layering organization of oil molecules adsorbed onto the pore surface, and then a series of equilibrium MD simulations were run to obtain the organization structures of the oil drop in the presence or absence of surfactant molecules. These simulated results showed that the surfactant could disturb the layering organization of the oil drop, since the hydrophobic chains of surfactant molecules could penetrate into the oil phase. And around the polar head of the surfactant, water molecules could form water channels between the oil phase and solid surface, which is vital to the displacement process. Finally, we used steered molecular dynamics (SMD) method to mimic the displacement and migration process of an oil drop inside the pore. From SMD calculations, detailed information about the process was obtained and the free energy of the process was calculated using the WHAM method. Through analysis of the free energy, we demonstrated the mechanism of surfactants aiding in the oil recovery at a molecular level. Our study provided information on

ACS Paragon Plus Environment

1


1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

the oil displacement within a nanoscale pore at molecular-level, which is expected to provide useful information for EOR experiments.

1. INTRODUCTION After water flooding, trapped residual oil is distributed in the pores of the reservoir rock in the shape of static oil droplets due to capillary action, which makes the oil recovery process more difficult.1 In recent years, enhanced oil recovery (EOR) techniques have been widely used to increase crude oil production. Many EOR techniques are applied in the oil recovery process, such as chemical injection and gas injection.2 In these techniques, the surfactant is widely used in chemical injection, which makes it feasible to recover the oil. It is known that the surfactant molecules can be adsorbed onto the oil-fluid interface, which could reduce the oil/water interfacial tension (IFT) and cause wettability alteration on the fluid/solid interface.3,4 The wettability-altering behavior of the solid/fluid system is a major factor in the EOR process, which can be controlled by surfactants due to their surface activity. Many experimental studies were performed to understand the mechanism of surfactants on residual oil by micro models of oil displacement over the past few decades.5-11 In these experiments, a single cylindrical capillary is usually employed as the micro model to construct a simplified system for multiphase flow in porous media, which is convenient to mimic the micro situations in the oil recovery process. Tiberg et al.8 used an oil-filled capillary tube to study the imbibition behavior of a surfactant solution. They found that it was difficult for the surfactant solution to invade into the capillary tube even in the case of high surfactant concentrations.8 On the basis of their observations, Hammond and Unsal9 modified the model for forced imbibitions. Their results showed that under additional forced conditions, the meniscus would move faster. Later, they reported a further investigation with diffusion of surfactant solutions in front of the mobile oil/water


2

Page 2 of 30

Page 3 of 30

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


meniscus.10 They found that the adsorption of surfactant molecules onto the capillary surface is significant for the fluid flow in the capillary tubes. Recently, Wasan et al.11 investigated the displacement of hexadecane by sodium dodecyl sulfate (SDS) micellar solution using glass capillaries. They found that the cylindrical droplets were finally detached from the capillary surface, forming spherical drops within the capillary. Based on their observations, a theoretical model to assess the dewetting velocity for the immiscible liquid inside a cylindrical tube was developed. These experimental investigations help us to understand the mechanism of oil displacement in a capillary, but the microscopic insights into the detachment process are difficult to be obtained through experimental technologies. Moreover, experimental measurements of thinner capillaries are challenging because of the relatively small size, especially for the measurements of thinner capillaries of a nanoscale diameter. A more adequate approach to gain more microscopic insights into the oil displacement process, which is considered supplemental to experimental observations, is to carry out computational simulations. Over the past decades, molecular simulations have become an adequate approach to investigate the behavior of liquids confined in nanoscale geometries.12-16 Simulations of surfactants or liquid systems adsorbed on a slab of solid surface are now relatively common,17-20 however, to our knowledge, molecular simulation studies on the oil displacement process inside capillaries are scarce. Early simulation studies have focused on the spreading drop dynamics on a porous surface or liquid penetration into a cylindrical pore using molecular kinetic theory (MKT) dynamics.21-23 More recently dissipative particle dynamics (DPD) simulations have been applied to research fluids in nanoscale geometries. Millan and Laradji used the DPD method to understand structural and transport properties of driven polymer solutions in nanoscale


3


1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

channels.24 Chen et al. carried out many-body DPD simulations to study the water-oil displacement process in capillaries under an external force.25 They demonstrated that both strengthening the interaction between water and capillaries and weakening the applied force can displace the entire oil drop from the capillaries. In this work, we performed a series of simulation studies of oil displacement within a cylindrical silica pore using the all-atom molecular dynamics (MD) method. Silica was selected since it is a major composition of glass capillaries or rock minerals (e.g. quartz sandstone) in many geological environments. Our attentions focused on the effect of the surfactant, cetyltrimethylammonium bromide (CTAB), on the adsorption of oil drops on the pore wall with equilibrium MD simulations. Using the steered MD method, the displacement and migration process of an oil drop under an externally applied force was investigated. We concluded with some observations from the simulation that may be useful in understanding the mechanism of oil detachment for EOR experiments. 2. SIMULATION METHOD 2.1. Model Systems. For crude oil, a model proposed by Matsuoka et al. was adopted, which consists of eight kinds of hydrocarbon molecules, including hexane (HEX), heptane (HEP), octane (OCT), nonane (NON), cyclohexane (CHEX), cycloheptane (CHEP), toluene (TOL), and benzene (BEN) molecules.3 A model silica nanopore was obtained according to previous publications.14, 15 First, a cylindrical hole of diameter d = ~30 Å was carved from an amorphous silica block of dimensions Lx = Ly = 98.26 Å, Lz = 108.10 Å, by removing all the atoms lying along the z axis within the diameter (Figure 1). The pore was set to be large dimensions to provide sufficient space for the subsequent pulling simulations to take place. Then, silanol groups were used to


4

Page 4 of 30

Page 5 of 30

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


saturate the bare Si atoms on the internal pore walls of the silica block, resulting two kinds of silanol groups (Si-(OH)2 and Si-OH), adsorbed on the internal surface (Figure S1, Supporting Information), with a total surface density of ~7.0 –OH/nm2, in accordance with previous simulation studies.16 (Figure 1) In order to get a model of the residual oil inside the nanopore, System I was performed with oil molecules fully filling the inside of the cylindrical pore. The oil molecules included the eight hydrocarbon types, and the proportion of components was in accordance with previous studies (Table 1).3 Systems II and III were constructed to investigate the effect of pure water and surfactant solution on the residual oil by truncating an oil cylinder with a length of ~30 Å, in which the oil cylinder was selected from the final configuration of System I after a 10 ns MD run. We note that the proportion of components of the oil drop was very similar to that in System I, which suggested that the hydrocarbon species were mixed uniformly along the pore in System I. In constructing the model of System III, six CTAB molecules were placed at each side of the oil cylinder in the pore. After water molecules were added to fill the rest space of pores besides the oil cylinder in System II and III, the box length Lz was increased to 113.10 Å to set the local water density to the bulk value under ambient conditions ( ρbulk = 33.3 nm-3) in the subsequent MD simulations. In System III, bromine anions were inserted to keep electrically neutral. (Table 1)

2.2. Computational Details. GROMACS 4.6.3 software package26 was employed to carry out all the MD simulations. The all-atom optimized performance for the liquid systems (OPLS-AA) force field27 was adopted for all of the potential function terms to calculate the interatomic interactions. The simulation


5


1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

parameters for CTAB and oil molecules used in this study were derived from OPLS-AA force field, which is similar to recent studies on the adsorption of CTAB on carbon nanotubes.28 The partial atomic charges of CTAB and oil molecules were assigned according to the OPLS-AA force field because of the force field’s good representation of small organic molecules. The force field parameters of amorphous SiO2 were referred to the work of Lorenz,29 which gave good results in the dynamics properties of amorphous SiO2. The particle charge of the surface atoms of the pore was derived from the previous report,30 which gave a good prediction of adsorbed water layers on the silicon oxide surface. Water molecules were described by the simple point charge/extend (SPC/E) model31. More detailed information on the force field parameters could be found in Table S1 in the Supporting Information. Each of the systems was first minimized using the steepest descent method. Following the minimization, MD simulations under canonical ensemble (NVT) were carried out for each system with periodic boundary conditions applied in the x, y and z directions. During the minimization and the equilibration MD simulation, position restraints were applied to the silica crystal except for surface atoms. The V-rescale thermostat algorithm32 was adopted to keep a constant temperature of 298 K. Bond lengths were constrained by the LINCS algorithm33. The particle mesh Ewald (PME) method34 was used to compute electrostatics interactions. The Lennard-Jones interactions were cut off at 1.4 nm for the non-bonded potential. For each system, the simulations were performed for 10 ns with a time step of 1 fs. Trajectories were collected at intervals of 0.1 ns and were viewed by VMD 1.9.1.35 The last 5 ns trajectories were used for further analysis. To study the displacement process of the residual oil inside the pore, we performed steered molecular dynamics (SMD)36 (i.e., pulling simulations) using the pull code of the GROMACS


6

Page 6 of 30

Page 7 of 30

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


package.26 Structures derived from the equilibrated MD trajectories of System II and III were used to perform the pulling simulations. For System II, the driving force was applied on the whole oil drop (i.e., the pulled group), while for System III it was on the oil drop together with CTAB molecules which were solubilized into the oil phase. A point on the central axis of the nanopore was selected as a fixed reference using in the pulling simulations. Figure S2 in the Supporting Information shows a schematic representation of pulling direction. For each system, the external force would draw the pulled group from its original position along the z axis, using a potential with K = 1000 kJ mol-1 nm-2. Four pulling rates (0.0025, 0.005, 0.01, and 0.025 nm ps-1) were used to investigate the influence of pulling rate on the SMD simulations. The target distance of the pulling simulations was set to be 4 nm. The potential of mean force (PMF) as a function of the center of mass (COM) of the oil drop displacement ξ along the pore was calculated through a series of biased simulations with various ξ spanning the path of oil drop displacement. The 20 starting configurations which correspond to the 20 sampling windows were derived from trajectories of the pulling simulation with a window spacing of 0.2 nm, for systems with and without CTAB. In each window, a short equilibration MD run (500 ps) has been performed before the sampling simulations. Then, umbrella sampling37 simulations were performed for 5 ns in each window. The details of other methodology for pulling and sampling simulations are the same as the equilibrated MD simulations. The weighted histogram analysis method (WHAM)38 was used to analyze the results from the sampling simulations.

3. RESULTS AND DISCUSSION 3.1. Structural Features of Crude Oil in the Silica Nanopore. To confirm the equilibrium of System I, the potential energy with time evolution was checked


7


1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 8 of 30

(Supporting Information, Figure S3). It was found that the energy remained steady during the last 5 ns, suggesting that the system was well equilibrated. Moreover, considering the system is complex, we performed one more simulations for System I with different initial configuration to guarantee equilibration. The structural properties of oil molecules inside the nanopore were investigated by exploring the average number density profiles of relevant species in the pore:

ρα (r ) =

1 2π rLZ

∑ 〈δ (r − Rα )〉 i

(1)

i

where Riα denotes the distance between the position of the molecule of species α and the pore axis, and LZ denotes the length of the oil cylinder. The relevant species including carbon atoms of oil molecules and oxygen atoms of the silanol groups were computed, to represent the density distribution of oil molecules and the silanol groups on the nanopore wall. Figure 2a shows the density profiles along the cylindrical radial direction with respect to the central axis of the pore (i.e., z axis of the box) for System I. (Figure 2) Three evident peaks can be observed from the density profiles, which exhibited a layering organization of oil molecules inside the nanopore. The width of a peak in the density profile was about 0.5 nm, which is approximately equal to a diameter of one n-alkane molecule, suggesting that monolayers had formed. Figure 2b shows a view of the layering structure. From the figure, three distinct layers, which correspond to the three peaks shown in the density profile, can be observed. The three peaks have different intensities, which suggests that the three layers formed with different organization aggregations. The results from the repeated simulations for System I can be found in Figure S4 in Supporting Information, which also confirmed the final equilibration and reproducibility of System I. A detailed view of the three layers inside the pore


8

Page 9 of 30

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


was shown in Figure S5. We observed that the oil molecules in the outer layer covered the pore surface, resulting in well-organized aggregates, while in the inner layers, there are more vacancies compared with that in the outer layer. This is probably due to the different molecular orientation for oil molecules with respect to the axis of the pore. The orientation of oil molecules inside the pore was calculated, which is described by the angle θ between vectors of oil molecules and a vector that extends from the pore axis to the pore wall through the COM of a given molecule. The orientation was considered for three kinds of the oil molecules, i.e., n-alkanes, cycloalkanes, and aromatic hydrocarbons. The molecular vector of the n-alkane molecules was defined by connecting the two terminal carbon atoms and a schematic illustrating the orientation angle of n-alkane kind molecules between the two identified vectors is shown in Figure 2c. Note that when the value of θ is 90°, the corresponding molecular vector is parallel to the pore wall. If the value is 0° or 180°, the molecular vector is vertical to the pore wall. The orientation probability distribution for n-alkane molecules is illustrated in Figure 2d. The distribution of cycloalkanes and aromatic hydrocarbons was provided in Figure S6. It can be seen that the n-alkane and aromatic molecules have a preference for orienting parallel to the surface of the pore wall at a distance of ~1.5 nm from the pore axis. This means that these species in the outer layer tend to lie on the pore wall. In the inner layers, the oil molecules tend to lie at various angles to the pore surface. There is no preferential orientation of cycloalkanes inside the pore. This is because the structure of these molecules is neither linear nor planar. When they adsorbed onto the pore wall, part of the carbon atoms of the cycloalkanes anchored onto the pore. The fully covered and well-ordered layer of oil molecules adsorbed on the pore wall is considered to be crucial for the process of oil detachment from the silica surface.


9


1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 10 of 30

3.2. Structural Features of Residual Oil in the Silica Nanopore. The structural properties of residual oil in the presence of displacing fluid, i.e., pure water, (System II) and CTAB solution (System III) were investigated by analyzing the corresponding equilibrium MD trajectories. The equilibration of these systems was determined by monitoring the thickness of the oil cylinder with time evolution. The plots were shown in Figure S7 in the Supporting Information. It can be seen that the thickness achieved a stable equilibrium after the first several ps of simulation, suggesting that these systems were well equilibrated. (Figure 3) Figure 3 shows the structures at the beginning and the end of the simulation. It can be found that the oil drop did not migrate after the MD simulations. CTAB molecules were found to be solubilized into the oil phase from the water phase, with their hydrophobic chains penetrating the oil phase while hydrophilic headgroups stayed at the oil/water interface. To guarantee equilibration and reproducibility, two more simulations from different initial configurations were performed. Configurations from these simulations are presented in the Supporting Information, Figure S8. It can be noted that the same observations were obtained from these results. Note that, in practice the in-situ EOR is normally carried out in the region of the ground under high pressures. To ensure reasonability of the subsequent pulling simulations, we also performed simulations using NPT ensemble on the silica slab under high pressures. Figure S9 in the Supporting Information shows a configuration of the NPT MD run at a pressure of 100 atm. We can observe that CTAB molecules can also enter into different layers of the oil aggregations with their polar headgroups toward the water phase. The location of the oil drop inside the pore was characterized by the average number density distributions of components along the axial direction of the silica pore (Supporting Information,


10

Page 11 of 30

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


Figure S10). The local density of water inside the pore in each system was found to be equal to the bulk value under ambient conditions, ρbulk = 33.3 nm-3, and the oil molecules were found to be located in the same region which extended from ~2.0 to ~5.0 nm in the two systems. Since the simulation models of these systems were all derived from System I, the unchanged location region reflects the immobile oil drop in the silica pore. For the systems with CTAB, the distribution of the oil molecules expanded due to the surfactant molecules entering into the oil phase, which resulted in an increase of the thickness of the oil cylinder. (Figure 4) The structures of the oil cylinder in the two systems were also checked by the number density profile of oil molecules along the cylindrical radial direction, as shown in Figure S11. It can be seen that the layering organization still existed in the residual oil drops. However, the intensity of the three peaks was significantly reduced, especially for systems with CTAB. This indicates that the layering structures in the residual oil drop have been changed, after the displacing fluid filled in the pore. For the system without CTAB, the distributions of water molecules were observed between the oil droplet and the pore surface. Figure 4a shows that water molecules moved into the vicinity of the oil aggregates and formed a molecular channel. Interestingly, the water channel formed on the surface with Si-(OH)2 groups. The water channel formed along with the mobility of oil molecules, which left a gap for water molecules to enter. We considered that the formation of the channel was attributed to the competitive adsorption of water molecules on different densities of the hydroxyl groups on the pore surface. The competitive adsorption was investigated by calculating the corresponding binding energies using the density functional theory (DFT). Figure S12 shows that the interaction between Si-(OH)2 groups and water molecules was stronger than that between Si-OH groups and water molecules.


11


1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Figure 4b shows that the CTAB surfactant molecules entered into different layers of the oil aggregations, under the dispersive interaction of the chain with the oil molecules. In the case of CTAB adsorbed onto the silica surface, and water molecules aggregated around the polar head of CTAB, water channels would form between the oil phase and silica surface. With the perturbation the alkyl chains of the surfactant into the oil layer, the arrangement of the oil molecules changed. Moreover, the thickness of the oil cylinder increased more than in the absence of surfactants (Figure S7) due to the solubilization of surfactants into the droplet. Thus, we can conclude that the surfactant molecules are more efficient than pure water in changing the microstructure of the oil drop adsorbed on the pore wall. Although both the structures in the two systems with and without CTAB changed, the oil drop was still undissociated from the silica surface. The molecular orientation of n-alkane molecules in the two systems was also investigated by calculating the corresponding orientation angle with respect to the pore surface, as shown in Figure S13. It can be seen that the linear alkane molecules that adsorbed on the pore wall still tended to lie on the pore wall. To further mimic the oil drop detach from the silica surface, pulling simulations were further performed.

3.3. Dynamical Process of Pulling Simulations. In SMD simulations, an external force was applied on the oil drop to simulate the driving force from the displacing fluid in the process of oil recovery. By pulling on the COM of the oil drop, the oil molecules were forced to move along the silica pore. To illustrate the movement of the oil drop inside the pore, the distance between the COMs of the oil drop and the immobile reference with time evolution was calculated (Supporting Information, Figure S14). It can be seen that the oil drop was pulled along the pore for about 3.5 nm. To verify the reproducibility of these non-equilibrium MD simulations, the pulling simulations of each system were performed


12

Page 12 of 30

Page 13 of 30

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


for three times. The profiles of the external force against the reaction coordinate of the repeated SMD simulations are shown in Figure S15. From the figure, we found that the variation trend of the force for the parallel simulations was very similar, which confirmed the reproducibility of the SMD method. The curves of the external force were also plotted with time evolution of the SMD simulations, as shown in Figure S16. We found that the force rapidly increases during the first ~300 ps. Based on the variation of the COM distance (Figure S14), we found that the movement of the oil drop is much slower during this period than that in the subsequent period. Take System II as an example, a series of configurations with time evolution in the SMD simulations are shown in Figure 5. It can be seen that the oil molecules in the inner layers were found to be first pulled by the external force. (Figure 5) After reaching the maximum at ~300 ps, the force remained relatively constant. By checking the SMD trajectories, we found that from ~300 ps, the oil molecules in the outer layer were pulled and detached from the pore surface. To confirm the observed phenomenon, the number of carbon atoms in the out layer (i.e., the carbon atoms located beyond a distance of 1.3 nm from the central axis of the pore) was plotted with respect to the reaction coordinate, as shown in Figure S17. It can be found that the number of carbon atoms began to decrease at 1.5 nm (i.e., ~300 ps), which means that the oil molecules began to be detached from the silica surface. At that time, the force needed to overcome the adsorption action between oil molecules and the surface. With the oil molecules adsorbed on the pore wall gradually decreasing, less force was required to pull the oil drop and the force was reduced (at ~600 ps). Until the entire oil drop has been pulled along the pore for about 3.5 nm (at the end of the SMD simulation), the oil drop was


13


1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 14 of 30

considered to be detached from the silica pore. At the moment, a few oil molecules were detained on the pore wall. For the systems with CTAB, similar observations were obtained from SMD simulations (Supporting Information, Figures S14-18). To check the reproducibility of the observations from the pulling simulations, we performed one more simulation for System II by truncating the oil cylinder in a different position of the silica pore from System I. Some results including the views of SMD simulations and the curves of the applied force and carbon atoms in the out layer can be found in Figure S19 in the Supporting Information. These results showed good agreement with those derived from System II. The detachment processes of the oil droplet under the other three pulling rates (i.e., 0.0025, 0.01 and 0.025 nm ps-1) in the SMD simulations were also investigated. The curves of the external force as well as configurations with time evolution of the SMD simulations are shown in Figure S20. The force and number of carbon atoms in the out layer plotted with the reaction coordinate were also included as shown in Figure S21. We may see that the structural change is following the same sequence under each of the pulling rates. The variation trends of the external force are also generalized and the maximum is related to the structural change. Thus, the used pulling rates in a certain range do not affect the process of the SMD simulations. Actually, the rates of the pulling simulations may not be in accord with the velocities in the real situation. The biased MD simulations were used here to achieve a particular phenomenon (i.e. the migration of oil drop along the pore) that might be difficult to occur using the conventional MD methods. (Figure 6) To further identify the effect of the surfactant on the oil detachment process from the silica surface, a more detailed analysis was conducted based on the SMD trajectories. We focused on a


14

Page 15 of 30

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


region of the oil drop where there was a CTAB molecule inserted into the layers of the oil organization. Figure 6a shows a configuration of System III at 300 ps of the SMD simulation. As discussed above, the oil molecules in the outer layer were pulled in order to be detached from the pore surface at about 300 ps. From the figure, we note that water molecules can aggregate around the headgroup forming a water channel at the solid-oil interface. Because the hydrophilic head group of CTAB can strongly interact with water molecules. When the oil drop moved along the pore, the water channel expanded which facilitated the detachment of oil molecules from the pore wall (Figure 6b). However, this phenomenon was not observed system without CTAB, even when the water channel existed in the initial structure of the SMD simulation. Figure S22 in Supporting Information shows the configurations of the water channel in the system without CTAB. It can be seen that the channel did not expand with the displacement of the droplet. Conversely, it became smaller due to the adsorption of oil molecules onto the pore wall.

3.4. Potential of Mean Force. The potential of mean force (PMF) of the oil detachment process inside the silica pore was determined using umbrella sampling method. The umbrella sampling method is an useful method to assemble the PMF.39 To perform the umbrella sampling, a series of configurations is generated by the SMD simulations along a reaction coordinate ξ between two groups first (e.g. the oil drop and the immobile reference in this work). A selected increasing COM distance between the two groups represents the sampling windows. Then, the position of the pulled group was maintained by a biased potential and independent sampling simulations were performed in each window. The PMF curves were plotted with the entire reaction coordinate using WHAM38 method, which led to the free energy ∆G of dissociation of the oil drop from the silica surface and movement along the pore.


15


1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 16 of 30

We have performed sampling simulations from initial configurations derived from pulling simulations with different rates (0.01, 0.005 and 0.0025 nm ps-1) for systems with and without CTAB. These curves are shown in the Supporting Information, Figure S23. As illustrated in these figures, a convergence was achieved by discarding the first 2 or 3 ns of the trajectories. Thus, the last 2ns of the simulation time was used for WHAM in all the sampling simulation systems. The PMF profiles for the detachment of the oil drop through the pore in System II and III under each of the pulling rates are shown in Figure 7. The histograms of the configurations within the umbrella sampling windows were plotted in Figure S24 in Supporting Information. The histograms show reasonable overlap between windows of the COM spacing, suggesting the entire reaction coordinate was properly sampled. (Figure 7) Figure 7 shows clearly that the free energy change is consistently smaller in systems with CTAB than those without a surfactant. Although there are certain differences among the PMF curves when using different pulling rates, the overall trend of the energy change is consistent. In the beginning of the reaction coordinate (i.e., the first 0.75 and 1.50 nm for System II and III respectively), the energy profiles show a slow-growth. This is because that the migration of the oil drop is very slow at the beginning of the SMD simulation, as discussed above. In this stage, we tend to believe that both the adsorption interaction of the drop onto the wall and dispersion interactions among the hydrocarbon molecules make the main contributions to the free energy, which suggests that the oil drop tended to keep their positions inside the pore when the external force was small. And the force could not counteract the adsorption energy between the oil drop and the pore wall. Thus, extra energy was required to make the droplet migrate, yielding an increment of ∆G in the two systems. The overall positive value of ∆G indicates a non-


16

Page 17 of 30

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


spontaneous process of the movement of the oil drop along the pore. In the subsequent stage, the energy profiles show a faster-growth. According to the migration displacement of the drop and the number of carbon atoms in the out layer (Supporting Information, Figure S14 and S17), we can see that the oil molecules in the outer layer of the drop detached from the pore wall under the external force at the moment. In the systems in the presence of CTAB, water channels that formed around the headgroups of the surfactant molecules expanded on the silica surface which caused a wettability alteration of the pore surface and facilitated the detachment process. It made the force overcome the interactions between the droplet and the pore wall easily. Thus, less extra energy was required to make the oil drop migrate in the presence of the surfactant solution. In addition, the increased thickness of the oil cylinder was another factor for the lower energy in the presence of CTAB, which resulted in a loose pack of oil molecules. We can conclude that in contrast with water, surfactant solutions may facilitate the process of oil detachment. Therefore, in order to enhance oil-displacement efficiency, it is important to reasonably determine the role of surfactant molecules on the inversion of wettability, which helps to reduce the interfacial tension at the oil-solid interface.

4. CONCLUSIONS A series of MD simulations were performed to study the effect of a surfactant on the oil detachment process inside a silica nanopore. With equilibrium MD simulations the microscopic structures of crude oil inside the pore were revealed. The results showed that there is a layering organization in the crude oil inside the silica pore, with one well-ordered layer adsorbed on the pore surface. To assess the mechanism of displacing fluid on residual oil, two systems with and without a surfactant were studied. The MD results showed that with the surfactant solution, the layering structures were evidently disturbed by the solubilization of surfactant molecules into the


17


1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 18 of 30

oil phase. In the case of CTAB adsorbed onto the pore wall, water molecules around the polar head formed water channels between the oil phase and silica surface, which is vital to the oil detachment. Steered MD method was used to study the displacement and migration process of the residual oil drop inside the silica pore. Detailed information about the process was achieved from SMD simulations. Briefly, oil molecules in the central part of the oil drop moved first under the applied force. At that time, the external force counteracted interactions of oil molecules with the pore wall. The structural change of the outer layer on the pore wall is important to oil detachment. Due to expansion of the water channel around the headgroup of surfactant molecules, oil molecules in the outer layer are more easily pulled away from the surface. The free energy profiles showed that less energy is required for the fluid to displace the oil droplet along the pore with surfactants present. These results are expected to provide molecular-level insight into the mechanism of surfactants aiding in displacing oil drops in enhanced oil recovery.

AUTHOR INFORMATION Corresponding Authors *E-mail: [email protected] Supporting Information Detailed view of the pore surface; Force field parameters; Illustration of pulling direction; Potential energy with time of System I; Layering structures inside the pore; Molecular orientation angles; Time evolution of the thickness of the oil cylinder; Number density profiles; COM distance between the oil drop and the reference; Results from repeated simulations for pulling simulation; Number of carbon atoms adsorbed on the pore wall; SMD configurations of System III; Pulling simulation using different rates and configuration of the water channel in a


18

Page 19 of 30

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


system without CTAB; PMF curves using different pulling rates and histograms of the configurations within the umbrella sampling windows. This material is available free of charge via the Internet at http://pubs.acs.org.

ACKNOWLEDGMENT This work was financially supported by the National Natural Science Foundation of China (21203084). Thanks to Dr. Edward C. Mignot, Shandong University, for linguistic advice.

REFERENCES (1) Salathiel, R. A. Oil recovery by surface film drainage in mixed-wettability rock. J. Pet. Technol. 1973, 25, 1216-1224.

(2) Petersen, K. S.; Christensen, P. L. Phase Behavior of Petroleum Reservior Fluids; CRC Press: Boca Raton, London and New York, 2007. (3) Kunieda, M.; Nakaoka, K.; Liang, Y.; Miranda, C. R.; Ueda, A.; Takahashi, S.; Okabe, H.; Matsuoka, T. Self-accumulation of aromatics at the oil−water interface through weak hydrogen bonding. J. Am. Chem. Soc. 2010, 132, 18281-18286. (4) Ghatee, M. H.; Fotouhabadi, Z.; Zolghadr, A. R.; Borousan, F.; Ghanavati, F. Structural and phase behavior studies of pyridine and alkyl pyridine at the interface of oil/water by molecular dynamics simulation. Ind. Eng. Chem. Res. 2013, 52, 13384-13392. (5) Kumar, N.; Varnasi, K.; Tilton, R.; Garoff, S. Surfactant self-assembly ahead of the contact line on a hydrophobic surface and its implications for wetting. Langmuir 2003, 19, 53665373. (6) Cachile, M.; Cazabat, M. Spontaneous spreading of surfactant solutions on hydrophilic surfaces: CnEm in ethylene and diethylene glycol. Langmuir 1999, 15, 1515-1521. (7) Starov, V.; Kosvintsev, S. R.; Velarde, M. G. Spreading of surfactant solutions over hydrophobic substrates. J. Colloid Interface Sci. 2000, 227, 185-190. (8) Tiberg, F.; Zhmud, B.; Hallstensson, K.; Von Bahr, M. Capillary rise of surfactant solutions. Phys. Chem. Chem. Phys. 2000, 2, 5189-5196.

(9) Hammond, P. S.; Unsal, E. Spontaneous and forced imbibition of aqueous wettability altering


19


1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 20 of 30

surfactant solution into an initially oil-wet capillary. Langmuir 2009, 25, 12591-12603. (10) Hammond, P. S.; Unsal, E. Forced and spontaneous imbibition of surfactant solution into an oil-wet capillary: The effects of surfactant diffusion ahead of the advancing meniscus. Langmuir 2010, 26, 6206-6221.

(11) Zhang, H.; Nikolov, A.; Wasan, D. Dewetting film dynamics inside a capillary using a micellar nanofluid. Langmuir 2014, 30, 9430-9435. (12) Lindsey, R. K.; Rafferty, J. L.; Eggimann, B. L.; Siepmann, J. I; Schure, M. R. Molecular simulation studies of reversed-phase liquid chromatography. J. Chromatogr. A 2013, 1287, 60-82. (13) Bródka, A.; Zerda, T. W. Properties of liquid acetone in silica pores: molecular dynamics simulation. J. Chem.Phys. 1996, 104, 6319-6326. (14) Elola, M. D.; Rodriguez, J.; Laria, D. Structure and dynamics of liquid methanol confined within functionalized silica nanopores. J. Chem.Phys. 2010, 133, 154707. (15) Elola, M. D.; Rodriguez, J.; Laria, D. Liquid methanol confined within functionalized silica nanopores. 2. Solvation dynamics of coumarin 153. J. Phys. Chem. B 2011, 115, 1285912867. (16) Ghoufi, A.; Hureau, I.; Morineau, D.; Renou, R.; Szymczyk, A. Confinement of tert-butanol nanoclusters in hydrophilic and hydrophobic silica nanopores. J. Phys. Chem. C 2013, 117, 15203-15212. (17) Liu, Q.; Yuan, S. L.; Yan, H.; Zhao, X. Mechanism of oil detachment from a silica surface in aqueous surfactant solutions: molecular dynamics simulations. J. Phys. Chem. B 2012, 116, 2867-2875. (18) Zehl, T.; Wahab, M.; Schiller. P.; Mögel. H. J. Monte carlo simulation of surfactant adsorption on hydrophilic surfaces. Langmuir 2009, 25, 2090-2100. (19) Seveno, D.; Blake, T. D.; Goossens, S.; De Coninck, J. Predicting the wetting dynamics of a two-liquid system. Langmuir 2011, 27, 14958-14967. (20) Zhang, X. F.; Lu, G. W.; Wen, X. M.; Yang, H. Molecular dynamic investigation into the adsorption of oil-water-surfactant mixture on quartz. Appl. Surf. Sci. 2009, 255, 6493-6498. (21) De Ruijter, M. J.; De Coninck. J.; Oshanin, G. Droplet spreading: partial wetting regime revisited. Langmuir 1999, 15, 2209-2216. (22) Seveno, D.; Ledauphin, V.; Martic, G.; Voué, M.; De Coninck. J. Spreading drop dynamics


20

Page 21 of 30

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


on porous surfaces. Langmuir 2002, 18, 7496-7502. (23) Martic, G.; Blake, T. D.; De Coninck, J. Dynamics of imbibition into a pore with a heterogeneous surface. Langmuir 2005, 21, 11201-11207. (24) Millan, J. A.; Laradji, M. Cross-stream migration of driven polymer solutions in nanoscale channels: A numerical study with generalized dissipative particle dynamics. Macromolecules

2009, 42, 803-810. (25) Chen, C.; Zhuang, L.; Li, X. F.; Dong, J. F.; Lu, J. T. A many-body dissipative particle dynamics study of forced water–oil displacement in capillary. Langmuir 2012, 28, 13301336. (26) van der Spoel, D.; Lindahl, E.; Hess, B and the GROMCAS development team, GROMACS User Manual version 4.6.3; Royal Institute of Technology and Uppsala University: Stockholm and Uppsala, 2013.

(27) Jorgensen, W. L.; Maxwell, D. S.; Tirada-Rives, J. Development and testing of the opls allatom force field on conformational energetics and properties of organic liquids. J. Am. Chem. Soc. 1996, 118, 11225-11236.

(28) Lin, S.; Hilmer, A. J.; Mendenhall, J. D.; Strano, M. S.; Blankschtein, D. Functionalization of single-walled carbon nanotubes in aqueous surfactant solutions. J. Am. Chem. Soc. 2012, 134, 8194-8204.

(29) Lorenz, C. D.; Webb III, E. B.; Stevens, M. J.; Chandross, M.; Grest, G. S. Frictional dynamics of perfluorinated self-assembled monolayers on amorphous SiO2. Tribol. Lett.

2005, 19, 93-98. (30) Wensink, E. J. W.; Hoffmann, A. C.; Apol, M. E. F.; Berendsen, H. J. C. Properties of adsorbed water layers and the effect of adsorbed layers on interparticle forces by liquid bridging. Langmuir 2000, 16, 7392-7400. (31) Berendsen, H. J. C.; Grigera, J. R.; Straatsma, T. P. The missing term in effective pair potentials. J. Phys. Chem. 1987, 91, 6269-6271. (32) Bussi, G.; Donadio, D.; Parrinello, M. Canonical sampling through velocity rescaling. J. Chem. Phys. 2007, 126, 014101.

(33) Hess, B.; Bekker, H.; Berendsen, H. J. C.; Fraaije, J. G. E. M. LINCS: A linear constraint solver for molecular simulations. J. Comput. Chem. 1997, 18, 1463-1472. (34) Essman, U.; Perera, L.; Berkowitz, M. L.; Darden, T.; Lee, H.; Pedersen, L. G. A smooth


21


1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 22 of 30

particle-mesh ewald method. J. Chem. Phys. 1995, 103, 8577-8593. (35) Humphrey, W.; Dalke, A.; Schulten, K. VMD: Visual Molecular Dynamics. J. Molec. Graphics 1996, 14, 33-38.

(36) Izrailec, S.; Stepaniants, S.; Isralewitz, B.; Kosztin, D. Lu, H.; Molnar, F.; Wriggers, W.; Schulten, K. Steered Molecular Dynamics, Springer-Verlag, Berlin, 1998. (37) Torrie, G. M.; Valleau, J. P. Nonphysical sampling distributions in monte carlo free-energy estimation: umbrella sampling. J. Comput. Phys. 1977, 23, 187-199. (38) Kumbia, S.; Bouzida, D.; Swendsen, R. H.; Kollman, P. A.; Rosenberg, J. M. The weighted histogram analysis method for free-energy calculations on biomolecules. I. The method. J. Comput. Chem. 1992, 13, 1011-1021.

(39) Lemkul, J. A.; Bevan, D. R. Assessing the stability of alzheimer’s amyloid protofibrils using molecular dynamics. J. Phys. Chem. B 2010, 114, 1652-1660.


22

Page 23 of 30

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


Figure 1. Representation of the silica nanopore model (left panel) and details of the pore wall functionalization (right panel). The axis label on the left panel shows the orientation of the nanopore in the simulation box relative to the z-axis, and d is the diameter of the silica pore.

Table 1. Details of the Simulation Systems Components of Oil system

CTA+

Br-

Water

HEX

HEP

OCT

NON

CHEX

CHEP

TOL

BEN

I

54

54

60

72

36

60

60

24

-

-

-

II

13

12

13

21

7

13

7

6

-

-

2431

III

13

12

13

21

7

13

7

6

6

6

2312


23


1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Figure 2. (a) Number density profiles along the cylindrical radial direction with respect to the zaxis of the pore; (b) Snapshots of the configurations of oil molecules inside the silica nanopore from the cylindrical axial view; (c) Scheme of orientation angle θ between vectors of the alkane molecules and a vector that extends from the pore axis to the pore wall through the COM of a given alkane molecule; (d) Probability distribution of the orientation angle for alkane molecules at each distance from the pore axis.


24

Page 24 of 30

Page 25 of 30

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


Figure 3. Configurations at the beginning (left) and the end (right) of the simulations for Systems II and III. Color scheme: Oil molecules, cyan; CTA+, violet, and Br-, gray.


25


1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Figure 4. Detailed views of the partial structures for Systems II (a) and III (b).


26

Page 26 of 30

Page 27 of 30

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


Figure 5. Configurations at different periods of the pulling simulation for System II.


27


1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Figure 6. Configuration of the water channel at 300 ps (a) and 400 ps (b) of the SMD simulation


28

Page 28 of 30

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


Potential of Mean Force / kcal/mol

Page 29 of 30

60 -1

v = 0.005 nm ps

without CTAB -1

v = 0.01 nm ps

40

-1

v = 0.0025 nm ps

20 with CTAB

0 0.5 1.0 1.5 2.0 2.5 3.0 3.5 4.0 4.5 Reaction coordinate ξ / nm

Figure 7. Potential of mean force (PMF) profiles as a function of COM distance between the oil drop and reference for the two systems in this study.


29


1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

TOC Graphics


30

Page 30 of 30

Molecular Dynamics Simulation of the Oil ... - ACS Publications

Recommend Documents