Functionalized Silica Nanoparticles within ... - ACS Publications

Subscriber access provided by MONASH UNIVERSITY

Article

Functionalized Silica Nanoparticles Within Multicomponent Oil/Brine Interfaces: A Study in Molecular Dynamics Lucas Stori de Lara, Vagner Alexandre Rigo, and Caetano R. Miranda J. Phys. Chem. C, Just Accepted Manuscript • DOI: 10.1021/acs.jpcc.5b11225 • Publication Date (Web): 15 Mar 2016 Downloaded from http://pubs.acs.org on March 18, 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 33

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

Functionalized Silica Nanoparticles Within Multicomponent Oil/Brine Interfaces: A Study in Molecular Dynamics Lucas S. de Lara1†, Vagner A. Rigo1, ‡ and Caetano R. Miranda*1,ǂ (1)

Centro de Ciências Naturais e Humanas (CCNH) - Universidade Federal do ABC (UFABC) – Santo André – SP – Brazil

ABSTRACT

Stability control of nanoparticles within polar and non-polar liquid interfaces can be influenced by surface effects and molecular-level interactions. This study uses fully atomistic molecular dynamics to investigate the behavior of functionalized silica nanoparticles (NPs) at crude oil/brine interfaces as a function of salt, brine concentration, temperature, and NP surface functionalization. The light crude oil model used in this study comprises aromatics, alkanes and cycloalkanes. Silanized (H-passivated), PEGlyated and sulfonated functionalized NPs are used to account for hydrophilicity variations. The size effect of the functional groups is evaluated for PEGlyated NPs. The highest contact angles (NP moves towards the oil phase) are observed for monovalent (NaCl) solutions, at higher salt concentrations, and for PEGlyated NPs. The findings indicate that the Young-Laplace equation is still valid at nanoscale for spherically symmetric

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

Page 2 of 33

nanoparticles. The mobility of all NPs indicates that the self-diffusion coefficient is ten times faster along than across the interface. The results also show that aromatic molecules adsorb to the NP surface, even on the face located in the brine phase, where they form patchy domains on each nanoparticle evaluated. As this result can greatly affect the stability of NPs at oil/brine interfaces, it has implications for several technological applications.

KEYWORDS: Interface, Nanoparticle, Oil, Water, Molecular Dynamics INTRODUCTION The presence of nanoparticles (NPs) at fluid-fluid interfaces has been subject to increasing interest in recent years,1-4 where the understanding and control of these interfaces are at the core of a broad range of technological applications such as biomedical,5 enhanced oil recovery (EOR) in petroleum fields,6-10 and emulsion and foam stabilization in oil leaks,11 as well as in the food and pharmaceutical industries.1 Silica NPs are of particular interest for advanced EOR processes due to their ability to act as a cheap additive at the oil/brine interface, which can reduce the oil/brine interfacial tension and viscosity,4,6,7,12-14 and stabilize CO2-water foams, as demonstrated experimentally.10 The interfacial tension of the three-phase interface (oil/NP/brine) defines a contact angle between the NP tangent at the interface and the liquid-liquid interface line, as described in previous studies.4,7,15 The oil/NP and brine/NP interfacial tension (and the contact angle) can be modified by using NP surface modifications,4,7,14,15 which can lead to enhanced oil recovery from porous rock. To address these effects, previous studies have considered silica NPs functionalized by molecular groups with different degrees of hydrophilicity,9-14,16 which are candidates to control the oil/brine interface properties and improve the EOR.


2

Page 3 of 33

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 stability, ordering, and dynamics of NPs at interfaces also represent a challenge from a theoretical point of view.17-24 For NPs at interfaces between a non-polar and an aqueous electrolyte media, it has been established that NPs located in the aqueous phase can ionize and form (and subsequently be surrounded by) charged layers of ions. According to the standard Derjaguin, Landau, Verwey, and Overbeek (DLVO) theory for particles within aqueous electrolyte media,1,2 this leads to a repulsive dipole-dipole interaction between NPs, which is the central stabilization mechanism, as discussed in Refs.1,2,23 In an interesting work, Chen et al.23 used optic and atomic force microscopy experiments to show that charged polystyrene latex spheres display patchy domains at a water/air interface. The authors demonstrated that the surface inhomogeneity generates a dipole contribution parallel to the interface, which accounts for an attractive force between NPs. The balance of the parallel and perpendicular components of the NP dipole provides the control of the suspension stability and NP concentration at the interface,22 as observed by Aubry et al.17 Colloidal NPs are expected to have some surface inhomogeneities, which is the case for functionalized and/or amorphous nanoparticles, and particularly for systems that adsorb molecules from the liquid phase to an NP surface. Consistent with these characteristics, crude oil (or other oil emulsions) contains a variety of hydrocarbon components, including alkanes, cycloalkanes and aromatics.8 Previous studies have used molecular dynamics simulations to show that aromatics accumulate on a hydroxylated silica surface25,26 and at a water interface,8,27 primarily due to van der Waals interactions. One study in particular, de Lara et al.,8 employed molecular dynamics simulations using a realistic light crude oil model which contained aromatic, alkane and cycloalkane molecules. The authors showed that these oil components are not equally dispersed when oil is in contact with another liquid phase. Silica and colloidal nanoparticles at an


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

