Composition Analysis and Viscosity Prediction of Complex Fuel

Department of Chemical Engineering, Imperial College London, South Kensington Campus, London SW7 2AZ, United Kingdom. § IFP Energies Nouvelles ...
9 downloads 0 Views 2MB Size
Article pubs.acs.org/EF

Composition Analysis and Viscosity Prediction of Complex Fuel Mixtures Using a Molecular-Based Approach Michelle Aquing,†,∥ Fausto Ciotta,‡ Benoit Creton,† Christophe Féjean,† Annabelle Pina,† Cyril Dartiguelongue,§ J. P. Martin Trusler,‡ Romain Vignais,‡ Rafael Lugo,† Philippe Ungerer,†,⊥ and Carlos Nieto-Draghi*,† †

IFP Energies Nouvelles, 1 et 4, Avenue de Bois-Préau, 92852 Rueil-Malmaison, France Department of Chemical Engineering, Imperial College London, South Kensington Campus, London SW7 2AZ, United Kingdom § IFP Energies Nouvelles, Rond point de l’échangeur de Solaize, BP3, 69360 Solaize, France ‡

S Supporting Information *

ABSTRACT: The automobile industry currently faces the challenge of developing a new generation of diesel motor engines that satisfy both increasingly stringent emission regulations and reduces specific fuel consumption. The performance of diesel engines, seen in terms of emissions and specific fuel consumption, generally improves with increasing fuel-injection pressure. The design of the next generation of diesel fuel injection systems requires the knowledge of the thermophysical properties, in particular viscosity, of a wide-type of diesel fuels at pressures up to 300 MPa or more. The objective of the present work is to demonstrate that it is possible to predict the viscosity of any petroleum-based diesel fuel, using, exclusively, its molar fraction distribution as provided by multidimensional gas chromatography techniques. The precise knowledge of the fuel chemical constituents allows the understanding of the influence of the different hydrocarbon families on the fluid viscosity by means of molecular dynamics simulations. The accuracy of the Anisotropic United Atom force-field was tested and was found to be in agreement with experimental viscosities obtained with a new vibrating wire device at different temperatures and pressures up to 300 MPa. Finally, the experimental and simulated viscosities have been compared with improved group contribution method that has been coupled with gas chromatography experimental measurements for a viscosity prediction that was estimated to be of less than 18% of mean absolute deviation.

1. INTRODUCTION A number of factors motivate the development of a new generation of high-performance diesel engines for automotive applications. These factors include stricter regulations limiting toxic and noxious vehicle emissions, the need to reduce specific fuel consumption, and the increasing diversity of fuels resulting from changes in petroleum refining practice and the widespread availability of biodiesel. Accordingly, there are increased efforts to optimize direct-injection diesel engine performance, primarily through the improved atomization of the fuel that results from higher injection pressures. The degree of droplet breakup inside fuel injectors is strongly related to the viscosity of the fuels; that is, the higher the viscosity, the bigger the size of the droplets are.1 The size of droplets is then strongly related to the pollutant emission, since small droplets favor vaporization and combustion. These efforts require better understanding of the processes involved in the injection, atomization, and combustion of fuels. The improvement in diesel engine modeling certainly requires detailed chemical kinetics models and reaction mechanisms for the different chemical species involved. 2 Such efforts require accurate knowledge of thermophysical properties, in particular, the viscosity and its dependence upon pressure, temperature, and fuel composition. Experimental studies of fuel properties under extreme injection conditions are expensive and time-consuming, and alternatives are, therefore, of interest. In particular, modeling approaches are sought that are capable of predicting the viscosity, and other © 2012 American Chemical Society

key thermophysical properties, at extreme conditions on the basis of the chemical composition of the fuel. Traditional correlative models for viscosity, mostly based on reference equations,3−8 corresponding states (based on general physical constants such as critical properties and acentric factor),9−11 semiempirical,12,13 hard spheres,14,15 free-volume,16−19 or reaction rate theory20,14 are not suited to this application because they have not been validated for high pressures and are not easily applicable for multicomponent systems such as the one treated in this work. Herein, we describe the approach followed to predict the viscosity and the rheological behavior of different diesel fuels under the relevant conditions of temperature and pressure. We show that the chemical composition can be determined through a new multidimensional gas chromatography technique that is able to identify and quantify, by carbon number and chemical family, in a single analysis all the significant chemical species present in a fuel. The chemical composition so determined is then used to estimate the viscosity through two different approaches: (i) molecular dynamics simulations and (ii) an optimized group contribution (GC) method.21−24 To validate our approach, a special measurement device was developed to obtain experimentally the density and viscosity for different diesel Received: January 18, 2012 Revised: March 15, 2012 Published: March 20, 2012 2220

dx.doi.org/10.1021/ef300106z | Energy Fuels 2012, 26, 2220−2230

Energy & Fuels

Article

Figure 1. Schematic representation of the multidimensional gas chromatography device: (a) The sample is separated by a modulator, which traps, refocuses, and reinjects the effluent from the first column into the second column. The combination of a first nonpolar column with a second (semi)polar column allows a separation of components according to their boiling points on the first dimension and according to their polarity on the second dimension in a single analysis. (b) The resulting polarity versus volatility separation leads to structured chromatograms where hydrocarbons are arranged according to their chemical group and to their number of carbon atoms (illustration: 2D chromatogram obtained for a diesel cut).

