Ignition in an Atomistic Model of Hydrogen Oxidation - ACS Publications

Feb 7, 2017 - technologies. Microscopic details of the combustion processes, ephemeral transient species, and complex reaction networks are necessary ...
1 downloads 14 Views 1MB Size
Subscriber access provided by YORK UNIV

Article

Ignition in an Atomistic Model of Hydrogen Oxidation Mohammad Alaghemandi, Lucas B Newcomb, and Jason R Green J. Phys. Chem. A, Just Accepted Manuscript • DOI: 10.1021/acs.jpca.7b00249 • Publication Date (Web): 07 Feb 2017 Downloaded from http://pubs.acs.org on February 8, 2017

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 A 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 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

The Journal of Physical Chemistry

Ignition in an Atomistic Model of Hydrogen Oxidation Mohammad Alaghemandi,† Lucas B. Newcomb,† and Jason R. Green∗,†,‡,§ † ‡ §

Department of Chemistry, University of Massachusetts Boston, Boston, MA 02125 Department of Physics, University of Massachusetts Boston, Boston, MA 02125 Center for Quantum and Nonequilibrium Systems, University of Massachusetts Boston, Boston, MA 02125

ABSTRACT: Hydrogen is a potential substitute for fossil fuels that would reduce the combustive emission of carbon dioxide. However, the low ignition energy needed to initiate oxidation imposes constraints on the efficiency and safety of hydrogen-based technologies. Microscopic details of the combustion processes, ephemeral transient species, and complex reaction networks are necessary to control and optimize the use of hydrogen as a commercial fuel. Here, we report estimates of the ignition time of hydrogen-oxygen mixtures over a wide range of equivalence ratios from extensive reactive molecular dynamics simulations. These data show that the shortest ignition time corresponds to a fuel-lean mixture with an equivalence ratio of 0.5, where the number of hydrogen and oxygen molecules in the initial mixture are identical, in good agreement with a recent chemical kinetic model. We find two signatures in the simulation data precede ignition at pressures above 200 MPa. First, there is a peak in hydrogen peroxide that signals ignition is imminent in about 100 ps. Second, we find a strong anti-correlation between the ignition time and the rate of energy dissipation, suggesting the role of thermal feedback in stimulating ignition.

∎ INTRODUCTION

pounding factors point to the need for a better understanding of the microscopic origins of ignition.

Using hydrogen for internal combustion engines or in fuelcell based electric vehicles holds promise as an effective way to mitigate pollution in the atmosphere. 1–3 However, one of the primary concerns of using hydrogen as a fuel is the potential for uncontrolled fire and explosion. 4 The root causes of this instability include the low ignition energy of hydrogen relative to hydrocarbons, the high flame temperature, and the wide flammability range. 5 Leaks from highly pressurized containers, for example, can ignite because the Joule-Thompson coefficient of hydrogen is negative above room temperature. 6 Fundamental predictions of the autoignition conditions of hydrogen from the detailed molecular dynamics are important to the development of reliable means of controlling ignition, as well as the broader use of hydrogen as a clean fuel. 1,7 Practical applications of combustion often rely on explosion limits, which are the temperature and pressure boundaries that divide explosive behavior from comparatively slower burning. 8 For hydrogen oxidation, at sufficiently high pressures and temperatures, previous experimental and theoretical evidence implicates HO2 and H2 O2 as critical precursors to ignition. The elementary reactions create and annihilate these species are part of a chain-branching chemical-feedback mechanism that can lead to radical pool growth and, ultimately, ignition. 9,10 A thermal feedback mechanism can also contribute to ignition when the heat released by the exothermic, elementary reactions does not efficiently dissipate. Phenomenological models of the hydrogen combustion mechanism also indicate that above ≈ 1150 K and ≈ 3.5 MPa, spontaneous ignition becomes pressure independent. 9 These comE-mail: [email protected]