Page 4 of 33

oil/water interface have been studied using molecular simulations4 and experiments.24 However, authors have adopted a homogeneous model for oil (composed solely of decane molecules), which did not include aromatics. A question arises from these previous findings: does the presence of aromatic molecules in the oil phase affect the properties of a brine/silica NP/oil interface? The possible accumulation of aromatics on silica NPs may form patchy domains in those systems. If this effect occurs with NPs at oil/brine interfaces, it can increase the dipole component parallel to the surface, and consequently decrease the NP stability at the oil/brine interface, as discussed in Refs.2,23 Consequently, this may modify the surface tension and contact angle of an NP at oil/brine interfaces. This paper presents a study of silica NPs at oil/brine interfaces, using molecular dynamics simulations. The light crude oil model includes alkane, cycloalkane and aromatic molecules. Overall neutral NPs are used, and deprotonated sites are not included (pH variations are not allowed in the system), while different chemical groups are used for NP functionalization. To the best of our knowledge, this study marks the first time that the accumulation of the different molecular components at an interface and surrounding the NPs has been determined. Wettability phenomena were investigated by using two different methodologies to calculate the contact angle of NPs at the interfaces. Finally, the mobility of the NPs was monitored across and along the interface for each system. METHODS Following is a description of how the molecular models consisting of brine, oil, and nanoparticles were created. Brine is composed of water and NaCl or CaCl2 salts, with the concentration ranging from 0 to 1 wt%. Crude oil is modeled using a combination of aromatics


4

Page 5 of 33

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


(156 toluene and 60 benzene molecules) and alkanes (144 hexane, 132 heptane, 156 octane, 180 nonane, 96 cyclohexane, and 156 cicloheptane molecules), consistent with our previous studies.8,27 In addition to hydroxylated NPs (silanol groups), NPs functionalized by (poly)ethylene and sulphonic (NP-SA) groups are used to account for hydrophilicity variations. These particles are overall neutral in order to assess the accumulation (or not) of oil on silica NPs. The atomistic model of a 3 nm-diameter hydroxylated NP is obtained through the method described in Ref.18 The functionalization is carried out by grafting the silica NP surface with Si(OH)3-(CH2)3-SO3H (SA) and Si(OH)3-(CH2-CH2-O)n-H (PEG) molecular groups, which represent low and high hydrophilic coverages, respectively. For PEGylated nanoparticles, chains with n = 1 and 2 ethylene glycol monomers are evaluated. Further details about the functionalization process and graft density appears in a previous study conducted by the authors of this paper.12 The

molecular

dynamics

calculations

are

undertaken

through

the

Large

Atomic/Molecular Massively Parallel Simulator (LAMMPS) package.28 The hydrocarbons are described by the CHARMM29 interatomic potential, whereas the parameterization of Cruz-Chu et al.30 is used to describe the SiO2 systems, and the use of the SPCE/FH force-field31 for water molecules. The description of ions includes the van der Waals (with a cutoff of 10 Å) and Coulomb interactions. The potential parameters for water, Na, Ca, and Cl were taken from Refs.31,32 The reciprocal space Particle-Particle, Particle-Mesh (PPPM) method28,33 is employed for long-range electrostatic interactions. The liquid systems (brine and crude oil only) were originally set as a cubic box with 10 nm in side under periodic boundary conditions. The box contains a total of 32,000 water molecules, while the salt ions are randomly distributed at concentrations between 0 and 1 wt%.


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

Page 6 of 33

A sequence of NVE, NVT, and NPT calculations are performed in order to obtain the equilibrium densities for each fluid, at a pressure of 1 atm and temperatures of 300 and 350 K. The nanoparticles are inserted into a spherical void at the oil/brine interface. The concentration of NPs is estimated to be ~1×106 particles/µm2 at the interfaces. For each system the simulations run along 1.0 ps in NVE, 10.0 ps in NVT, and 2.0 ns in the NPT ensemble, using a timestep of 0.5 fs. After equilibration, production occurs in NVT along 8 ns. The Nosé-Hoover thermostat34,35 and Andersen barostat36 are applied to control the temperature and pressure, respectively. The density profiles for each molecular component are determined by considering the linear density dl (z) perpendicular to the oil/brine interface (z-axis), and the radial density profile dr (r) with respect to the center of mass of the nanoparticles. The oil/brine interface is defined by the Gibbs dividing surface (GDS), consistent with the authors’ previous study.8 To quantify the amount of aromatic molecules accumulated on the surfaces of the nanoparticles, their contribution to radial density components dr (r) is monitored separately at the oil and brine phases. The interfacial tension for spherical fluids/nanoparticles interfaces is calculated by means of the pressure tensor in the Gibbs formulation:27-39  (P − P )2 ∞ dP (r )  γ = − l g ∫ r 3 N dr  dr 8   0