diaromatics (naphthalenes, naphtheno-diaromatics), and polyaromatic compounds (phenanthrenes, pyrenes). Within each family, the distribution of isomers according to carbon number is obtained. Although other techniques, such as advanceddistillation curve (ADC), are available for fuel characterization, these are less adapted for heavy fluids having hydrocarbons with carbon numbers over 30.26 For the analysis of samples containing polynaphthenes and olefins, a specific prefractionation technique was developed using supercritical fluid chromatography (SFC) linked to twin 2D gas chromatography systems. Through this powerful technique, saturated and unsaturated compounds are first separated by SFC and transferred on two different 2D gas phase chromatography columns sets (twin-two-dimensional gas phase chromatography) placed in the same oven. In this way, polynaphthenes and olefins are specifically quantified by families and by groups of isomers.27 2.2. Density Measurements. At atmospheric pressure, they have been measured with an oscillating U-tube density meter (DMA 4500, Anton Paar) with a method based on NF EN ISO 12185. As specified by Anton Paar society, the calibration is done at 20 °C with air and water and checked at the temperature of measurement. The precision of this method obtained by statistical examination of tests results is as follows in the case of fuels: repeatability (the difference between

samples at temperatures from 273 to 423 K and at pressures from 0.1 to 300 MPa. We show that our methodology is able to predict the viscosity, both at equilibrium and under shear, of any diesel sample under different thermodynamic conditions solely on the basis of the detailed characterization of the fluid obtained by gas chromatography.

2. METHODS 2.1. Multidimensional Gas Chromatography. Comprehensive two-dimensional (2D) gas phase chromatography consists in the association of two capillary columns with different polarities, separated by a modulator, which traps, refocuses, and reinjects the effluent from the first column into the second column. The combination of a first nonpolar column with a second (semi)-polar column allows a separation of components according to their boiling points on the first dimension and according to their polarity on the second dimension, as illustrated by Figure 1. 2D chromatograms are integrated using dedicated software developed by IFPEN (2DChrom). 2D gas phase chromatography hyphenated to a flame ionization detector (FID) allows quantifying about 250 groups of compounds25 in petroleum middle distillates (kerosenes and diesels). Chemical families identified include linear and branched paraffins, naphthenes (1 and 2 rings), monoaromatics (alkylbenzenes, indanes and tetralins, indenes), 2221

dx.doi.org/10.1021/ef300106z | Energy Fuels 2012, 26, 2220−2230

Energy & Fuels

Article

Figure 2. Schematic representation of the approach followed in this work to directly link fuel viscosity with its molecular composition.

Diesel densities and viscosities were computed by means of equilibrium molecular dynamics (EMD) as well as nonequilibrium molecular dynamics (NEMD) simulations using the Newton code.31 Dynamic viscosities have been computed by EMD simulations using the Einstein formula, whereas the shear dependent viscosities were computed using a modified version of the reverse nonequilibrium molecular dynamics (RNEDM) simulation of Müller-Plathe.32,33 In all simulations the temperature and pressure applied were constrained by the Berendsen weak coupling bath.34 The equations of motion were integrated through the velocity Verlet algorithm with constrained bonds using the Rattle algorithm35 using a time step of 2 fs (1 fs = 10−15 s) and a cutoff radius of 0.12 nm. All EMD simulations were performed with 500 molecules in cubic boxes with periodic boundary conditions. The number of molecules for each species corresponds to a surrogate model obtained form the gas chromatography experiments (section 4.3). Each system was equilibrated for at least 1 ns (1 ns = 10−10 s), followed by 10 ns of a production run. The production step additionally recorded the pressure tensor elements, which were used in determining the system’s viscosity. Each simulation point for viscosity was the resulting average of four independent viscosity calculations. R-NEMD simulations were performed with 800 molecules placed in a parallelepiped box with Lz = 3 × Lx = 3 × Ly, where Li is the size of the simulation box in the different ith directions of the space (Li are chosen to match the system density extracted from the EMD simulations). An equilibration run was then employed in order to relax the system at constant temperature using the Nose− Hoover NVT integration algorithm. Then, a second equilibration step was applied to generate the velocity gradient in a run of 2 ns, using the R-NEMD scheme. In this case, we have used different exchange frequencies ν to obtain different shear rates. Shear viscosity was then computed as the average of 4 different runs of 2 ns using the methodology previously described. The integration of the equations of motion was performed with a time step of 1 fs.

successive test results obtained by the same operator with the same apparatus under constant operating conditions on identical material) is 2.0 × 10−4 g/cm3, reproducibility (the difference between two single and independent results obtained by different operators working in different laboratories on identical tests material) is equal to 5.0 × 10−4 g/cm3 and the uncertainty of the fuels density is ±3.535 × 10−4 g/cm3. At high pressures, the density was obtained by thermodynamic integration of sound-speed measurements (see below). 2.3. Viscosity Measurements. The measurements at high pressures were made with a vibrating-wire sensor, described by Peleties,28 fitted with a tungsten wire of nominal diameter 0.15 mm. The radius of the vibrating wire was calibrated with noctane at T = 298.15 K and P = 0.1 MPa, and the density and viscosity of the calibration fluid were obtained from NIST Standard Reference Database 69.29 To facilitate the measurements under high pressure, the vibrating-wire sensor and associated magnet assembly were mounted inside a suitable steel pressure vessel immersed in an oil bath for purposes of temperature control. The temperature was measured by means of a calibrated platinum resistance thermometer inserted into a thermowell in the cap of the pressure vessel, and the pressure was measured with a calibrated pressure transducer in the external pipework. The determination of viscosity with the vibrating-wire viscometer requires knowledge of the density ρ of the fluid at each state point. To obtain this information, the ρ(T,P) surface of each fluid studied was obtained by thermodynamic integration from additional experimental measurements of the speed of sound u(T,P) over the same region of T and P, subject to initial conditions (also derived from experimental measurements) at P = 0.1 MPa. The experimental technique and method of analysis have been described in detail elsewhere.30 The experimental uncertainties obtained with this methodology were ±2% for η and ±0.05% for ρ. 2.4. Molecular Dynamics Simulations. The dynamical behavior of the different fuels were analyzed through molecular dynamics simulations at different temperatures and pressures. 2222

dx.doi.org/10.1021/ef300106z | Energy Fuels 2012, 26, 2220−2230

Energy & Fuels

Article

Table 1. General Description of the Different Fluids Analyzed in This Work nature of the sample

sample

σNCa

ρ (g/cm3) @ 288 K

ν (mm2/s)b @ 293/313 K

1

2.73

0.8550

7.780/4.490

diesel

atmospheric distillation (SR)

Middle East SR

2

3.20

0.8730

6.430/3.820

diesel

diesel hydrotreatment

highly naphthenic

3

2.01

0.9520

4.440/2.720

diesel

olefinic and aromatic diesel

4

1.76

0.8428

3.261/2.170

diesel

5

3.10

0.8629

5.090/3.110

diesel

6

2.84

0.8710

5.818/3.488

diesel