A hydrogen-oxygen mixture spontaneously (auto)ignites when initially slow thermal reactions suddenly give way to chain-branching reactions that proliferate radicals, sustaining and accelerating the overall oxidation reaction. 11 Massivelyparameterized chemical kinetic models form the foundation of our current understanding of chain-branching reactions, but, even with extensive experimental validation, these models may still be inaccurate, incomplete, or only work for a specific range of temperatures and pressures. 12–28 Reactive molecular dynamics (MD) simulations, however, provide a comparatively unbiased description of complex chemical reactions in that they do not require any prior knowledge of the possible mechanisms or intermediates. 29–33 Parameterized models are also necessary for a classical treatment of atomistic systems and simulation of the molecular dynamics. However, one can readily validate the interatomic potentials of these models against ab initio electronic structure calculations. 30,31,33,34 At the molecular level, the ignition of a gaseous mixtures of hydrogen and oxygen quickly and irreversibly produce water and large amounts of energy. Molecular dynamics simulations have been used to study hydrogen combustion. 31,35 Agrawalla et al. 31 used optimized ReaxFF potential functions to run MD simulations of hydrogen combustion at fixed number of atoms, volume, and temperature, NV T . Their analysis of the reactions produced by the dynamics suggested hydroperoxyl radical, HO2 , plays a key part in the reaction kinetics at T > 3000 K and P > 400 atm. More recently, 35 two of us investigated hydrogen oxidation at NV T , over wide range of temperature and density at NV T conditions, and analyzed the chemistry with tools from symbolic dynamics. This formalism avoids ele1

ACS Paragon Plus Environment

The Journal of Physical Chemistry

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 ACS Paragon Plus Environment

Page 2 of 9

Page 3 of 9

100

∎ MODEL AND METHODS

Counts

(c)

An objective here is to collect statistics for the ignition of hydrogen-oxygen mixtures directly from atomistic dynamics for comparison to chemical kinetic models. We accomplish this through extensive molecular simulations over a range of mixture compositions. To quantitatively characterize the initial composition, we use the equivalence ratio, φ . The equivalence ratio is the fraction of fuel-oxidizer ratio in a mixture of interest to the fuel-oxidizer ratio in the stoichiometric mixture. With this definition, a mixture with φ = 1 is stoichiometric, with φ < 1 is fuel lean, and with φ > 1 is fuel rich. We consider a range of φ from 0.08 to 3.31. For the stoichiometric mixture, 66 H2 and 33 O2 molecules occupy a cubic simulation box with volume V = 8 nm3 and periodic boundary conditions. For other equivalence ratios, a simulation box with identical dimensions was initially filled with the corresponding number of hydrogen and oxygen molecules such that there are always 99 molecules. At high pressures, the chemical kinetic model from Li et al. 37 suggests the dominant initiation reaction in hydrogen oxidation is H2 + O2 → HO2 + H, which is followed by the fast reaction of H + O2 → O + OH. 10 In our simulations, we seed the chemistry with a single OH radical, which neglects the effect of the O radical, but reduces the simulation time (cost) to initiate the reaction through this sequence of reactions. All molecular dynamics (MD) simulations were performed with the PuReMD-GPU simulation package 38 using the ReaxFF potential. 30–33,39 The MD simulations were carried out at constant number of atoms, volume, and energy (NV E) with a time step of 0.1 fs and total simulation time of 3 ns. This time was sufficient to convert more than 80% of the reactants to water. Under these conditions, both temperature and pressure vary during the course of the reaction (Figure 1 and Figure S1). For each φ , in total 11 unique mixtures, we ran 100 independent simulations at an initial temperature, Ti , of 1500 K and a density, ρi , of 250 kg m−3 . Note that subscripts indicate initial conditions. A separate set of 1000 parallel simulations of stoichiometric mixtures with Ti =1000 K and ρi =250 kg m−3 was also simulated.

∎ RESULTS AND DISCUSSION Kinetic Temperature and the Onset of Ignition. Hydrogen oxidation is an exothermic reaction producing large amounts of heat and water. If the heat generated does not or cannot dissipate to an external bath, thermal energy will accumulate, increase the temperature of the mixture, and, consequently, the rate of the overall reaction. This thermal feedback mechanism is apparent in our MD simulations. Figure 1 (a) shows the time profile for the total kinetic energy, Ekin , from 8 out of 1000 trajectories, selected to show the range of

0

200

3.0 Ignition time, tig (ns)

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

400

0 800

600

(a)

3.0

2.5

2.5

2.0

2.0

1.5

1.5

1.0

1.0

0.5

0.5 0

Counts 0

100

(b)

200 400 600 800 Dissipation rate, α (kcal/mol ns)

Figure 2 (a) The dissipation rate, the rate at which Ekin increases prior to ignition, anti-correlates with the ignition time, tig . The kinetic energy increases more rapidly for trajectories with shorter ignition time. Longer ignition times correspond to slower rates of energy dissipation. (b) Empirical distribution of the ignition times over a sample of 1000 trajectories. (c) Empirical distribution of the dissipation rate, α .