1/ 3

,

(1)

where Pl and Pg represent the internal and external pressures acting on the nanoparticle, and PN (r) is the normal component of pressure, given by:


6

Page 7 of 33

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


2 PN ( r ) = 2 r

∞

∫ rP (r )dr ,

(2)

T

0

where PT (r) is the tangential component of the pressure tensor. For the oil/brine interface, the interfacial tension is calculated as: L

1 a γ = ∫ dz ( p ab ( z ) − pT ( z )) , 2 − Lb

(3)

where the distances La and Lb are the limits of the interfacial region, and pT (r) and pab (r) are the tangential and normal components of the pressure tensor, respectively. The oil/NP/brine contact angle is determined by two different methods. The first adopts a modified version of the Young-Laplace equation in order to obtain the contact angle from the interfacial tensions, as described by Pieranski21 for a spherical particle:40

cos θ γ =

γ o − np − γ b − np γ o −b

,

(4)

where γo,b-np represents the interfacial tension (eq 1 and 3). The labels o, b, and np are the oil, brine, and nanoparticle values, respectively. The second method is based on a direct measurement of the contact angle by considering the average NP position in the interface. The angle is obtained from the distance between the nanoparticle center of mass and the interface (l) and the NP radii (R), as described in Ref.4,7,15

cos θ e =

l R

(5)

The self-diffusion coefficient for evaluated NPs is obtained by means of mean square displacement (MSD) in the Einstein equation:


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

1 1 N MSD (t ) D= lim ∑ , t 2d t →∞ N i =1

Page 8 of 33

(6)

where d represents the dimension of the system. The interface is set up to be perpendicular to the z direction, allowing the total diffusion (D) and its components along (parallel to) (Dxy) and across (Dz) the interface to be calculated. RESULTS Figures 1a-d present snapshot images of the four evaluated silica NPs in the oil/brine interfaces. The analysis of these systems begins with the density profiles across the interface (linear density) and around the nanoparticle (radial density), followed by a calculation of the interfacial tension and contact angle. Finally, the transport properties of nanoparticles in the interfaces are presented.


8

Page 9 of 33

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. Snapshot images of (a) NP-H, (b) NP-SA, (c) NP-EG, and (d) NP-PEG2 nanoparticles within the oil/brine interface, for a brine solution composed of 1 wt% of NaCl, at 300 K and 1 atm. The aqueous phase is on the right of each image. Density Profiles of oil/nanoparticles/brine interfaces Figures 2a-d present the linear density profile dl (z) for NP-H, NP-SA, NP-EG, and NPPEG2 in the oil/brine (NaCl) interfaces, respectively. For all systems, the nanoparticle center of mass is located within the aqueous phase. A similar situation is observed for a NP-H in the interface between an ionic liquid and water.41 However, for NP-SA and NP-PEG2, the average nanoparticle position is found to be located closer to the interface than NP-H and NP-EG. This


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 33

information is central to understand the robustness of the nanoparticle suspension, principally because nanoparticles are subject to Coulombic and van der Waals interactions through the brine, but primarily a van der Waals interaction through the oil. The dipole-dipole interaction between NPs can be enhanced by taking the NPs even closer within the brine phase. In Figure 2a-d, a structuring of the ions density (layers of charge) in the aqueous phase is observed for all interfaces, even for overall neutral NPs, as reported by de Lara et al.13 for silica NPs within an electrolyte solution. The dl (z) results show an oil accumulation at the brine-oil interface for all nanoparticles studied. Figures 2a-d demonstrate that this accumulation is primarily due to aromatics, and cycloalkanes to a lesser degree. In the absence of nanoparticles, this accumulation of aromatics is observed in the brine-oil interface, as previously reported results.8 For NP-H, NP-EG, and NPPEG2 (Figure 2a, c, and d, respectively) the results indicate a weak accumulation of aromatic molecules on the nanoparticle surface. However, this does not occur with NP-SA (Figure 2b). A similar trend is found for cycloalkanes.


10

Page 11 of 33

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. Linear density profile (dl) for (a) NP-H, (b) NP-SA, (c) NP-EG, and (d) NP-PEG2 within the oil/brine interface, with the brine solution composed of 1 wt% of NaCl for 300 K and 1 atm. The right-hand scale of the dl-axis refers to the nanoparticle and ions (Na and Cl) densities.


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

Page 12 of 33