fluid catalytic cracking (LCO) hydrotreatment residue (fixed bed) hydroconversion residue (ebullated bed) blend straight-run + LCO

7

3.06

0.8700

7.818/4.264

diesel

atmospheric distillation (SR)

8

1.11

0.7927

3.378/1.448

kerosene

atmospheric distillation (SR)

SR from heavy european crude oil Middle East SR

9

1.40

0.7758

3.899/1.617

kerosene

formulated kerosene

10

2.83

0.8244

3.567/2.522

diesel

atmospheric distillation (SR)

African SR

11

2.95

0.8523

6.974/4.077

diesel

atmospheric distillation (SR)

African SR

12

3.23

0.8864

8.187/4.547

diesel

13

3.87

0.8758

8.249/4.725

diesel

mild hydrocracking of vacuum gas oil (VGO) atmospheric distillation (SR)

blend of 6 SR diesels

14

2.28

0.8495

3.75/2.44

diesel

a

process/origin

VGO hydroconversion (ebullated bed)

characteristics

few olefins

analysis used 2D gas chromatography SFC-twin-2D gas chromatography SFC-twin-2D gas chromatography 2D gas chromatography 2D gas chromatography 2D gas chromatography 2D gas chromatography 2D gas chromatography 2D gas chromatography 2D gas chromatography 2D gas chromatography 2D gas chromatography 2D gas chromatography 2D gas chromatography

σNC is the dispersion in the distribution of carbon number of the sample. bKinematic viscosity measured using NF EN ISO 3103 norm.

250 component clusters, represented by a unique dieseltype molecular compound database (section 4.1). (iii) Viscosity and density measurements were performed on the selected fuel samples at temperatures between 273 and 423 K and pressures up to 400 MPa with specially designed vibrating-wire equipment. The thermodynamic properties of the liquids were obtained by combining the sound-speed data with measurements of density and isobaric specific heat capacity at ambient pressure (section 4.2). (iv) Viscosity was modeled using two different approaches: (a) The first approach reduced fuel samples to typically five pseudocomponents (surrogates) by means of an adapted Lumping method designed to mimic the bulk physical properties of the initial fuel (section 4.3).43 This reduction allowed the use of Molecular Dynamics simulations at equilibrium (EMD) and out of equilibrium (R-NEMD) with the modified Anisotropic United Atom (AUA4m) intermolecular potential to predict thermodynamic properties, equilibrium and shear viscosity (section 4.4). (b) In the second approach molecular-based group contribution methods were used to estimate the viscosity of pure compounds and the detail composition of fluids was used to calculate the mixture viscosity through a new corrected mixing rule (section 4.5). Viscosity at different temperatures was computed with standard ASTM-D1250 and ASTMD341 methods and a new corrected pressure law was used to extend the viscosity from ambient conditions up to the desired pressure.

The Anisotropic United Atom (AUA4m) model for hydrocarbons was used to describe fuel compounds, including n-alkanes,36 iso-alkanes,37 cyclic-alkanes,38 olefins,39 and aromatic compounds.40 According to the AUA approach, certain groups of atoms were treated together simplifying complex molecules and reducing computational time. For hydrocarbons, the effect of hydrogen atoms was not directly accounted for; the carbon atom is considered as the group center on which the group mass was assigned (i.e., CH3, CH2, CH, etc.). In this study, we also used the AUA4m model with a few adjustments to the torsion parameters, notably with regards to improve transport properties41,42 and to better describe aromatic rings, where torsion parameters were adjusted to introduce more rigidity in the molecules. Simulation parameters, a description of the new R-NEDM algorithm, and the AUA4m force-field intermolecular parameters are described in detail in the Supporting Information.

3. WORKING STRATEGY A schematic representation of the approach used in this work can be observed in Figure 2 and is described as follows: (i) A representative set of fuel samples were selected according to its diversity in composition (fuels with different refining processes, composition, oil origins, and distillation cuts and without additives, see Table 1 for details). (ii) A detailed characterization of these fuels was obtained by comprehensive multidimensional gas chromatography and supercritical fluid chromatography (SFC). These techniques enable the characterization of diesel fuels by both chemical family and carbon atom number, leading, in this case, to the molar quantification of approximately

4. RESULTS AND DISCUSSION 4.1. Multidimensional Gas Chromatography Composition Analysis. Seven different diesels (designated 1 to 7) 2223

dx.doi.org/10.1021/ef300106z | Energy Fuels 2012, 26, 2220−2230

Energy & Fuels

Article

Figure 3. Molar distribution by number of carbon atoms and chemical families of hydrocarbons for (a) fluid 1 (2D gas chromatography analysis) and (b) fluid 2 (SFC-twin-2D gas chromatography analysis).

and two aviation-like kerosenes (8 and 9) were selected for the analysis and evaluation of group contribution methods and molecular simulations. Additional diesel samples (10 to 14) were finally used to validate the correlation methods proposed. A brief description of samples 1 to 14, in addition to densities and kinematic viscosities at some thermodynamic conditions, is given in Table 1. The properties, compositions, and gas chromatograms of diesel fluids are shown in the Supporting Information. Diesels 2 and 3, which contain multi-ring naphthenes and olefins, respectively, were analyzed by SFCtwin-2D gas chromatorgraphy, while conventional 2D gas chromatography was applied to the analysis of other samples. A special molecular database have been created to represent all possible molecules identified with the different chromatograms. To illustrate the gain offered by these techniques, detailed characterizations of samples 1 and 2 by chemical family and carbon numbers, with resolution of isomers, are show in Figure 3. Composition analysis for the other diesels are available in the Supporting Information. Significant differences were observed in the compositions of fluids 1 and 2. Fluid 1 was a conventional straight-run (SR) fuel distilled from Middle Eastern crude oil, whereas fluid 2 was hydrotreated after distillation to hydrogenate aromatics, leading to a fuel rich in naphthenic compounds. In general, we can consider both fluids as very close to commercial diesels but without additives. 4.2. Viscosity Measurements. The dynamic viscosities η for fluids 1, 2, and 3 were measured at temperatures of 273 K, 313 K, 373 K, and 423 K in a range of pressures from 0.1 to 400.0 MPa, or until wax formation was observed. The kinematic viscosities for fluids 4−14 were obtained only at ambient pressure and two temperatures. For extrapolation at other temperatures, viscosity was computed using two standards: ASTM-D1250 44 and ASTM-D341.45 According to our calculations and measurements, the mean deviation of this methodology is 5%. The discussion about the behavior of experimental viscosity for fluids 1 and 2 is addressed with the analysis of the molecular simulation results and correlative models (sections 4.4 and 4.5). 4.3. Surrogate Models. There are many different ways to obtain surrogate models for fuel samples,46 using for example a combination of gas chromatography and mass spectroscopy,47 combustion experimental information,48 reproducing distillation curves by discrete multicomponent fuel models.49 In this work, we have used the compositional analysis extracted from the gas chromatography and the physical property database