temperature profiles. The kinetic energy (or the kinetic temperature) of the isolated system gradually increases before a sudden jump, which signals an ignition event. The ignition delay time is of paramount importance in practical applications of hydrogen combustion. The data in Figure 1 (a) are an opportunity for comparison to macroscopic estimates of the ignition time scale. An initial challenge to this comparison is that the ignition time varies across the trajectory ensemble sample: the discrete jump in kinetic energy occurs at a different time for each simulation. The temporal location depends sensitively on the initial phase point. Numerically, it also depends on the integration algorithm and the hardware architecture. A second challenge in comparing the predictions from microscopic dynamics to bulk models is the definition of the zero of time. In this work, we set the origin of time, t = 0 at the simulation time of the first reaction in each simulation. We can then define the ignition time with a criterion analogous to the “temperature inflection criterion”. 18 The ignition time, tig , is the time associated with the kinetic energy midpoint indicated by a horizontal dashed line in Figure 1 (a); max min )/2. The maximum that is, when Ekin is equal to (Ekin + Ekin and minimum kinetic energy are essentially fixed quantities for all trajectories. As the overall reaction progresses, the total kinetic energy of each mixture increases steadily prior to ignition but with a slope that depends on the particular deterministic trajectory. An interesting question is whether the dynamical observable, 3

ACS Paragon Plus Environment

The Journal of Physical Chemistry

1.00

Ekin , has any prescient signatures of ignition. The one thousand trajectory set (for an initially stoichiometric mixture) is sufficient to begin addressing this question. A chemical signature of an upcoming ignition is the average number of hydrogen peroxide molecules over time, shown in Figure 1 (b) for a select set of trajectories. Each time series is an average over 11 trajectories with a similar Ekin profile for each of the trajectories in Figure 1 (a). The dispersion of tig in the set of trajectories is less than 10 ps. The mole fraction of hydrogen peroxide gradually increases before a sudden drop at tig . An apparent thermal indicator of mixtures with a pending ignition event is the initial rate of increase in kinetic temperature. Figure 2 shows ignition time estimates against the slope, α, of the kinetic energy versus time for 1000 trajectories from an initially stoichiometric composition. The dissipation rate, as measured by the slope α, is from a linear fit of kinetic energy data points up to 0.1 ns prior to ignition. The statistics exhibit a crescent-like region and a clear nonlinear correlation between tig and α. While shallower growth of the kinetic energy prior to ignition corresponds to a longer delay before ignition, a rapid dissipation of kinetic energy translates into short ignition times. We note that the energy is dissipated by individual reactions but remains trapped within the simulation cell. There is trajectory-to-trajectory variability in the sample data set, but there are coarse trends. For example, in the upper region, the ignition time decreases gradually with increasing α, while the lower region the decay of the ignition time tig is nearly exponential with α. The distribution of ignition times and dissipation rates, α, are shown in Figure 2 (b) and (c). The mean dissipation rate is 321.131 ± 96.834 kcal/mol ns and the mean ignition time is t¯ig = 1.778 ± 0.372 ns. Both distributions appear to be asymmetric with a tail at larger values of tig and α, though the tail is more pronounced for the distribution of dissipation rates.

Figure 3 Time series for the average number of (a) reactants and product and (b) intermediate species. In these data, the time, τ = t −tig , is the simulation time, t, shifted by the corresponding ignition time, tig > 0, of each trajectory estimated by the temperature-inflection criterion. The resulting zero of time is at tig , with negative (positive) values of τ prior to (following) ignition. Results are averaged over one thousand simulations with Ti =1000 K, ρi =250 kg m−3 , and φ = 1.

Species Evolution and the Ignition Time. The simulated combusting systems of hydrogen and oxygen evolve away from, but relax towards, equilibrium. Over the course of the reaction a population of intermediate species (including OH, H2 O2 , HO2 , H, and O) is generated and consumed upon formation of the final product, water. Figure 3 (a) shows the time series for the mole fraction, χ, of reactants and products. All data in this figure is averaged over 1000 trajectories and the time axis is shifted so that the origin, τ = 0, denotes the ignition time of each trajectory. Times τ < 0 are prior to ignition and times τ > 0 are after ignition. During the course of the reaction, the mole fractions of H2 and O2 steadily reduce before a sharp crossover at the ignition time when χH2 O increases sharply. A few picoseconds after ignition, the mole fractions are stationary. The final mole fractions χH2 and χO2

are nonzero, and hence χH2 O ≠ 1. This result is a consequence of the total energy being fixed, which means the thermal energy dissipated by the elementary reactions is not removed, leaving the mixture at an stationary state with a high temperature, typically more than 4000 K (Figure 1). Time profiles for the mole fractions χ of the transient intermediate species are shown in Figure 3 (b). The mole fraction of H is near zero before ignition, but increase sharply at τ = 0. The magnitude of its crossover, however, is small in comparison to the jump in χOH . After ignition, the mole fraction of these species reach nearly constant values, χOH ≈ 0.06 and χH ≈ 0.01. The mole fraction of HO2 , with a value of about 0.01 at τ < 0, shows a sharp peak at τ = 0. χO is close to zero at τ < 0 and increases slightly after ignition, as expected

(a)

Mole fraction, χ

0.80

H2 O H2

0.60

0.40 O2 0.20

0.00 -0.4

-0.2

0.0

0.2

0.4

0.10

(b) 0.08 Mole fraction, χ

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 9

H2 O2

OH

0.06

0.04

0.02

H

HO2 0.00

O -0.4

4

ACS Paragon Plus Environment

-0.2

0.0 0.2 Time, τ (ns)

0.4

Page 5 of 9

(a)

(b) 3.0

3.0 2.5

2.5

2.0

Ignition time, t¯ig (ns)

2.5 Ignition time, t¯ig (ns)

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

1.5

2.0

1.0

0.1

1.0

10.0

1.5

1.5

1.0

1.0

0.5 0.0

2.0

0.5

1.0 1.5 2.0 2.5 Equivalent ratio, φ

3.0

3.5

0.5 0.0

0.2

0.4 0.6 0.8 Hydrogen fraction, χH2

1.0

Figure 4 Average ignition time versus (a) equivalence ratio, φ , and (b) initial hydrogen mole fraction, χH2 . Each data point is a sample average over one hundred simulations of the ignition time using the temperature inflection criterion. The initial temperature and density were Ti =1500 K and ρi =250 kg m−3 . Inset is a linear-log plot of the same data. The solid blue line in (b) is the quadratic function in the main text. The green line is the inverse Shannon entropy in the main text. The gray line is the prediction of the ignition time from the work of S´anchez et al. 18 at a constant temperature of 3090 K and a constant pressure of 600 bar; the temperature is close to the inflection point of the stoichiometric mixture and the pressure was adjusted to fit the simulation data.

for a short-lived intermediate that only participates in one of the branching reactions. The kinetics of H2 O2 is particularly interesting in the context of the simulated ignition process. At τ < 0, χH2 O2 increases gradually up to τ ≈ −100 ps with χH2 O2 ≈ 0.08, then it suddenly decreases. The sudden drop in χH2 O2 continues until a few picoseconds after ignition, where it slightly increases and remains constant with a value of 0.002. The peak in the mole fraction of H2 O2 is an apparent signature that precedes ignition by 100 ps. However, as Figure 3 (b) indicates, prior to ignition, χH2 O2 is always the most probable intermediate species before the system ignites. Thus, there might be the possibility of monitoring the H2 O2 concentration and inferring its peak to predict the onset of ignition. The results of our atomistic simulations are in a good agreement with the species profiles predicted by a kinetic model using CHEMKIN 40 and GRI-Mech 3.0 41 (for a direct comparison see Figures S2, S3, and S5 in the SI). Comparing the profiles of H2 , O2 , OH, HO2 , and H2 O for atomistic simulations and numerical solutions of mass-action rate equations under adiabatic conditions, we see strikingly similar results. There is one discrepancy worth further discussion. In the average profiles from molecular dynamics simulations, there is a high concentration of H2 O2 (relative to other intermediates) preceding ignition. There is also a distinct H2 O2 peak in the kinetic model concentration profiles, but the peak is a smaller compared to that from molecular dynamics, the concentration even stays below that of HO2 , and occurs at the ignition time.

One other minor discrepancy is the larger concentration of H and O radicals suggested by the GRI-Mech 3.0 than molecular dynamics. It remains to be tested whether these discrepancies are a consequence of the finite-size of the simulations or the quality of the kinetic model at these temperatures and pressures. It is now possible to measure the concentration of shortlived chemical species at high flame temperatures using ultrafast spectroscopy. 42,43 Such experiments, if they can be done at high pressures, might be capable of monitoring the concentration of H2 O2 and verifying this characteristic temporal profile. Although the mole fraction of HO2 is much smaller than that of H2 O2 (Figure 3 (b)), the elementary reactions that produce and consume HO2 have the highest probability among the observed reactions in our simulations (SI Table S1). The most probable elementary reactions are HO2 → O2 + H and H2 + O2 → HO2 + H. Therefore, manipulating χHO2 might also be another control parameter for the ignition delay time. Effect of Fuel-Oxidizer Ratio on the Ignition Time. From the macroscopic chemistry perspective, the temperature, pressure, and stoichiometric ratio of reactants determines whether the reaction will complete. Having more oxygen in the mixture, although it ensures the complete burning of the fuel, might change the ignition time and hence have practical implications (e.g., the efficiency of a combustion engine). It is known from experiments that fuel-oxidizer ratios influence en5

ACS Paragon Plus Environment

The Journal of Physical Chemistry

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

gine spark timings and also determine flame speeds. 44 Driven by these practical issues, we ran simulations over a range of equivalence ratios. Figure 4 (a) shows the average ignition time, t¯ig , as a function of the equivalence ratio. For equivalence ratios less than 1, the mixtures have a higher load of oxygen relative to a stoichiometric mixture. Starting from a stoichiometric mixture and adding more oxygen reduces t¯ig until the minimum at φ = 0.5. Further increasing the amount of oxygen also leads to an increase in t¯ig . Mixtures with lower levels of oxygen, φ > 1, also produce larger ignition times. The shape of t¯ig in Figure 4 (a) is partly a result of using the equivalence ratio φ as the dependent variable. Figure 4 (b) shows the same data instead using the hydrogen mole fraction in the initial mixture, χH2 . This representation of the data shows that a mixture with χH2 = 0.5 has the shortest ignition time. Any deviation from this composition increases t¯ig . Under the chosen simulation conditions, any change in the oxygen load of the mixture away from φ = 0.5 increases the ignition time. We considered fits of these data to several common functions. For example, a parabolic function of the form 2 2π (χH2 − 21 ) +1.0 matches well with the simulation data near the minimum (blue line in Figure 4 (b)); however, at small or large values of χH2 , i.e., χH2 → 0 or χH2 → 1, this function has a finite values for t¯ig of about 2.6 ns (dashed blue lines). Another function, the inverse of the Shannon entropy, −1 t¯ig = c(−χH2 log(χH2 ) − (1 − χH2 )log(1 − χH2 )) ,

(1)

also agrees with these data. The constant c ≈ 2/e and e is Euler’s number and c is a T and P (or a E/N, and N/V ) dependent quantity. Although in the range of our data, Eq. (1) does not fit as well as the parabolic function (green line in Figure 4 (b)), at the limits of χH2 it does exhibit some features we expect: t¯ig → ∞ when χH2 → 0 or 1, dashed green line in Figure 4 (b). These results suggest the qualitative features of the ignition time are a consequence of the number of collisions between hydrogen and oxygen in the initial mixture. In a combusting mixture with initially identical proportions of H2 and O2 , χH2 = 0.5, there is a high probability for hydrogen and oxygen collisions to facilitate the production of the transient intermediate species. However, if the number of oxygen (hydrogen) molecules is disproportionately large, such that it dilutes the mixture, this lowers the probability of reactive collisions among the intermediate species. A recent study examines the collisions and initial reactions more closely. 34 Since the intermediates are essential for the chain-branching mechanism, a longer ignition time is observed. In Figure 4, we also show the result of a recent chemical kinetic estimate of the ignition time. S´anchez et al. have derived a formula for the ignition delay time, tig , of hydrogen-air

mixtures for a reduced mechanism of five reactions. 18 The tig is shown (solid grey line) as a function of the mixture composition at a fixed pressure of 600 bar and a fixed temperature T = 3090 K. This temperature is near the inflection in MD results (see Figure S1). The pressure, chosen to fit the MD data, is lower than the simulation pressures we estimate through the ideal gas equation of state. Although our simulations are at a constant volume and energy, not fixed T and P, the analytical form is within the statistical error of our ignition time estimates. The ignition delays predicted by CHEMKIN 40 and GRI-Mech 3.0 41 exhibit similar trends at constant temperature and pressure conditions (T = 3090 K and P = 600 bar), but also at constant volume and energy conditions (with Ti = 3000 K and Pi = 1000 bar), although the ignition times are further offset from our MD results (Figure S7). In part, the deviations between the chemical kinetic model predictions and simulation data are the result of their differing assumptions and conditions. The conditions of an adiabatic, homogeneous reactor seem to most closely match our simulation conditions. Of course, while the system boundary is adiabatic in the NV E simulations, we impose no homogeneity and molecules are free to diffuse within the available space. An important difference is that MD simulation box consists of only 200 atoms while the rate coefficients and kinetic models only strictly hold in an infinite-system size limit. The finite size of our simulations, however, also have the beneficial effect of faster mixing. Another difference is the elementary reactions and rate coefficients in kinetic models are validated at much lower temperatures and pressures than those in our MD simulations (T varies from 1000 − 5000 K and pressures P > 1000 bar). Despite these factors, these data are preliminary evidence that the ignition delay time, as defined from a time zero of the first reaction, is a reasonable definition. Though computationally demanding, these data open up the possibility of using ignition times from MD simulations to compare against more complicated fuels and as a potential independent reference point for kinetic modelling. In sum, the differences await further constant energy simulations to assess qualitative and quantitative features caused by finite-size effects.

∎ CONCLUSIONS Using reactive molecular dynamics simulations, we studied ignition in an atomistic model of hydrogen combustion with a variety of fuel-oxidizer ratios over temperatures and pressures that span the second extended explosion limit. These simulations show that mixtures with a higher dissipation rate initially will ignite on shorter time scales. Our findings also indicate that hydrogen-oxygen mixtures with an equivalence ratio of 0.5 (fuel-lean conditions) have the shortest ignition delay; increasing or decreasing the oxygen load uniformly enhances the time for ignition. By tracking the mole fraction

6

ACS Paragon Plus Environment

Page 6 of 9

Page 7 of 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

The Journal of Physical Chemistry

of transient intermediate species over time, features of these time series are possible indicators of the onset of ignition. Among five metastable intermediate species in the hydrogen combustion, H2 O2 and HO2 profiles are promising signatures for forecasting ignition at elevated temperatures and pressures. Looking ahead, manipulating the magnitude of this component could be a means of controlling the ignition delay in high-temperature and pressure hydrogen combustion.

∎ ASSOCIATED CONTENT Supporting Information I Figure S1 (a), temperature and pressure profiles during the hydrogen oxidation of a stoichiometric mixture; Figure S1 (b), instantaneous temperature T and pressure P estimates from a single simulation compared to the second explosion limit associated with the S´anchez et al. kinetic model; Figure S2, profiles of the chemical species generated by CHEMKIN and GRI-Mech 3.0 using the T and P from Figure S1 (a); Figure S3, profile of the chemical species generated by CHEMKIN and GRI-Mech 3.0 with Ti = 2000 K and Pi = 2000; Figure S4, temperature and pressure profiles generated by CHEMKIN and GRI-Mech 3.0 with Ti = 2000 K and Pi = 2000; Figure S5, profiles of chemical species generated by CHEMKIN and GRI-Mech 3.0 with Ti = 2000 K and Pi = 6000; Figure S6, temperature and pressure profiles generated by CHEMKIN and GRI-Mech 3.0 with Ti = 2000 K and Pi = 6000; Figure S7, estimated ignition time using two kinetic models and varying macroscopic conditions; Table S1, elementary reactions observed in the atomistic simulation for a stoichiometric mixture.

Supporting Information II A movie of a molecular dynamics simulation that shows a mixture undergoing an ignition event as a result of hydrogen combustion.

∎ ACKNOWLEDGMENTS This material is based upon work supported by the U.S. Army Research Laboratory and the U.S. Army Research Office under grant number W911NF-14-1-0359. We acknowledge the use of the supercomputing facilities managed by the Research Computing Group at the University of Massachusetts Boston as well as the University of Massachusetts Green High Performance Computing Cluster.

∎ REFERENCES 1 Chu, S.; Majumdar, A. Opportunities and challenges for a sustainable energy future. Nature 2012, 488, 294–303. 2 Schlapbach, L. Hydrogen-fuelled vehicles. Nature 2009, 460, 809–811. 3 Schlapbach, L.; Z¨uttel, a. Hydrogen-storage materials for mobile applications. Nature 2001, 414, 353–358. 4 Dryer, F. L.; Chaos, M.; Zhao, Z.; Stein, J. N.; Alpert, J. Y.; Homer, C. J. Spontaneous Ignition of Pressurized Releases of Hydrogen and Natural Gas into Air. Combust. Sci. Technol. 2007, 179, 663–694. 5 Alcock, J. L.; Shirvill, L. C.; Cracknell, R. F. Compilation of Existing Safety Data on Hydrogen and Comparative Fuels; 2001; pp 1–15. 6 Yamada, E.; Kitabayashi, N.; Hayashi, A. K.; Tsuboi, N. Mechanism of high-pressure hydrogen auto-ignition when spouting into air. Internation Journal of Hydrogen Energy 2011, 36, 2560–2566. 7 Bates, J.; Hill, N. The Role of Physics in Renewable Energy RD & D. A report produced for the Institute of Physics by Future Energy Solutions; 2005. 8 Lewis, B.; Elbe, G. V. Combustion, Flames and Explosions of Gases, 3rd ed.; Acedemic Press, 1987. 9 Cain, T. Autoignition of hydrogen at high pressure. Combustion and Flame 1997, 111, 124–132. 10 Burke, M. P.; Chaos, M.; Ju, Y.; Dryer, F. L.; Klippenstein, S. J. Comprehensive H2 /O2 kinetic model for highpressure combustion. Int. J. Chem. Kin. 2011, 44, 444– 474. 11 Yang, C. Theory of ignition and auto-ignition. Combust. Flame 1962, 6, 215–225. 12 Moore, J. W.; Pearson, R. G. Kinetics and Mechanism, 3rd ed.; John Wiley & Sons: New York, 1981; p 480. 13 Vlachos, D. G. Stochastic modeling of chemical microreactors with detailed kinetics - Induction times and ignitions of H2 in air. Chem. Eng. Sci. 1998, 53, 157–168. 14 Law, C. K. Combustion Physics, 1st ed.; Cambridge University Press: Cambridge, 2006. 15 Mittal, G.; Sung, C. J.; Yetter, R. A. Autoignition of H2/CO at Elevated Pressures in a Rapid Compression Machine. Int. J. Chem. Kinet. 2006, 38, 516–529. 16 Marin, G.; Yablonsky, G. S. Kinetics of Chemical Reactions, 1st ed.; WILEY-VCH: Singapore, 2011; p 446. 17 Zhao, Z.; Chen, Z.; Chen, S. Correlations for the ignition delay times of hydrogen/air mixtures. Chinese Sci. Bull. 2011, 56, 215–221. 18 S´anchez, A. L.; Fern´andez-Tarrazo, E.; Boivin, P.; Li˜na´ n, A.; Williams, F. A. Ignition time of hydrogen-air diffusion flames. Comptes Rendus - Mec. 2012, 340, 882– 893. 7

ACS Paragon Plus Environment

The Journal of Physical Chemistry

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

19 Zhang, Y.; Huang, Z.; Wei, L.; Zhang, J.; Law, C. K. Experimental and modeling study on ignition delays of lean mixtures of methane, hydrogen, oxygen, and argon at elevated pressures. Combust. Flame 2012, 159, 918–931. 20 Fern´andez-Tarrazo, E.; S´anchez, A. L.; Williams, F. A. Hydrogen-air mixing-layer ignition at temperatures below crossover. Combust. Flame 2013, 160, 1981–1989. 21 Wang, X.; Law, C. K. An analysis of the explosion limits of hydrogen-oxygen mixtures. J. Chem. Phys. 2013, 138, 134305–134317. 22 K´eromn`es, A.; Metcalfe, W.; Heufer, K.; Donohoe, N.; Das, A.; Sung, C.-J.; Herzler, J.; Naumann, C.; Griebel, P.; Mathieu, O.; et al, An experimental and detailed chemical kinetic modeling study of hydrogen and syngas mixture oxidation at elevated pressures. Combust. Flame 2013, 160, 995–1011. 23 de Vijver, R. V.; Vandewiele, N. M.; Bhoorasingh, P. L.; Slakman, B. L.; Khanshan, F. S.; Carstensen, H.-H.; Reyniers, M. F.; Marin, G. B.; West, R. H.; Geem, K. M. V. Automatic Mechanism and Kinetic Model Generation for Gas- and Solution-Phase Processes: A Perspective on Best Practices, Recent Advances, and Future Challenges. Int. J. Chem. Kin. 2015, 47, 199–231. 24 Gao, C. W.; Allen, J. W.; Green, W. H.; West, R. H. Reaction Mechanism Generator: Automatic construction of chemical kinetic mechanisms. Comput. Phys. Commun. 2015, 203, 212–225. 25 Burnham, A. K.; Zhou, X.; Broadbelt, L. J. Critical Review of the Global Chemical Kinetics of Cellulose Thermal Decomposition. Energy & Fuels 2015, 29, 2906– 2918. 26 Donohoe, N.; Heufer, K. A.; Aul, C. J.; Petersen, E. L.; Bourque, G.; Gordon, R.; Curran, H. J. Influence of steam dilution on the ignition of hydrogen, syngas and natural gas blends at elevated pressures. Combust. Flame 2015, 162, 1126–1135. 27 Varga, T.; Nagy, T.; Olm, C.; Zs´ely, I. G.; P´alv¨olgyi, R.; Valk´o,; Vincze, G.; Cserh´ati, M.; Curran, H. J.; Tur´anyi, T. Optimization of a hydrogen combustion mechanism using both direct and indirect measurements. Proc. Combust. Inst. 2015, 35, 589–596. 28 Wagner, A. L.; Yelvington, P. E.; Cai, J.; Green, W. H. Combustion of Synthetic Jet Fuel: Chemical Kinetic Modeling and Uncertainty Analysis. J. Propuls. Power 2017, in press. 29 Bal, K. M.; Neyts, E. C. Direct observation of realistictemperature fuel combustion mechanisms in atomistic simulations. Chem. Sci. 2016, 7, 5280–5286. 30 Duin, A. V.; Kimberly, C. ReaxFF Potential Functions. J. Phys. Chem. A 2001, 1–10. 31 Agrawalla, S.; van Duin, A. C. T. Development and appli-

32

33

34

35

36

37

38

39

40 41

42

43

44

cation of a ReaxFF reactive force field for hydrogen combustion. J. Phys. Chem. A 2011, 115, 960–972. Cheng, T.; Jaramillo-Botero, A.; Goddard, W. A.; Sun, H. Adaptive accelerated ReaxFF reactive dynamics with validation from simulating hydrogen combustion. J. Am. Chem. Soc. 2014, 136, 9434–9442. Chenoweth, K.; van Duin, A. C. T.; Goddard, W. A. ReaxFF reactive force field for molecular dynamics simulations of hydrocarbon oxidation. J. Phys. Chem. A 2008, 112, 1040–1053. Monge-Palacios, M.; Rafatijo, H. On the role of the termolecular reactions 2O2 + H2 → 2HO2 and 2O2 + H2 → H + HO2 + O2 in formation of the first radicals in hydrogen combustion: ab initio predictions of energy barriers. Phys. Chem. Chem. Phys. 2017, 19, 2175–2185. Alaghemandi, M.; Green, J. R. Reactive symbol sequences for a model of hydrogen combustion. Phys. Chem. Chem. Phys. 2016, 18, 2810–2817. Nicholson, S.; Alaghemandi, M.; Green, J. Learning the mechanisms of chemical disequilibria. J. Chem. Phys. 2016, 145, 084112–084119. Li, J.; Zhao, Z.; Kazakov, A.; Dryer, F. L. An updated comprehensive kinetic model of hydrogen combustion. Int. J. Chem. Kinet. 2004, 36, 566–575. Kylasa, S.; Aktulga, H.; Grama, A. PuReMD-GPU: A reactive molecular dynamics simulation package for GPUs. J. Comput. Phys. 2014, 272, 343–359. Senftle, T. P.; Hong, S.; Islam, M.; Kylasa, S. B.; Zheng, Y.; Shin, Y. K.; Junkermeier, C.; Engel-herbert, R.; Janik, M. J.; Aktulga, H. M.; et al, The ReaxFF reactive force-field : development, applications and future directions. Comput. Mater. 2016, 2, 15011–15025. CHEMKIN-PRO 15131, Reaction Design: San Diego, 2013. Smith, G. P.; Golden, D. M.; Frenklach, M.; Moriarty, N. W.; Eiteneer, B.; Goldenberg, M.; Bowman, T. C.; Hanson, R. K.; Song, S.; Jr., W. C. G.; et al, GRI-MECH 3.0. http://www.me.berkeley.edu/ gri_mech/. Gord, J. R.; Meyer, T. R.; Roy, S. Applications of Ultrafast Lasers for Optical Measurements in Combusting Flows. Annu. Rev. Anal. Chem. 2008, 1, 663–687. Wrzesinski, P. J.; Stauffer, H. U.; Schmidt, J.; Roy, S.; Gord, J. R. Single-Shot Thermometry via Femtosecond, Fully Resonant, Electronically Enhanced CARS (FREECARS). 32nd AIAA Aerodyn. Meas. Technol. Gr. Test. Conf. 2016, 41, 2021–2024. Mohammadi, A.; Shioji, M.; Nakai, Y.; Ishikura, W.; Tabo, E. Performace and combustion characteristics of a direct injection SI hydrogen engine. Hydrogen Energy 2006, 32, 296–304.

8

ACS Paragon Plus Environment

Page 8 of 9

Page 9 of 9

The Journal of Physical Chemistry

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 ACS Paragon Plus Environment