Figure 3a-d present the dr (r) curves for NP-H, NP-SA, NP-EG, and NP-PEG2, respectively. The results show that a greater number of aromatics is concentrated on the NP-H (Figure 3a) and NP-EG (Figure 3c) surfaces, compared with NP-SA (Figure 3b) and NP-PEG2 (Figure 3d). For NP-SA, the dr (r) peak of the aromatics indicates the presence of these molecules within the shell NP region. This is due to the weak hydrogen bonding between the aromatic rings and the silanol groups on the nanoparticle surface,25,26 and the functional group spacer. Additionally, the dr (r) peak of aromatics for NP-H (Figure 3a) is clearly located on the silanol groups. The same results occur for nanoparticles within the oil/CaCl2 interface, as presented in Figure S1 in the supplemental material. Interestingly, Figure 2a-d also show that an oil accumulation can occur on a nanoparticle surface, even on the nanoparticle side which is within the brine. This result is not found to the same degree for NP-H in an interface between hexane and an ionic liquid.41 To better analyze these results, the time average position of aromatics is presented in Figure 4a-d, for NP-H, NPSA, NP-EG, and NP-PEG2, respectively. The accumulation of aromatics on a nanoparticle surface, even within the brine, is qualitatively observed for those systems. Added to this, the accumulation is not symmetric. These aromatic patchy domains can have considerable effects on the stability of those systems,1,2,24 and must be taken into account when NPs are at an interface in the presence of aromatic molecules. Based on experiments conducted by Chen et al.,24 charged polystyrene latex spheres display patchy domains at the water/air interface, with drastic consequences for NP stability at the interfaces. The present study reports the first time that a similar process has been observed for fluid (oil)-particle intermolecular interactions.


12

Page 13 of 33

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. Radial density profile (dr) for (a) NP-H, (b) NP-SA, (c) NP-EG, and (d) NP-PEG2 in the oil/brine interface, with the brine solution composed of 1 wt% of NaCl for 300 K and 1 atm. The right-hand scale of the dr-axis refers to the nanoparticle and ions (Na and Cl) densities.


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 33

Figure 4. Density map of the time averaged distribution of aromatics molecules (gray surfaces) for (a) NP-H, (b) NP-SA, (c) NP-EG, and (d) NP-PEG2, with the brine solution composed of 1 wt% of NaCl. The panel on the left shows the perpendicular (z-axis), while the parallel (xyplane) view along the interface appears in panel on the right. The isosurfaces represent values of 0.18 g/cm3.

To quantify the accumulation of aromatics on the surface of the NPs at both the oil and brine phases, the components of dr (r) are presented separately in Figure 5a-d for NP-H, NP-SA, NP-EG, and NP-PEG2, respectively, for a brine solution consisting of 1 wt% NaCl. Comparing all systems, while the NP-H presents the largest accumulation of aromatics in the oil phase, the NP-SA shows the smallest accumulation within the brine phase.


14

Page 15 of 33

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 5b shows that aromatic molecules may penetrate the region of the functional groups for NP-SA. For the same NP models, a similar effect was also observed recently for water molecules.12 Interestingly, NP-EG and NP-PEG2 both accumulate a similar quantity of aromatics for both regions (in brine and oil). Compared with the results for NaCl (Figure 5), the aromatic accumulation on the surface of the NPs is not significantly affected within a CaCl2 brine for each fluid phase of an oil/CaCl2 brine interface, as demonstrated in Figure S3 of the supplemental material.

Figure 5. Components of the radial density profile of aromatics for (a) NP-H, (b) NP-SA, (c) NP-EG, and (d) NP-PEG2 in the oil/brine interface, with the brine solution composed of 1 wt% of NaCl for 300 K and 1 atm. Brine phase and oil phase each refer only to the radial density


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 33

profile of aromatics surrounding the NP face within the brine and oil. The average NP profile is indicated in each case. Table 1 presents the mass of aromatics near the NP surface (up to 28 Å from the NP center of mass) at each phase (oil and brine). For the NP-H at NaCl brine, the total amount of aromatics increases by approximately 29% when compared to the aromatics in the oil model without nanoparticles. The greatest difference in the mass of aromatics between the oil and brine phases (5,205.2 u) appears with the NP-SA at NaCl brine. Results for CaCl2 show similar behavior.

Table 1. Mass of aromatics surrounding the nanoparticle systems for oil/brine at 300 K, up to 28 Å from the nanoparticle center of mass. The brine and oil phases only represent the mass on the nanoparticle face within brine and oil. Total refers to the sum of both hemisphere components.

Nanoparticle NP-H NP-SA NP-EG NP-PEG2 NP-H NP-SA NP-EG NP-PEG2

NaCl brine Brine Phase (u) Oil Phase (u) 4,319.1 8,527.5 2,214.9 7,420.1 4,596.0 5,814.2 4,651.4 5,205.1 CaCl2 brine 3,986.9 8,195.3 1,993.4 7,087.8 4,540.6 5,648.1 4,706.8 6,035.7

Total (u) 12,846.6 9,635.0 10,410.2 9,856.5 12,182.2 9,081.2 10,188.7 10,742.5