Figure 4. Surrogate molecules generated by the lumping procedure for (a) model A and (b) model B to represent fluid 1.

previously developed to determine surrogate compounds in order to model fluids 1 and 2. Diesels 1 and 2 were represented with two different sets of surrogates (model A and B), containing five components issued by the lumping procedure proposed by Montel and Goul.43 The lumping method employed in this work to determine the surrogate compounds preserved the macroscopic thermo-physical properties, such as the phase envelope (critical temperature Tc, critical pressure Pc, acentric factor ω, molecular weight Mw, number of carbons NC, carbon-to-hydrogen ratio C/H), of the original fluid.50 This approach was successfully applied to determine surrogate molecules for molecular simulation of physical properties of gasoline. 42 The molecules employed in our molecular simulations are shown in Figures 4 and 5. Analysis of the surrogate molecules generated for fluids 1 and 2 indicated that structurally the surrogates represent the 2224

dx.doi.org/10.1021/ef300106z | Energy Fuels 2012, 26, 2220−2230

Energy & Fuels

Article

Figure 6. Comparison between simulation results (models A and B) and experimental measurements of viscosity for (a) fluid 1 and (b) fluid 2.

Figure 5. Surrogate molecules generated by the lumping procedure for (a) model A and (b) model B to represent fluid 2.

chemical family distribution shown in Figure 3. Though models A and B have different surrogate components, both tended to mirror the chemical families represented in the original fluids: (i) N-paraffins (23.0%), iso-paraffins (25.8%) and one-ringed naphthenes (15.4%) had the highest concentrations in fluid 1 according to gas chromatography data. N-paraffins represented 30.0% of the surrogate components and iso-paraffins represented 20.0% (see Figure 3a). No one-ringed naphthene was present among the surrogate components. It should be noted that fluid 1 had a relatively sharp chemical family distribution; that is, the hydrocarbons present in 1 correspond to a small range of carbon numbers. (ii) Fluid 2 gas chromatography data indicated that its highest concentrations were in naphthenes (a total of 46.2% for one, two- and threeringed naphthenes), indanes and tetralins (15.9%), and alkylbenzenes (9.6%) (see Figure 3b). Its surrogate components were 30.0% naphthenes and contained an equal number of tetralines and alkylbenzenes (20%) (see Figure 5). 4.4. Molecular Simulation. The zero shear rate viscosities of surrogate models A and B were calculated for both fluids 1 and 2 at specific temperatures of 323, 373, and 423 K and at pressures of 0.1, 100.0, 200.0, and 300.0 MPa. The comparison between experimental data and simulation results are shown in Figure 6 for fluids 1 and 2. The mean absolute deviation (MAD) between the experimental and simulated viscosities for fluids 1 and 2 at different thermodynamic conditions is shown in Table 2. The viscosity calculated using model A was very similar to that obtained for model B, and the viscosity predicted by molecular simulation differed by approximately 14−23% from the experimental measurements. MD simulation was also able to predict closely the experimental densities of the samples with

Table 2. Mean Absolute Deviation (MAD) of Viscosity Predicted by Molecular Simulation at All Thermodynamic Conditions from Experimental Viscosities fluid

model

MAD (%)

1

A B A B

23.2 20.6 13.7 16.1

2

a global deviation of less than 5% (see the Supporting Information). Additional viscosity simulations on fluid 1 revealed that the results obtained are not affected by the number of surrogates compounds employed, since we obtained identical results within statistical uncertainties using 10 surrogates as with 5 surrogates (Figure S-19 in the Supporting Information). To evaluate the accuracy of the simulation in predicting the variation of the viscosity with pressure for each fluid, the experimental and simulation data were analyzed with the Barus equation:51,52 ln(η) = ln(η0P ) + αP ·P

(1)

where η0P is the viscosity extrapolated at the limit of zero pressure and αP is the Barus coefficient that accounts for the increment of the viscosity of each diesel with the pressure P. Experimental and simulation data were fitted with eq 1, and the results are shown in Table 3. This table shows that the Barus coefficient was different for each fluid and also is a reverse function of temperature. This relationship with temperature was expected, since an increment in temperature reduces the 2225

dx.doi.org/10.1021/ef300106z | Energy Fuels 2012, 26, 2220−2230

Energy & Fuels

Article

despite their difference in composition, behave very similarly. The comparisons between simulations and experiments suggested that the component database adequately describes the constituents of diesel fuels and that the lumping methodology leads to a reasonable description of the fluid properties even with a small number of surrogates. We computed the shear viscosity as a function of the shear rate for fluid 1 at 313 K and pressures of 0.1 and 200 MPa and at 373 K and 200 MPa in Figure 6. We observed that the shear viscosity obtained at higher shear rates was reduced with respect to the limiting value at zero shear rate (extrapolated from EMD). This phenomenon is known as the shear thinning effect, where viscosity decreases with increasing rate of shear stress. Materials that exhibit shear thinning are called pseudoplastic. Fluid 1 shows a clear shear thinning, and the critical shear rates varied with the thermodynamic condition used in the calculation. Our interest in the rheological behavior of diesels is motivated by the phenomenon of cavitation occurring at the entrance orifice of the diesel injector due to a process of excessive shearing. The drop of viscosity observed in Figure 7

Table 3. Comparison of the Barus Coefficient Obtained from Fitting of eq 1 to Experimental Measurements (exp.) and Simulation Results with Different Surrogate Models αP (×103 MPa−1) fluid 1

