Article pubs.acs.org/IECR
Dissipative Particle Dynamics Study on the Aggregation Behavior of Asphaltenes under Shear Fields Xianyu Song,†,‡ Peng Shi,‡,§ Shuangliang Zhao,∥ Ming Duan,*,‡,§,⊥ Chengjie Wang,‡ and Yongzhang Ma‡ †
Department of Mechanical and Electrical Engineering, Dazhou Vocational and Technical College, Dazhou, Sichuan 635000, P. R. China ‡ College of Chemistry and Chemical Engineering, Southwest Petroleum University, Chengdu 610500, P. R. China § Oil & Gas Field Applied Chemistry Key Laboratory of Sichuan Province, Chengdu 610500, P. R. China ∥ State Key Laboratory of Chemical Engineering, East China University of Science and Technology, Shanghai 200237, China S Supporting Information *
ABSTRACT: In the present work, the effects of shear fields on the aggregation of asphaltene molecules in heptane were investigated by means of dissipative particle dynamics simulations. The geometries of asphaltene aggregates without shear fields were studied, and the simulation results provide an interpretation of the experimental results on the microscopic level. The effects of shear fields on asphaltene aggregates were also investigated by accessing the radial distribution functions, spatial orientation correlation functions, and the radii of gyrations. We show that the shear fields can destroy the conformational order of the aggregates by damaging the organized structure and isolating the asphaltenes. As the radius of gyration results show, the asphaltene molecules are elongated to be alike-polymers by shear fields. Moreover, the reason why the viscosity decreases under shear fields is that the shear fields lead to the increase of dimerization free energies.
1. INTRODUCTION Asphaltenes are usually classified as the most polar fraction of crude oil, which are soluble in some solvents, such as aromatics (benzene, toluene, etc.), but insoluble in others, such as normal alkanes (paraffins, heptane, etc.).1,2 Asphaltenes have been identified as the molecules containing polyaromatic and polycyclic rings with short aliphatic chains and heteroatoms, such as nitrogen, oxygen, and sulfur.3 Because of these special molecular structures, asphaltenes can easily self-associate with themselves, thus leading to precipitation and deposition in crude oils.4,5 It is well-known that the precipitation and deposition of asphaltenes may cause many serious problems in petroleum industries during the processes of production, transport, storage, and refinement of crude oil,6−8 and to resolve these problems a thorough interpretation of the asphaltene deposition is needed. This aggregation process is variously known as agglomeration, coagulation, or flocculation.9 Most flocculators are operated under shear to ensure complete mixing, high collision frequency between particles, and therefore rapid floc growth.10 As floc grows larger, fluid shear stresses can also break the floc into smaller fragments, and the breakage tends to occur at the weakest points in the floc structure. Thus, these produced fragments are stronger and more compact but smaller than the parent floc.11 When oil mixtures are under shear, both aggregation and break age processes are important and shall © XXXX American Chemical Society
be investigated for understanding the aggregation structure and aggregate size distribution of asphaltenes.3 These problems have been studied for the last 50 years; however, the aggregation behavior and the aggregate structure of asphaltenes under shear fields are still been poorly understood because of the strong dependence on the experimental conditions and also the absence of effective approaches.12−14 Fortunately, different simulation methods, such as Monte Carlo simulation (MC), stochastic rotation dynamics (SRD), molecular dynamics (MD), and dissipative particle dynamics (DPD) are applicable to address these problems. Recently, dynamic evolution of asphaltene aggregate models based on Monte Carlo simulation have been developed to predict the dynamic evolution of asphaltene aggregate size distribution.15 In addition, the aggregation and breakage processes of asphaltene aggregates under shear-induced petroleum mixtures are also investigated by using Monte Carlo simulation.16 Barcenas et al.17 employed an off-lattice MC simulation to study the irreversible cluster−cluster aggregation in the presence of associative inhibitors. The stochastic rotation dynamics (SRD) method was used by Boek Received: June 22, 2016 Revised: July 27, 2016 Accepted: August 1, 2016
A
DOI: 10.1021/acs.iecr.6b02400 Ind. Eng. Chem. Res. XXXX, XXX, XXX−XXX
Article
Industrial & Engineering Chemistry Research
2. METHODOLOGY 2.1. DPD Method. DPD simulation, as a stochastic simulation technique, was first introduced by Hoogerbrugge27 and Koelman.28 DPD is often employed for studying the dynamical behaviors of complex fluid systems. In the DPD method, each bead represents a group of molecules or atoms, of which the size is large on the atomistic scale but still macroscopically small. All beads in a system interact with each other through three kinds of forces, and they are conservative repulsive forces, representing excluded volume; dissipative forces, representing viscous drag; and random forces. An extensive description on DPD simulation can be found elsewhere.22−24 Some general aspects of DPD technique are also presented in the Supporting Information.Through the introduction of bonded potentials including bond potential, angle potential, and inversion angle potential, the island and archipelago architecture asphaltene models are constructed in this paper. 2.2. Coarse-Graining Method. Though it is well-known that each asphaltene molecule consists of condensed aromatic rings and peripheral alkane side chains with some heteroatoms (like N, O, and S atom) in the form of functional groups such as acid or base, the degree and the way of condensed aromatic rings of asphaltenes are in dispute at present.29,30 It was common to take the idea that asphaltene molecules are island architecture and archipelago architecture. The island architecture is usually shaped “like your hand” with a single aromatic core (palm), which has alkane substitution (fingers). This is also termed the “like your hand” model. Another widely accepted model for asphaltenes is an archipelago architecture, which was thought to be composed of several small fused aromatic rings connected by bridge chains of small molecular weights. The archipelago architecture model is more applicable to light crude oil rather than to heavy or ultra-heavy oil. In our study, two typical asphaltene architectures are designed. As shown in Figure 1a, a coarse-grained model of the fused aromatic rings is constructed by creating the rigid sheet of the hexaparticle ring. To effectively coarse-grain the asphaltenes and heptane, three typical beads have been identified as the building blocks of the different molecular components, as shown in Figure 1e. For the three typical beads, B represents the moiety of aromatic rings, which is denoted by benzene molecule; H corresponds to the alkyl chain; and T is the functional group containing heteroatoms. The fused aromatic ring core presented in Figure 1a contains 22 benzenoid rings and which were kept as “square shape”. In this paper, the heptane is clustered into two identical coarse-grained particles which are represented as two butane molecules. The thiourea is selected as the heteroatom group in asphaltenes. Two kinds of island architecture and archipelago architecture coarse-grained models of asphaltenes are presented in Figure 1c,d. In the DPD simulation, in general, a bead corresponds to Nm water molecules. The number Nm (degree of coarse-graining) can be viewed as a real-space renormalization factor.25 In our work, we group three water molecules as one DPD particle, and this treatment has already been confirmed to produce ideal crude oil systems.24,25 The length scale RC in angstroms, mass scale m, and the time scale τ in picoseconds can be evaluated, as RC = 3.107(ρNm)1/3 Å, m = Nm·mwater amu, τ = (1.41 ± 0.1) Nm5/3 ps,25 where ρ is the DPD number density and mwater is the mass of the water molecule. In practice, because the number of bead−bead interactions increases with density,22 the
et al. to study the deposition of colloidal asphaltenes in capillary flow.18 It has been shown that the simulation results on the properties of asphaltene colloid aggregates agree well with experimental measurements. MC and stochastic rotation dynamics simulations, however, hinder their applications to obtain further information on the aggregation structure of asphaltenes at the molecular or mesoscopic scale. In MC simulation, the asphaltene molecule is treated with a discotic seven-center model and the resins are represented as singlesphere particles with selected volumes. Molecular dynamics can provide specific information on the behavior of asphaltene aggregates without losing atomistic features. By focusing on the effect of the asphaltene molecule structure, organic solvent, and temperature, the association and aggregation behaviors of asphaltenes were studied by using molecular dynamics.19−21 The dissipative particle dynamics is a widely used mesoscale simulation approach, and it can deal with systems with much larger size and time scale compared with MD simulation. On one hand, by grouping, for example, two or three water molecules into a single coarse-grained bead, DPD provides an ideal tool to investigate the asphaltenes properties in crude oil. On the other hand, because of the introduction of dissipative force and the coup of friction coefficient and noise amplitude, DPD is an excellent methodical method for the simulation of coarse-grained systems over long length and time scales.22 Alvarez et al. simplified one asphaltene molecule as a core and aliphatic regions. In this coarse-graining process, only one bead represents each core region and one bead denotes each aliphatic chain. They used these larger coarse-grained models to simulate crude oil−water emulsions in the presence of a functionalized copolymer.23 By introducing a rotational algorithm for the motion of rigid body into dissipative particle dynamics, Zhang et al.24 studied the aggregation behaviors of asphaltenes in heavy crude oil. Ruiz-Morales and Mullins proposed a different “mapping” of the fused aromatic ring. The fused aromatic ring region of asphaltene was coarse-grained as an ovalene core: each fused benzenoid ring is represented by one bead. Ovalene is composed of ten fused aromatic hexagonal rings. Based on this coarse-grained model, the orientation at the oil−water interface of asphaltenes with peripheral oxygen moieties was investigated.25 Skartlien et al. considered asphaltene as two types of beads and constructed six hexagonal rings in the core structure of asphaltene. They studied asphaltene adsorption on hydrophilic substrates by changing the polar groups and solubility.26 Though the theoretical investigations about the aggregation behavior of asphaltenes has been the focus of these recent works, aggregation behavior of asphaltenes under shear conditions has not been explored. In this work, two types of coarse-grained asphaltene models are adopted. One is the island architecture model, and the other is the archipelago architecture model, which can be interpreted as the integration of three island architecture models. The bonded potentials including bond potential, angle potential, and inversion angle potential were introduced to construct the asphaltene with the rigid planes of aromatic rings. DPD simulations are performed to investigate the aggregation behaviors of asphaltenes with or without shear fields. The variations of viscosity under shear fields have also been discussed in this study. B
DOI: 10.1021/acs.iecr.6b02400 Ind. Eng. Chem. Res. XXXX, XXX, XXX−XXX
Article
Industrial & Engineering Chemistry Research
Table 1. Hansen Solubility Parameters (J/cm3)1/2 and Molar Volume (cm3/mol) at 298 K34 material
Hansen solubility parameter δ at 298 K (J/cm3)1/2
molar volume at 298 K (cm3/mol)
benzene butane thiourea
18.51 14.10 33.01
89.4 101.4 72.8
calculated by Blends method, and the calculated value of the interlayer distance of asphaltenes is distributed in the range of 3.75−4.05 Å, which is close to the experimental data of 3.35 Å.24 Therefore, aHH = 81.0 and aTT = 47.0. Consequently, all the parameters of conservative force are listed in Table 2. Table 2. Parameters of Conservative Force bead
B
H
T
B (benzene) H (butane) T (thiourea)
62.05 80.4 100.5
80.4 81.05 119.1
100.5 119.1 47.05
The DPD simulations were carried out in NVT ensembles using the Mesocite module embedded in the Materials Studio 6.0 package.31 Cui et al. and Liu and Zhong used this package for studying the microphase separation of diblock copolymer melts in the presence of a steady shear flow, showing that this package is validated to access nonequilibrium properties.35,36 In our work, all DPD simulations were accomplished in a cubic box with the size of 60 × 60 × 60Rc3. The periodic boundary conditions were applied at three directions,and the system temperature was fixed at 298 K. The enclosed simulation cells are composed of 20%, 30% asphaltenes, the remaining fractions are heptanes which are selected as organic solvent. The total number of beads is 2378 when the density of all systems is set to 3.0 in reduced units. A total of 6 × 106 DPD simulation steps were carried out with a time step Δt = 0.005τ. The scales used in the Mesocite DPD units were as follows: length scale, 6.46 Å; mass scale, 54 amu; energy scale, 0.59191 kcal/mol; time scale, 3.0158 ps. The total real dynamic time ttotal = Nsteps·Δt·τtime scale, namely, 90.5 ns.
Figure 1. Coarse-grained model molecules of (a) fused aromatic rings, (b) solvent heptane, (c) archipelago architecture model of asphaltenes, and (d) island architecture model of asphaltenes; (e) corresponding schematic representation of coarse-grained beads in simulations.
DPD algorithm is most efficient when the density ρ is set to 3.0.22The length and time scales in physical units are RC = 6.46 Å, m = 54 amu, and τ = 8.8 ps with Nm = 3, the rigid sheet of the hexaparticle ring = 3.0. 2.3. DPD Parameters and Simulation Details. The results of DPD simulations are usually determined by two parameters, namely, the mode for coarse-grained of molecules and the interactions of DPD particles. The bead−bead interaction parameters are determined by the following equation:22 aij = aii + 3.27χij (1) To satisfy the compressibility of water, the mapping of three water molecules per bead leads to aii = 78.0.22,32 The values of χij can be calculated from the solubility parameters by the equation vij (δi − δj)2 χij = (2) RT where vij is the average of molar volumes of two beads; δi and δj are the solubility parameters of beads i and j, respectively. The Hansen solubility parameters are a good choice to calculate the χij. Shi et al. employed the Hansen solubility parameters as solubility parameters δ, and this method has already been confirmed to produce ideal interfacial tensions in the water/ benzene/caprolactam system.33 Ruiz-Morales et al. used the Hansen solubility parameters to investigate the properties of asphaltenes at the oil−water interface.25 The Hansen solubility parameters in this work used for different components are given in Table 1.34The parameter aBB is roughly determined as aBB = 62.0 because the beads of type B are placed in the rigid sheets and the bead density must be larger than the average density value of the system.24 The particles of type H and T have a particle density similar to that presented in ref 24, which were
3. RESULTS AND DISCUSSION To study the effects of shear rates on the aggregation behaviors of asphaltenes, a set of DPD simulations with different shear rates (from 0 to 0.8 ps−1) at the y-direction for island and archipelago asphaltene models at concentrations of 20% and 30%, respectively, were performed. The shear fields were carried out after 6 × 106 steps system aggregates finished. All simulations with shear rates were also carried out with 6 × 106 steps. To examine whether the simulations are completely equilibrated, we checked the time evolution of pressure and total potential energy of the system. The results are shown in Figure S2 in the Supporting Information. Figure S2a shows that the pressure in terms of shear rate decreases dramatically and reaches equilibrium rapidly, and FigureS2b plots the total potential energies at different shear rates, and it shows that the total potential energy reaches the individual stable stage within 90.5 ns. It implies that the simulation time of 6 × 106 steps is sufficiently long for our systems to reach a steady state. As a first step, we validate our DPD calculation by accessing the diffusion coefficient and interlayer distance of asphaltenes. The diffusion coefficient and interlayer distance of asphaltenes C
DOI: 10.1021/acs.iecr.6b02400 Ind. Eng. Chem. Res. XXXX, XXX, XXX−XXX
Article
Industrial & Engineering Chemistry Research
intermolecular attractive forces, while the larger colloids are self-associates, formed from these small aggregates by weaker interaction forces.39 To get further information on aggregation structure, Figure 3a,b gives the morphology of asphaltene aggregates. Two
obtained by DPD simulations agree well with the existing experiment data.24,37 The diffusion coefficient values were calculated by mean square displacements; more details can be found in Figure S3. The diffusion coefficient of diluted asphaltene in toluene is 2.26−3.38 × 10−10 m2 s−1, which is similar to experimental data, 2.2−6.3 × 10−10 m2 s−1.37 The interlayer distance value calculated by DPD simulations is about 5.0 Å for the system without shear fields and slightly larger than the experimental data (∼3.35 Å) and DPD simulation results (3.75−4.05 Å).24 3.1. Aggregation Structure without Shear Fields. 3.1.1. Aggregation Structure. To study the effects of shear fields, the aggregation process with no shear fields was first carried out by DPD dynamics. The simulated results for the systems which are composed of continental and archipelago asphaltene models are shown in Figure 2. The simulated
Figure 3. Morphology of asphaltene aggregates obtained by DPD simulations (a, b), aggregate formation mechanism based on the modified Yen-Mullins model for (c−e),40 and micrographs of 5000 ppm asphaltenes in (f) 80:20 and (g, h) 70:30 heptane solutions imaged on Si3N4 grids.41 The same one color represents one asphatene molecular in panels a and b. Panels c−e, reprinted from ref 40. Copyright 2012 American Chemical Society. Panels f−h, reprinted from ref 41. Copyright 2014 American Chemical Society.
Figure 2. Nanoaggregate structure of island asphaltene model at concentrations of 20% (a) or 30% (b) and archipelago asphaltene model at concentrations of 20% (c) or 30% (d). The morphology of aggregates is obtained at around 6 × 106 steps (90.5 ns). Each color represents one asphaltene molecule.
typical structures of asphaltene aggregates are observed in simulations, including stacking in a face-to-face structure which is induced by π−π bonding and T-shaped geometry which is greatly influenced by the hexaparticle rings in the asphaltene model molecule. These two types of aggregate morphology have been reported by Zhang et al.24 Mullins et al.40 proposed the Yen−Mullins model. With sufficient concentration, asphaltene molecules usually form nanoaggregates with small (