The effect of divalent salts (CaCl2) on dl (z) curves is evaluated for NP-H, NP-SA, NPEG, and NP-PEG2 systems, as shown in Figure 5a-d, respectively. The results for all coverages indicate that the nanoparticle center of mass is even closer to the interface, compared with results for NaCl brine. The distribution of ions in brine also indicates a layered structure (similar to


16

Page 17 of 33

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


NaCl), as previously reported.13 With respect to NP-H (Figure 5a), a predominant accumulation of Cl- (anions) close to the NP surface is observed. This large accumulation of anions on NP-H, when compared to cations, is also reported by Frost et al.41 in a molecular dynamics investigation of NPs at an ionic liquid and water or hexane interface. With respect to oil partition, aromatics accumulate on a nanoparticle surface, particularly for NP-H, NP-EG, and NP-PEG2 (Figure 5a, c, and d, respectively). For NP-SA (Figure 5b), the accumulation of aromatics take place close to the interface (z ~0 Å), but reaches a minimum at the nanoparticle surface (z ~ -20 Å). This indicates that the SA groups do not accumulate aromatics, while the peak at ~0 Å is primarily due to accumulation on the water interface and silanol groups at the surface. The time average position of aromatics for NP-H, NP-SA, NP-EG, and NP-PEG2 in the CaCl2 brine is reported in Figure S2 in the supplemental material, where the accumulation of aromatics on NPs can be confirmed, principally in the case of NP-H.

Interfacial tension and contact angles at oil/nanoparticles/brine interfaces The interfacial tension profiles are evaluated for all nanoparticle coverages, and brine solutions composed of 0, 0.25, 0.50 and 1 wt% of NaCl, at 300 and 350 K. Compared to the brine/oil interface without NPs8 (under the same thermodynamic conditions), the inclusion of NPs promotes a reduction of interfacial tension (the values appear in Table S1 in the supplementary material). Using the interfacial tension data, the contact angle was obtained by means of eq 4, and using geometrical data, through eq 7. The contact angles from eq 4 and eq 7 are named θγ and θe, respectively. If the nanoparticle center of mass is located in the brine phase, the contact angle is less than 90º with the interface, whereas if the center of mass is inside the oil, the contact angle is greater than 90º.


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 33

Figure 5. Linear density profile (dl) for the (a) NP-H, (b) NP-SA, (c) NP-EG, and (d) NP-PEG2 in the oil/brine interface, for a brine solution composed of 1 wt% of CaCl2 for 300 K and 1 atm. The right-hand scale of the dl-axis refers to the nanoparticle and ions (Ca and Cl) densities. Figure 6a and b present the θγ results for the NaCl and CaCl2 solutions, respectively. For both salts, the θγ increases with increasing salt concentration. A comparison of θγ for NP-EG and


18

Page 19 of 33

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


NP-PEG2 reveals that including the ethylene monomers results in an increase of θγ for NP-PEG2 over NP-EG. The θe values confirm these general trends, as presented in Figure 7a and b for NaCl and CaCl2, respectively. This effect can be explained as a result of the hydrophobicity of the spacer in the ethylene glycol group. As the functional group increases from NP-EG to NPPEG2, the increasing of the spacer on NP-PEG2 leads to an increase in the contact angle for this coverage (where the NP moves toward the oil phase). The numerical values of each surface tension component and the contact angles for NPs in the NaCl and CaCl2 solutions are available as supplementary material in Tables S1 and S2. The θe and θγ values show a reasonable agreement (Figure 6a, b and 7a, b). Figure 8a and b present the θe against θγ data for both NaCl and CaCl2, and show a direct comparison between the Young-Laplace and geometric methods, respectively. As indicated by the ideal fitting (y = x), the results from both methods are consistent for the spherical silica NPs in this study.


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 33

Figure 6. Contact angles in oil/NP/brine interfaces from the interfacial tension values (YoungLaplace equation) for a brine solution composed of (a) NaCl and (b) CaCl2 salt. The θγ values are obtained by averaging a total of 2,000 configurations along each simulation at 300 K (circles) and 350 K (squares) and 1 atm. The margin of error for θγ is within ± 0.1º.


20

Page 21 of 33

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 7. Direct measurement of the contact angle in oil/NP/brine interfaces (according to eq 7) for a brine solution composed of (a) NaCl and (b) CaCl2 salts. The θe values are obtained by averaging a total of 2,000 configurations along each simulation at 300 (circles) and 350 K (squares) and 1 atm. The margin of error for θe is within ± 2º.


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 33

Figure 8. Linear correlation of the comparison between oil/NP/brine contact angles (θγ and θe) for all results of (a) NaCl and (b) CaCl2. The solid line in each represents an eye guide for θγ = θe.

Transport properties Figure 9 shows the diffusion coefficients obtained for each NP in the oil/brine interfaces. The total self-diffusion coefficients (D) for nanoparticles within brine solutions composed of NaCl and CaCl2 appear in Figure 9a and d, respectively. With increasing salt concentration a similar trend is seen for all nanoparticle coverages and the salts evaluated. When different salts are used, diffusion is lower for NPs in the oil-brine interface with divalent ions (CaCl2) than for the monovalent interface (NaCl). Another material result is that the D values are dependent of the nanoparticle functional groups. In particular, the NP-H presents the lowest D values among