2

condition exp. model model exp. model model

323 K

373 K

423 K

11.2 10.4 9.4 12.0 13.2 12.9

8.3 9.1 8.1 9.0 9.1 8.0

7.0 9.0 8.1 7.6 7.3 7.3

A B A B

energy barriers and causes the pressure evolution of the viscosity to be less dependent on the structural details of each fluid (topology of the molecules, rotational dynamics linked with the principal moments of inertia, etc). Hence αP is very similar for each fluid at higher temperatures. The molecular composition details of each fluid are important only at lower temperatures or at high pressures where the slope of viscosity with pressure is different for each sample. We can see that our simulation results described well the variation of viscosity with pressure for each fluid, since the obtained values for αP are in good agreement with experimental values (αP of fluid 2 > αP of fluid 1). The variation of viscosity with temperature for the different fluids at ambient pressure was in better agreement with experiments at higher temperatures. At lower temperatures we observed a small cross over for the computed viscosities between fluids 1 and 2 which was not observed experimentally. To better assess the temperature effect on our viscosity calculations, we determined the activation energy of the viscosity for the different fluids according to 1 ln(η) = ln(ηT 0 ) + ΔEa RT

(2)

where η0T is the viscosity at the limit of infinite temperature, ΔEa is the activation energy of the viscosity, R is the gas constant, and T is the temperature. The use of eq 2 over the data shown in Figure 6 produced the coefficients reported in Table 4.

Figure 7. Variation of the shear viscosity with the shear rate for fluid 1 at 313 K and 0.1 MPa (○) and 200.0 MPa (□) and 373 K at 200.0 MPa (Δ). Filled symbols correspond to model A, and open symbols correspond to model B. Dashed-lines are fitted with the Carreau model. Data at 10−5 ps−1 (1 ps = 10−12 s) are extrapolated from equilibrium MD simulations.

Table 4. Activation Energies (ΔEa) for the Viscosity of the Different Diesels Obtained with the Experimental and Simulation Data at 0.1 MPa diesel 1 2 a

condition exp sima exp sima

η0T(mPa·s) 8.30 1.40 7.90 1.90

× × × ×

−3

10 10−3 10−3 10−3

for fluid 1 is on the order of 70% at high shear rates with respect to the values obtained at equilibrium conditions at 313 K and 0.1 MPa. This difference can reach values up to 92% at 200 MPa. One can observe that both models A and B behaved similarly when the shear rate was increased. Data shown in Figure 7 can be fitted with a three parameters Carreau-like model53

ΔEa (kJ/mol) 15.97 20.5 15.68 19.34

Average results computed using models A and B.

η(γ̇) =

We observed similar activation energies for fluids 1 and 2 as far as the experimental data were concerned. This fact was also observed in our simulations; however, we obtained higher values for both fluids, compared with experiments. It is important to note that the differences in viscosity between fluids 1 and 2 at ambient pressure were on the order of 10% only, which is the same order of precision of our viscosity calculations. At higher temperatures, we observed that the viscosities computed for both fluids became closer. This fact is the consequence of thermal agitation where different fluids,

η0 [1 + (τγ̇)2 ]p

(3)

where η0 is the limiting Newtonian viscosity as the shear rate γ̇ → 0 and τ is a characteristic time constant, which has been found to be roughly equal to the rotational relaxation time for many systems. Simulation results adjusted with eq 3 are shown in Table 5. Our simulations revealed that the inferred rotational relaxation time τ was of the same order of magnitude for models A and B at 313 K and 0.1 MPa (0.18 to 0.24 ns, 2226

dx.doi.org/10.1021/ef300106z | Energy Fuels 2012, 26, 2220−2230

Energy & Fuels

Article

The transition shear rate (TSR ∼ 1/τ, as shown in Table 5), where the fluid becomes pseudoplastic, was shifted to lower values when pressure was increased at constant temperature. When temperature was increased at the same pressure, however, TSR seemed to increase. With this information, we observed that fluid 1 was still Newtonian under the injection conditions mentioned before (T ∼ 400 K and P > 150 MPa with shear rates of 6.5 × 106 s−1), and consequently, no shear thinning induced cavitation should exist for this fluid. It would be interesting to provide, for instance, by means of systematic R-NEMD simulations, a novel correlation of the TSR according to the operational diesel conditions (T, P) and fluid composition. This correlation might be used to screen the properties of new generations of fuels/biofuels in order to establish new operational parameters for engine design. 4.5. Group Contribution (GC) Model. In a similar way to the previous section 4.4, the first step in modeling the viscosity of fuels by GC is to represent all chemical families obtained by the gas chromatography method. We used the database containing 245 compounds to represent fuel constituents identified by the chromatographic analysis. The viscosity of fluids 1−7 were calculated with an original combination of GC methods (Sastri-Rao,21 Van Velzen,22 and Grunberg−Nissan mixing rules23) using the mole fraction distribution for each fluid, and the results are compared with the experimental measurements in Figure 9a. The values reported (empty symbols) show an average absolute deviation of about 20% with respect to experiments. A closer analysis of the results of each fluid revealed a relationship between the carbon−number distribution curve and the deviations from experiment. In fact, the wider the distribution is (presence of more dispersed compounds), the greater becomes the observed error. We used this information to improve the accuracy of the GC method by means of an empirical correction as follows:

Table 5. Regression of the Simulation Results for the Variation of the Shear Viscosity with the Shear Rate for Fluid 1 Using eq 3 model

T (K)

P (MPa)

τ (ns)

p

A B B B

313 313 313 373

0.1 0.1 200.0 200.0

0.24 0.18 2.28 1.28

0.29 0.32 0.24 0.48

respectively). This time should be proportional to the mean average relaxation time of the molecules in the system. This time became strongly dependent on the temperature and pressure, since we obtain values of 2.28 ns compared to 1.28 ns when fluid 1 is at the same pressure but at 313 and 373 K, respectively. The observed values of the exponent p were consistent with those obtained for normal molecular liquids.53 It is interesting to note that experimental devices (such as conic and cylindric viscosimeters) can reach maximum shear rates on the order of 10−7 ps−1 (105 s−1). This value is unreachable by actual nonequilibrium molecular dynamic simulations as a result of a poor signal-to-noise ratio; consequently, the Newtonian plateau can not be accessed in direct simulations.54 This fact prevents any possible validation of our results with standard rheometers. We can mention that multidimensional fluid dynamic modeling flows including the effect of cavitation inside a multipoint diesel fuel injection predicts shear rates of the order of 6.5 × 106 s−1 at ∼400 K and 150 MPa (this value was estimated with an injector of 0.165 mm of diameter with a maximum velocity profile of 540 m·s−1 in the middle of the stream).55 This shear rate is impossible to acquire by actual experimental as well as with our molecular simulation code due the small signal-to-noise ratio in the creation of the velocity profile.56,57 However, a future implementation of TTCF (transient time correlation function) may overcome this limitation.58,59 We have attempted to estimate the magnitude of the Newtonian plateau by means of eq 3, where we imposed the value of η0 as the one obtained in our EMD simulations. Finally, to analyze the Newtonian− pseudoplastic transition under different conditions, the reduced viscosity in function of the shear rate is shown in Figure 8.

