Fig. 3 Initial configurations of the seven PH slit pore models
of 5 nm width with varying water concentration. Cw stands for the water
concentration. The color codes of clay are as in Fig. 1. In the first
panel, the direction of the arrows indicates the direction of an imposed
acceleration in subsequent sections of this paper. The fluid-color
regimes: ethane (light green), H2O (red) and dodecane
(light blue).
2.2 Simulation Details
Our modeling relies on the use of
the ClayFF and OPLS All-Atom force field to describe
illite53,54 and organic components55respectively, while the flexible SPC model and the shake algorithm are
applied to model water. Lorentz-Berthelot mixing rules describe the
interactions between different atoms. We use LAMMPS (Large-scale
Atomic/Molecular Massively Parallel Simulator) 56 with
periodic boundary conditions applied in 3 directions. A more detailed
description is available in our previous work42,48.
Our workflow is as follows: In the initial set-up, fluid molecules are
placed randomly between the slit pore using the Packmol
package57. We then run equilibrium MD (EMD)
simulations using an NPT ensemble for 5 ns. After equilibrium of
isothermal–isobaric (NPT) ensemble, we further run the simulations for
another 5 ns with a canonical (NVT) ensemble. Temperature and pressure
are maintained at 400 atm and 350 K using Parrinello-Rahman
barostat58 and Nose Hoover thermostat59 respectively.
After EMD simulations, non-equilibrium MD (NEMD) simulations are
performed for another 10 ns to mimic hydrocarbon-water transport.
Several methods exist for inducing flow in molecular dynamics, such as
forced60–62, surface-induced, pressure
difference63, osmotic pressure64 and
gravitational approaches65–67. In our simulations, we
adopt the gravitational technique because other methods have been known
to cause a distortion of the velocity-field and density profile along
the flow direction67–69. A uniform external force
(along the x-direction as shown in Fig. 3), is imposed on all molecules
inside the illite nanopores70,71. Published literature
suggests using small accelerations on the order of
10-4 to 10-3 nm/ps2
72,73 to avoid velocity jumps in the system and therefore, in our
study, the acceleration ranges from 0.0005 to 0.002
nm/ps2.
It should be noted that in MD simulations, fluid temperature is
calculated from the kinetic energy74. However, in NEMD
simulations, the molecule velocities in the direction of the driving
force consists of both thermal velocities and the imposed center-of-mass
velocity. To ensure the temperature is maintained, our fluid temperature
calculations do not include the center-of-mass
velocity69. The temperature and pressure during the
NEMD simulations are shown in Fig. S2 in the Supporting Information.