22

Page 23 of 33

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


all salt concentrations evaluated. This effect is related to the lower contact angle for NP-H (most of the nanoparticles are in the water phase), as well as more significant aromatic accumulation on the nanoparticle surface. In a previous study, Frost et al.41 analyzed silanized silica nanoparticles within an interface between an ionic liquid and hexane (without aromatics). The authors showed that the mobility of an ionic liquid is reduced by the presence of an NP-H nanoparticle, indicating an interaction of silanol sites with the ionic liquid. Here, for the first time the adsorption of hydrocarbons on an NP-H surface appear as an additional contribution to reduce the diffusion of NP-H at the interface. Particles at interfaces have distinct self-diffusion values along (Dxy) and across (Dz) the interface plane,24,42 as shown in Figure 9b and c, and Figure 9e and f for NaCl and CaCl2, respectively. The increase in D values (Figsure 9a and d), as a function of the salt concentration, is due in great part to the in-plane contribution to diffusion (Dxy). Additionally, the ratio between the two diffusion components is ~10 times and remains almost constant with increasing salt concentration. A major difference between the in-plane and perpendicular diffusion has also been reported in other studies, such as in experiments for polystyrene particles at a water-decane interface,24 and theoretically for much simpler Lennard-Jones models of NPs within a fluid-fluid interface.42


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

Page 24 of 33

Figure 9. Total (a), xy (b), and z (c) diffusion coefficients for nanoparticles in oil/NP/brine interfaces with a brine solution composed of concentrations for NaCl. The related results for


24

Page 25 of 33

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


CaCl2 are in (d), (e), and (f), respectively. Circles and squares represent data at 300 and 350 K, respectively.

CONCLUSIONS In this study the authors use molecular dynamics simulations to analyze the density profile, diffusion, and contact angle between silica nanoparticles (NPs) in an oil/brine interface under different NP surface coverage, temperature, and salt conditions for the brine solution. For the first time a multicomponent model for light crude oil (including aromatics and alkanes) has been employed to study this brine/NP/oil interface. Results show that the highest contact angles occur at higher salt concentrations. The NP-PEG2 coverage gives the highest contact angle between the NPs evaluated, as it is displaced towards the interface and has greater penetration in the oil phase, making it a potential candidate for EOR applications (where these NPs can act as a chemical agent to modify the interface properties). Additionally, by determining the interface tension using two distinct methods, our results indicate that the Young-Laplace equation continues to be valid at a nanoscale for spherically symmetric nanoparticles. The density profile across interfaces shows that aromatic molecules accumulate on the NP surface, even in the brine phase. This effect occurs for all nanoparticles evaluated, due primarily to silanol groups on the NP surface. For the NP-SA, no aromatic accumulation is observed near the -SO3H head groups. These results indicate the formation of patchy domains of oil on the NP surfaces. This can lead to a reduction in NP stability within the oil/NP/brine interfaces. The effects on the diffusion coefficient of NPs at the interfaces show that the diffusion along the interfaces is 10 times greater than diffusion across the interface. This result is seen for


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

Page 26 of 33

all the systems evaluated, similar to those reported for suspensions of particles at low densities.24,42

AUTHOR INFORMATION

Corresponding Author *Corresponding author: C.R. Miranda, +55 11 30917009, e-mail: [email protected]

Present Addresses †Departamento de Física, Universidade Estadual de Ponta Grossa (UEPG) – Ponta Grossa – PR – Brazil ‡Universidade Tecnológica Federal do Paraná (UTFPR) – Campus Cornélio Procópio – PR – Brazil ǂInstituto de Física, Universidade de São Paulo (IF-USP) – São Paulo – SP – Brazil

Author Contributions This paper is the result of contributions from all the authors, each of whom has given their approval to the final version.

ACKNOWLEDGMENT


26

Page 27 of 33

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


We acknowledge the financial support of the Advanced Energy Consortium (AEC) and the Brazilian agencies CAPES, FAPESP and CNPq. The calculations have been partially performed at CENAPAD-SP, CESUP-RS, UFABC and UTFPR-Cornélio Procópio supercomputer facilities.

ASSOCIATED CONTENT

Supplementary Information Available. The components of the interface tension and contact angles for all NPs, as well as the radial density profile and the time averaged density map for NPs within CaCl2 solutions, are included as supplementary information. This information is available free of charge via the Internet at http://pubs.acs.org

REFERENCES (1) Hunter, T. N.; Pugh, R. J.; Franks, G. V.; Jameson, G. J. The Role of Particles in Stabilising Foams and Emulsions. Adv. Colloid Int. Sci. 2008, 137, 57-81. (2) Bresme, F.; Oettel, M. Nanoparticles at Fluid Interfaces. J. Phys.: Condens. Matter 2007, 19, 413101. (3) Binks, B. P. Colloidal Particles at Liquid Interfaces. Phys. Chem. Chem. Phys. 2007, 9, 62986299.


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