c ηGC (xi , Nic , Ti) =

ηGC(xi , Ni , ηi )

[Δη(xi , Nic , Ti) + 1]

(4)

c where ηGC and ηGC are the corrected viscosity (c-GC from now on) and the viscosity obtained by the direct evaluation of the original GC method. Δη is the correction of the mixture viscosity, which is given by

Δη(xi , Nic , Ti) = a(Ti) ·σNc + b(Ti)

(5)

where coefficients a and b are temperature dependent (given in the Supporting Information) and σNC represents the average dispersion of carbon-number obtained for each fluid using gas chromatography: N

σNC =

∑ xi|NiC − ⟨NC⟩| (6)

i=1 N

⟨NC⟩ =

∑ xiNiC (7)

i=1

where ⟨NC⟩ is the average fluid carbon-number, and xi are the carbon number and mole fraction respectively of constituent i and N is the total number of constituents. It is important to notice that evaluation of eq 4 only required knowledge of the carbon-number and mole fraction of the constituents in the fluid, which were obtained by 2D gas NiC

Figure 8. Variation of the reduced shear viscosity (η/η0) in function of the shear rate for fluid 1 at 313 K and 0.1 MPa (○), 200.0 MPa (□) and 373 K and 200.0 MPa (Δ). Filled symbols correspond to model A and open symbols to model B. Lines joining the data are fitting with the Carreau model. 2227

dx.doi.org/10.1021/ef300106z | Energy Fuels 2012, 26, 2220−2230

Energy & Fuels

Article

Figure 10. Comparison of the viscosity predicted with the c-GC method with the experimental data for fluids 1 and 2 in the pressure range from 0.1 to 340.0 MPa and for temperatures of 323 and 373 K. Empty symbols represent the estimation of the viscosity with the c-GC method and the original extrapolation law of Kanti, full symbols represent the viscosity obtained with the c-GC method and the corrected Kanti’s law.

high pressure using an approach similar to the one used in the optimization of the GC method, that is, introducing a correction factor that is a function of temperature and pressure: ηc (P , T ) =

ηK (P)Tref , Pref [1 − Δη(P , T )]

(8)

Figure 9. Comparison of the viscosity predicted by the group contribution methods with the experimental measurements at ambient pressure in the range from 288 to 373 K. (a) Viscosity of fluids 1−9 used in the optimization procedure of the group contribution method. (b) Evaluation of the optimized group contribution method for the calculation of the viscosity for fluids 10−14 between 288 and 423 K. Empty symbols represent the estimation of the viscosity with the original GC method, full symbols represent the viscosity obtained by the c-GC proposed here.

where ηc(P,T) is the corrected viscosity at a given pressure and temperature, ηK(P)Tref,Pref is the viscosity obtained with the evaluation of the Kanti’s law at a given pressure using the viscosity obtained with the c-GC method at temperature Tref and Pref as an input value. Δη(P,T) is the correction term given by

chromatography or SFC-twin-2D gas chromatography measurements. The corrected c-GC method for fluids 1−7 can be observed in Figure 9a (full symbols). In this case, the MAD observed was 9% with respect to the experimental value. The corrected c-GC method proposed here was tested in a predictive way for fluids 10 to 14 in Figure 9b. Here, we predicted the viscosities over a larger range of temperatures (288 to 423 K) and obtained an MAD of 14% and 23% for the c-GC and original GC methods, respectively. The high pressure limit of our approach was explored by using the Kanti’s law24 to extrapolate the viscosity computed for fluids 1 and 2 at ambient pressure obtained by the c-GC method. Comparison between the experimental data and the cGC method results, extrapolated in pressure using the Kanti’s law, in the range from 0.1 to 340.0 MPa for two different temperatures can be observed in Figure 10. The coefficients used in the evaluation of the Kanti’s law have been optimized using experimental viscosity of medium/ heavy crude oils.24 An overall MAD of 22% was observed for these two fluids in the temperature and pressure range evaluated (empty symbols). We can correct the viscosity at

with the pressure and temperature ranges 0.1 MPa < P < 340 MPa and 323 K < T < 373 K. Coefficients C(T) are described in the Supporting Information. The viscosity obtained with this correction can be observed in Figure 10 (full symbols). In general, we obtained an MAD of about 7.3% in the temperature and pressure ranges evaluated. At atmospheric pressure, a direct use of the proposed combination of different GC methods (combined with a mixing rule) and multidimensional gas chromatography analyses on nine fuel mixtures (fluids 1−9) gives an average absolute deviation of 20% close to the one given by the authors21 for pure compounds. In a first step, the use of a corrective term (determined from the 9 fuel mixtures database) improved the predictability of the c-GC method as verified on 5 new fuel mixtures (fluids 10 to 14, MAD = 14%). In a second step, to extend this predictive tool to high pressures, a correction factor was added to the pressure law.24 Finally, for pressures up to 300.0 MPa, the global MAD was reduced to 7%, which was slightly greater than the experimental uncertainty of the best equipment (2% for the vibrating wire viscosimeter used in this work).

Δη(P , T ) = C1(T )P 2 + C2(T )P + C3(T )

2228

(9)

dx.doi.org/10.1021/ef300106z | Energy Fuels 2012, 26, 2220−2230