Page 28 of 33

(4) Fan, H.; Resasco, D. E.; Striolo, A. Amphiphilic Silica Nanoparticles at the Decane-Water Interface: Insights from Atomistic Simulations. Langmuir 2011, 27, 5264-5274. (5) Oroskar, P. A.; Jameson, C. J.; Murad, S. Surface-Functionalized Nanoparticle Permeation Triggers Lipid Displacement and Water and Ion Leakage. Langmuir 2015, 31, 1074-1085. (6) Sakthivela, S.; Velusamya, S.; Gardasb, R. L.; Sangwaia, J. S. Adsorption of Aliphatic Ionic Liquids at Low Waxy Crude Oil–Water Interfaces and the Effect of Brine. Colloids and Surfaces A: Physicochem. Eng. Aspects 2015, 468, 62-75. (7) Garbin, V.; Crocker, J. C.; Stebe, K. J. Nanoparticles at Fluid Interfaces: Exploiting Capping Ligands to Control Adsorption, Stability and Synamics. J. of Colloid and Int. Sci. 2012, 387, 111. (8) de Lara, L. S.; Michelon, M. F.; Miranda,C. R. Molecular Dynamics Studies of Fluid/Oil Interfaces for Improved Oil Recovery Processes. J. Phys. Chem. B 2012, 116, 14667-14676. (9) de Lara, L. S.; Michelon, M. F.; Metin, C. O.; Nguyen, Q. P.; Miranda, C. R. Interface Tension of Silica Hydroxylated Nanoparticle with Brine: A Combined Experimental and Molecular Dynamics Study. J. Chem. Phys. 2012, 136, 164702. (10) Worthen, A. J.; Bagaria, H. G.; Chen, Y. S.; Bryant, S. L.; Huh, C.; Johnston, K. P. Nanoparticle-Stabilized Carbon Dioxide-in-Water Foams with Fine Texture. J. of Colloid and Int. Sci. 2013, 391, 142-151. (11) Worthen, A. J.; Foster, L. M.; Dong, J.; Bollinger, J. A.; Peterman, A. H.; Pastora, L. E.; Bryant, S. L.; Truskett, T. M.; Bielawski, C. W.; Johnston, K. P. Synergistic Formation and


28

Page 29 of 33

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


Stabilization of Oil-in-Water Emulsions by a Weakly Interacting Mixture of Zwitterionic Surfactant and Silica Nanoparticles. Langmuir 2014, 30, 984-994. (12) Rigo, V. A.; de Lara, L. S.; Miranda, C. R. Energetics of Formation and Hydration of Functionalized Silica Nanoparticles: An Atomistic Computational Study. Appl. Sur. Sci. 2014, 292, 742-749. (13) de Lara, L. S.; Rigo, V. A.; Michelon, M. F.; Metin, C. O.; Nguyen, Q. P.; Miranda, C. R. Molecular Dynamics Studies of Aqueous Silica Nanoparticle Dispersions: Salt Effects on the Double Layer Formation. J. Phys.: Condens. Matter. 2015, 27, 325101. (14) de Lara, L. S.; Rigo, V. A.; Miranda, C. R. The Stability and Interfacial Properties of Functionalized Silica Nanoparticles Dispersed in Brine Studied by Molecular Dynamics. Eur. Phys. J. B 2015, 88, 261. (15) Reguera, J.; Ponomarev, E.; Geue,T.; Stellacci, F.; Bresme, F.; Moglianetti, M. Contact Angle and Adsorption Energies of Nanoparticles at the Air–Liquid Interface Determined by Neutron Reflectivity and Molecular Dynamics. Nanoscale 2015, 7, 5665-5673. (16) Metin, C. O.; Baran Jr., J. R.; Nguyen, Q. P. Adsorption of Surface Functionalized Silica Nanoparticles Onto Mineral Surfaces and Decane/Water Interface. J Nanopart. Res. 2012, 14, 1246. (17) Aubry N.; Singh, P. Physics Underlying Controlled Self-Assembly of Micro- and Nanoparticles at a Two-Fluid Interface Using an Electric Field. Phys. Rev. E, 2008, 77, 056302. (18) Singh, P.; Joseph, D. D.; Fisher, I. S.; Dalal, B. Role of Particle Inertia in Adsorption at Fluid-Liquid Interfaces. Phys. Rev. E, 2011, 83, 041606.


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

Page 30 of 33