Energy & Fuels



5. CONCLUSION AND PERSPECTIVES

ACKNOWLEDGMENTS The authors thank Dr. Bernard Rousseau for the use of the Newton Molecular Dynamics code. The authors are indebted with Xavier Montagne, who contributed to focus the project.

The ability of molecular simulation to adequately predict the viscosity and density of diesel samples was studied. Multidimensional gas chromatography provided useful detailed characterization of diesel and kerosene samples, which were reasonably described with the aid of a unique diesel database containing a few hundred of common diesel compounds. Improved group contribution methods, mixing rules, and pressure extrapolation laws were used to model the viscosity with 9−18% mean absolute deviation with respect to experimental data. To the best of our knowledge, this is the first work to propose the use of complete molar composition of fuels as the unique input information to estimate their viscosity for a wide range of thermodynamic conditions. The methodology employed to reduce the large number of diesel constituents to five surrogate components proved to be satisfactory; the surrogate components retained the general distribution and physical properties (viscosity and density) displayed by the original samples when compared with the experimental data at different temperatures and pressures. The simulation reported densities closely mirrored experimentally reported values. Generally 18% and less than 5% deviations, with respect to experimentally determined viscosity and density, respectively, were observed over the range of considered temperatures, 273−423 K, and pressures, 0.1−300 MPa. These preliminary results on the rheological behavior of diesels at high pressure reveals that nonequilibrium MD methodology can be used to estimate whether or not a particular diesel undergoes shear thinning under different shear rates. This approach has the limitation that we are not able to reach the low shear rates values encountered inside normal diesel fuel injectors. The results reported in this work can be used to develop correlative models to aid in the design of high pressure injector engines. In particular, the fact that the viscosity of fuels at 200−300 MPa may increase as much as 1 order of magnitude with respect to viscosities at ambient conditions. This variation should be considered for the development of more accurate injection models, and eventually, they can be adapted or extended to predict the physical properties of diesel/biofuel blends for the development of new generations of fuels.





REFERENCES

(1) Joseph, D.; Belanger, J.; Beavers, G. Int. J. Multiphase Flow 1999, 25, 1263−1303. (2) Farrell, J. T.; Cernansky, N. P.; Dryer, F. L.; Friend, D. G.; Hergart, C. A.; Law, C. K.; McDavid, R. M.; Mueller, C. J.; Patel, A. K.; Pitsch, H. SAE J. 2007, 01, 0201−0227. (3) Hendl, S.; Millat, J.; Vogel, E.; Vesovic, V.; Wakeham, W. A.; Luettmer-Strathmann, J.; Sengers, J. V.; Assael, M. J. Int. J. Thermophys. 1994, 15, 1−31. (4) Reid, R. C.; Prausnitz, J. M.; Poling, B. E. The Properties of Gases and Liquids, 4th ed.; McGraw-Hill: New York, 1987. (5) Vogel, H. Phys. Z. 1921, 22, 645−646. (6) Mehrotra, A. K. Can. J. Chem. Eng. 1994, 72, 554−557. (7) Sandler, S. I.; Orbey, H. The viscosity of liquid hydrocarbons and their mixtures. In Advances in Engineering Fluid Mechanics: Multiphase Reactor and Polymerization System Hydrocarbons, 1st ed.; Gulf Professional Publishing: Oxford, U.K., 1996. (8) Orbey, H.; Sander, S. I. Can. J. Chem. Eng. 1993, 71, 437−446. (9) Aasberg-Petersen, K.; Knudsen, K.; Stenby, E. H. F. Phase Equilib. 1991, 70, 293−308. (10) Queimada, A. J.; Rolo, L. I.; Caço, A. I.; Marrucho, I. M.; Stenby, E. H.; Coutinho, J. A. P. Fuel 2006, 85, 874−877. (11) Srinivasan, C.; Krishna M., M. V. Int. J. Refrig. 1985, 8, 13−16. (12) Lohrenz, J.; Bray, B. G.; Clark, C. R. J. Pet. Tech. 1964, 1171− 1176. (13) Et-Tahir A. Détermination des variations de la viscosité de divers hydrocarbures en fonction de la pression et de la température. Etude critique de modèles représentatifs. PhD thesis; Université de Pau et des Pays de l’Adour, Pau, France, 1993. (14) Dymond, J. H.; Assael, M. J. Modified Hard-Spheres Schemes, in Transport Properties of Fluids, Their Correlation, Prediction, and Estimation, 1st ed.; Millat, J., Dymond, J. H., Nieto de Castro, C. A., Eds.; Cambridge University Press: Cambridge, U.K., 1996. (15) Ö zdogan, S.; Gürbüz Yücel, H. Fuel 2000, 79, 1209−1214. (16) Hildebrand, R. H.; Lamoreaux, J. H. J. Phys. Chem. 1973, 77, 1471−1473. (17) Przezdziecki, J. W.; Sridhar, T. AIChE J. 1985, 31, 333−335. (18) Marcus, Y. Fluid Phase Equilib. 1999, 154, 311−321. (19) Latini, G.; Cocci Grifoni, R.; Passerini, G. Transport Properties of Organic Fluids, 1st ed.; WITPress: Cambridge, U.K., 2006. (20) Marano, J. J.; Holder, G. D. Ind. Eng. Chem. Res. 1997, 36, 2409−2420. (21) Sastri, S. R. S.; Rao, K. K. Chem. Eng. J. 1992, 50, 9−25. (22) Van Velzen, D.; Cardozo, R. L.; Langenkamp, H. Ind. Eng. Chem. Fundam. 1972, 11, 20−25. (23) Grunberg, L.; Nissan, A. H. Nature 1949, 164, 799−800. (24) Kanti, M. Viscosité de mélanges d’alcanes et d’alkylbenzènes en fonction de la pression et de la température: application aux coupes pétrolières. PhD thesis,Université de Pau et des Pays de l’Adour, Pau, France, 1988. (25) Vendeuvre, C.; Ruiz-Guerrero, M.; Bertoncini, F.; Duval, L.; Thiébaut, D.; Hennion, M. C. J. Chromatogr., A 2005, 1086, 21−28. (26) Bruno, T. J.; Ott, L. S.; Lovestead, T. M.; Huber, M. L. J. Chromatogr., A 2010, 1217, 2703−2715. (27) Adam, F.; Thiebaut, D.; Bertoncini, F.; Courtiade, M.; Hennion, M. C. J. Chromatogr., A 2010, 1217, 1386−1394. (28) Peleties, F. Advanced fluid property measurement for oilfield applications. PhD Thesis, University of London, London, U.K., 2007. (29) Lemmon, E. W.; McLinden, M. O.; Friend, D. G. NIST Standard Reference Database Number 69: Thermophysical Properties of Fluid Systems; National Institute of Standards and Technology: Giathersburg, MD, 2011. Available online: http://webbook.nist.gov/ chemistry/.