(19) Cheung, D. L.; Bon, S. A. F. Interaction of Nanoparticles with Ideal Liquid-Liquid Interfaces. Phys. Rev. Lett. 2009, 102, 066103. (20) Garbin, V.; Jenkins, I.; Sinno, T.; Crocker, J. C.; Stebe, K. J. Interactions and Stress Relaxation in Monolayers of Soft Nanoparticles at Fluid-Fluid Interfaces. Phys. Rev. Lett. 2015, 114, 108301. (21) Pieranski, P. Two-Dimensional Interfacial Colloidal Crystals. Phys. Rev. Lett. 1980, 45, 569. (22) Oettel, M.; Domínguez, A.; Dietrich, S. Effective Capillary Interaction of Spherical Particles at Fluid Interfaces. Phys. Rev. E, 2005, 71, 051401. (23) Chen, W.; Tan, S.; Huang, Z.; T.-K. Ng,; Ford, W. T.; Tong, P. Measured Long-Ranged Attractive Interaction Between Charged Polystyrene Latex Spheres at a Water-Air Interface. Phys. Rev. E, 2006, 74, 021406. (24) Parolini, L.; Law, A. D.; Maestro, A.; Buzza, D. M. A.; Cicuta, P. Interaction Between Colloidal Particles on an Oil–Water Interface in Dilute and Dense Phases. J. Phys.: Condens. Matter 2015, 27, 194119. (25) Ledyastuti, M.; Liang, Y.; Kunieda, M.; Matsuoka, T. Asymmetric Orientation of Toluene Molecules at Oil-Silica Interfaces. J. Chem. Phys. 2012, 137, 064703. (26) Coasne, B.; Fourkas, J. T. Structure and Dynamics of Benzene Confined in Silica Nanopores. J. Phys. Chem. C 2011, 115, 15471-15479.


30

Page 31 of 33

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


(27) 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. (28) Thompson, A. P.; Plimpton, S. J.; Mattson, W. General Formulation of Pressure and Stress Tensor for Arbitrary Many-Body Interaction Potentials Under Periodic Boundary Conditions. J. Chem. Phys. 2009, 131, 154107. (29) Brooks, B. R.; Brooks III, C. L.; Mackerell Jr., A. D.; Nilsson, L.; Petrella, R. J.; Roux, B.; Won, Y.; Archontis, G.; Bartels, C.; Boresch, S. et al. CHARMM: The Biomolecular Simulation Program. J. Comp. Chem. 2009, 30, 1545-1614. (30) Cruz-Chu, E. R.; Aksimentiev, A.; Schulten, K. Water–Silica Force Field for Simulating Nanodevices. J. Phys. Chem. B 2006, 110, 21497-21508. (31) Alejandre, J.; Chapela, G. A.; Bresme, F.; Hansen, J.-P. The Short Range Anion-H Interaction is the Driving Force for Ccrystal Formation of Ions in Water. J. Chem. Phys. 2009, 130, 174505. (32) Gavryushov, S. Dielectric Saturation of the Ion Hydration Shell and Interaction between Two Double Helices of DNA in Mono- and Multivalent Electrolyte Solutions: Foundations of the ε-Modified Poisson-Boltzmann Theory. J. Phys. Chem. B 2007, 111, 5264. (33) Hockney, R. W.; Eastwood, J. W. Computer Simulation using Particles, Adam Hilger, New York, 1989. (34) Nosé, S. A Unified Formulation of the Constant Temperature Molecular Dynamics Methods. J. Chem. Phys. 1984, 81, 511.


31


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 32 of 33

(35) Hoover, W.G. Canonical dynamics: Equilibrium Phase-Space Distributions. Phys. Rev. A

1985, 31, 1695-1697. (36) Andersen, H. C. Molecular Dynamics Simulations at Constant Pressure and/or Temperature. J. Chem. Phys. 1980, 72, 2384. (37) Schweigert, I. V.; Lehtinen, K. E. J.; Carrier, M. J.; Zachariah, M. R. Structure and Properties of Silica Nanoclusters at High Temperatures. Phys. Rev. B 2002, 65, 235410. (38) Hawa, T.; Zachariah, M. R. Internal Pressure and Surface Tension of Bare and Hydrogen Coated Silicon Nanoparticles. J. Chem. Phys. 2004, 121, 9043. (39) Thompson, S. M.; Gubbins, K. E.; Walton, J. P. R. B.; Chantry, R. A. R.; Rawlinson, J. S. A Molecular Dynamics Study of Liquid Drops. J. Chem. Phys. 1984, 81, 530. (40) Udayana Ranatunga, R. J. K.; Kalescky, R. J. B.; Chiu, C.-C.; Nielsen, S. O. Molecular Dynamics Simulations of Surfactant Functionalized Nanoparticles in the Vicinity of an Oil/Water Interface. J. Phys. Chem. C 2010, 114, 12151-12157. (41) Frost, D. S.; Machas, M.; Dai, L. L. Molecular Dynamics Studies on the Adaptability of an Ionic Liquid in the Extraction of Solid Nanoparticles. Langmuir 2012, 28, 13924-13932. (42) Cheung, D. L. Molecular Simulation of Nanoparticle Diffusion at Fluid Interfaces. Chem. Phys. Lett. 2010, 495, 55-59.


32

Page 33 of 33

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


Table of Contents (TOC) Image


33

Functionalized Silica Nanoparticles within ... - ACS Publications

Recommend Documents