ASSOCIATED CONTENT

S Supporting Information *

Complete description of the experimental devices, methods, parameters, and the molecular fluid database. This material is available free of charge via the Internet at http://pubs.acs.org.



Article

AUTHOR INFORMATION

Corresponding Author

*E-mail: [email protected]. Present Addresses ∥

Atlantic LNG Company, Trinidad and Tobago Materials Design, Materials Design s.a.r.l., 18 rue de Saisset, Montrouge, France ⊥

Notes

The authors declare no competing financial interest. 2229

dx.doi.org/10.1021/ef300106z | Energy Fuels 2012, 26, 2220−2230

Energy & Fuels

Article

(30) Peleties, F.; Segovia, J. J.; Trusler, J. P. M.; Vega-Maza, D. J. Chem. Thermodyn. 2010, 42, 631−639. (31) Van-Oanh, N.-T.; Houriez, C.; Rousseau, B. Phys. Chem. Chem. Phys. 2010, 12, 930−936. (32) Müller-Plathe, F. J. Chem. Phys. 1997, 106, 6082−6085. (33) Bordat, P.; Müller-Plathe, F. J. Chem. Phys. 2002, 116, 3362− 3369. (34) Berendsen, H. J. C.; Postma, J. P. M; van Gunsteren, W. F.; Di Nola, A.; Haak, J. R. J. Chem. Phys. 1984, 81, 3684−3690. (35) Andersen, H. C. J. Comput. Phys. 1983, 52, 24−34. (36) Ungerer, P.; Beauvais, C.; Delhommellle, J.; Boutin, A.; Rousseau, B.; Fuchs, A. J. Chem. Phys. 2000, 112, 5499−5510. (37) Bourasseau, E.; Ungerer, P.; Boutin, A.; Fuchs, A. Mol. Simul. 2002, 28, 317−336. (38) Bourasseau, E.; Ungerer, P.; Boutin, A. J. Phys. Chem. B 2002, 106, 5483−5491. (39) Bourasseau, E.; Haboudou, M.; Boutin, A.; Fuchs; Ungerer, P. J. Chem. Phys. 2003, 118, 3020−3034. (40) Ahunbay, M. G.; Perez-Pellitero, J.; Contreras-Camacho, R. O.; Teuler, J. M.; Ungerer, P.; Mackie, A. D.; Lachet, V. J. Phys. Chem. B 2005, 109, 2970−2976. (41) Nieto-Draghi, C.; Ungerer, P.; Rousseau, B. J. Chem. Phys. 2006, 125, 044517. (42) Nieto-Draghi, C.; Bocahut, A.; Creton, B.; Have, P.; Ghoufi, A.; Wender, A.; Boutin, A.; Rousseau, B.; Normand, L. Mol. Sim. 2008, 34, 211−230. (43) Montel, F.; Gouel, F. A new lumping scheme of analytical data for compositional studies. SPE Annual Technical Conference and Exhibition, Houston, TX, Sept. 16−19, 1984; SPE 1984-13119-MS. (44) ASTM D1250, Standard Guide for Petroleum Measurement Tables, Table 53B, Generalized Products, Correction of Observed Density to Density at 15 °C; ASTM: West Conshohocken, PA, 2008. (45) ASTM D341, Standard Practice for Viscosity Temperature Charts for Liquid Petroleum Products; ASTM: West Conshohocken, PA, 2009. (46) Pitz, W. J.; Mueller, C. J. Prog. Energy Combust. Sci. 2011, 37, 330−350. (47) Huber, M. L.; Lemmon, E. W.; Diky, V.; Smith, B. L.; Bruno, T. J. Energy Fuels 2008, 22, 3249−3257. (48) Herbinet, O.; Pitz, W. J.; Westbrook, C. K. Combust. Flame 2010, 157, 893−908. (49) Anand, K.; Ra, Y.; Reitz, R. D.; Bunting, B. Energy Fuels 2011, 25, 1474−1484. (50) Lugo, R.; Ebrahimian, V.; Lefèbvre, C.; Habchi, C.; de Hemptinne, J. C. A compositional representative fuel model for biofuelsApplication to diesel engine modelling. SAE 2010 Powertrains Fuels and Lubricants Meeting, San Diego, CA, Oct. 2010; SAE 2010-01-2183. (51) Hamrock, B. J. Fundamentals of Fluid Film Lubrication; McGrawHill: New York, 1994. (52) Stachowiak, G. W.; Bachelor, A. W. Engineering Tribology; Elsevier: Amsterdam, The Netherlands, 1993. (53) Kelbar, M.; Rafferty, J.; Maginn, E.; Siepmann, I. J. Fluid Phase Equilib. 2007, 260, 218−231. (54) Evans, D. J.; Morris, G. P. Statistical Mechanics of Nonequilibrium Liquids; Academic Press: New York, 1990. (55) Habchi, C.; Dumont, N.; Simonin, O. Atomization Sprays 2008, 18, 129−162. (56) Müller-Plathe, F. Phys. Rev. E 1999, 59, 4894−4898. (57) Bordat, P.; Müller-Plathe, F. J. Chem. Phys. 2002, 116, 3362− 3369. (58) Pan, G.; McCabe, C. J. Chem. Phys. 2006, 125, 194527. (59) Mayzar, O.; Pan, G.; McCabe, C. Mol. Phys. 2009, 107, 1423− 1429.

2230

dx.doi.org/10.1021/ef300106z | Energy Fuels 2012, 26, 2220−2230