- Email: [email protected]
CVI processes
Carbon 40 (2002) 675–683
A reduced reaction model for carbon CVD/ CVI processes Narayana Birakayala, Edward A. Evans* Chemical Engineering Department, University of Akron, Akron, OH 44325 -3906, USA Received 26 June 2000; accepted 3 June 2001
Abstract Carbon–carbon composites are produced by chemical vapor deposition / chemical vapor infiltration (CVD/ CVI) processes. Models of carbon–carbon composite production processes will help reduce production costs. Reliable process models must, however, include details of the gas phase kinetics in order to identify optimal conditions. We have combined detailed gas phase kinetics, surface kinetics, and a pore closure model to predict pore geometry changes with respect to time. To determine the dominant gas phase kinetics, we reduced a large set of reactions to a minimal set using a sensitivity, rate, and dimensional analysis approach. These robust and relatively fast techniques can be used under a variety of conditions, including those within the pores of the composite. The process model shows that the deposition profile depends on the kinetic model chosen. Using the more realistic reaction model, conditions for uniform, or inside-out, densification can be suggested. 2002 Elsevier Science Ltd. All rights reserved. Keywords: A. Pyrolytic carbon; B. Chemical vapor deposition; Chemical vapor infiltration
1. Introduction Carbon–carbon (C–C) composites are a class of synthetic, pure carbon materials consisting of a carbon matrix reinforced by carbon fibers. C–C composites combine the unique properties of bulk carbon (e.g., low-density, high thermal conductivity, and low thermal expansion) with the excellent mechanical properties of carbon fibers [1]. C–C composites are used in a variety of high temperature structural applications. The high thermal and mechanical shock resistance and high stiffness of these composites contribute to their use as a material for advanced aerospace applications; applications include rocket nose cones and nozzles, heat shields for re-entry vehicles, and heat sinks and radiators for satellite applications. C–C composite materials are also used in jet engine rotors and stators, where higher operating temperatures increase the efficiency of the heat engines. C–C composites can be used to significantly reduce engine weight and size which results in lower fuel consumption [2]. The excellent friction and wear characteristics of these composites combined with the high thermal conductivity make them particularly useful for aircraft brake discs. The high thermal conductivity and thermal capacity of C–C *Corresponding author. Fax: 11-330-972-5856. E-mail address: [email protected] (E.A. Evans).
composites helps in absorbing and conducting away large quantities of heat. Production costs of C–C composites are high which limits their use in some markets. Our objective, therefore, is to develop models of the C–C composite production process that will help reduce production costs. We have focused our efforts on understanding the details of the kinetics involved during C–C composite manufacturing; identification of the dominant gas phase kinetics for C–C composite fabrication is required for improved modeling of the production process [1,3–8].
2. Background Carbon–carbon composites are produced by chemical vapor deposition / chemical vapor infiltration (CVD/ CVI) processes. In the CVD/ CVI process the densification of a porous carbon preform is achieved through carbon deposition. Model development for improving CVD/ CVI technology has received much attention in the literature over the last decade. Most of the published research on the modeling of C–C composites has emphasized the modeling of flow through porous media [2,10,11]. These models, which use sophisticated pore and flow models, typically represent the gas phase kinetics by a single overall reaction; more generally, the reactions have been modeled as reaching equilibrium. Simple reaction kinetics can be
0008-6223 / 02 / $ – see front matter 2002 Elsevier Science Ltd. All rights reserved. PII: S0008-6223( 01 )00184-1
676
N. Birakayala, E. A. Evans / Carbon 40 (2002) 675 – 683
used when the time evolution of the species distribution is not considered. The equilibrium assumption is valid for long residence times and a well-mixed reactor. To obtain ideal conditions for densification, however, the reaction times should be longer than the transport times within the system — a reaction-limited regime. Using a single gas phase reaction and surface reaction for a generic CVD/ CVI process, Middleman has shown that inside-out densification is possible under a reaction-limited regime [9]. In a reaction-limited regime, the growth species have time to diffuse before reacting which results in a more uniform deposition profile. In a reaction-limited regime the details of the kinetics must be known. Most CVD/ CVI processes are carried out under conditions closer to plug flow rather than well mixed [12] (Peclet number |1). Under these conditions, the species distribution should be predicted by the kinetics rather than by an equilibrium point. The CVD/ CVI modeling of titanium composites by Moene et al. [13] showed that the reaction rate profile within the porous region depended strongly on the kinetic model chosen. Accurate kinetics are therefore required to predict the species distribution within the porous structure of the carbon preform. To model the composite production process, a three-dimensional numerical model may be required. There are limitations to the extent of the reaction mechanism that can be used when three-dimensional numerical models are required. Most three-dimensional models can only handle a small set of reactions [20]. This set of reactions should be chosen in a methodical, robust fashion. Furthermore, multiple techniques for selecting the reactions should be used to verify the importance of those reactions. Finally, techniques for selecting reactions should be based upon both the characteristic reaction and transport times within the system. A model based on a reliable kinetic model can be used to improve process efficiency and predict the composite properties with more confidence. Huttinger and co-workers have recently published several articles which detail the kinetics involved during C–C composite production [3–8,14–18]. It is difficult, however, to understand the reaction pathways and the dominant kinetics from a large set of reactions. We have used sensitivity and rate analyses to identify the dominant reactions and then used them to model the deposition of carbon into a C–C composite. The analysis is based on the kinetic parameters in the literature and not on parameters estimated or measured by the authors.
3. Model development Three components of the model that have been developed are the gas phase kinetics, the surface kinetics, and the pore closure model. All three components are required to predict the pore geometry changes with respect to time. Our focus is on developing the gas phase kinetics; our
results can be applied to other flow and deposition models. Future work will focus on improving the surface kinetic model and the pore closure model.
4. Gas phase kinetics To determine the dominant gas phase kinetics, we reduced a large set of reactions to a minimal set using a sensitivity and rate analysis approach. We started with a set of gas phase reactions and their associated rate parameters from the Gas Research Institute (GRI) database; the NIST database was used to check the rate parameters [19,20]. We included all of the species in the mechanisms used by Becker and Huttinger [3] but used the reactions and rate parameters found in the GRI and NIST databases. Experimental results in the literature and our own thermodynamic calculations indicate the importance of hydrocarbon species containing up to six carbon atoms [21]; other aromatic species and higher hydrocarbons were left out. To include all of these species, the complete mechanism is made up of 47 reversible reactions with 19 species. The rate constant is expressed in a modified Arrhenius form in which the pre-exponential term includes a temperature dependence. Pressure dependency is included for some of the reactions and Troe expressions are used for the calculation of the pressure-dependent rate coefficients [19]. Thermodynamic properties are taken from the CHEMKIN thermodynamic database [22]. The predictions of the general gas phase mechanism at long times were checked against equilibrium predictions and the predictions at short times were assumed to be correct. Reaction sets of this size (greater than 20) cannot be efficiently used with current finite difference modeling packages. Furthermore, it is difficult to visualize the reaction pathways with such a large number of reactions. The dominant reactions must be identified. To identify the dominant reactions we have applied sensitivity, rate, and dimensional analyses; sensitivity and rate analysis techniques have been used previously to identify the dominant reactions in a large set [23–26]. Briefly, a sensitivity analysis approach selects reactions based on the response of the model solution to changes in the kinetic parameters (e.g., activation energy); rate analysis selects reactions based on their contribution to the overall production rate of the important species. We have also used a dimensional analysis approach based on the Damkohler number to determine the important reactions within the porous region. The Damkohler number is the ratio of the reaction rate constant to the diffusion coefficient of the reactants. Reactions with the largest Damkohler numbers are selected because they are typically the fastest reactions in a parallel set and would be the dominant terms in a species balance equation. To reduce the large reaction mechanism to a smaller set, we established criteria for a good comparison of the
N. Birakayala, E. A. Evans / Carbon 40 (2002) 675 – 683
reduced mechanism to the large mechanism. The chosen reactions represent the general mechanism only when the deviations of the predictions of the reduced mechanism from the complete mechanism are below 10% for the dominant species. The dominant species depend on the reaction conditions chosen. We have chosen the same conditions as those used in the work of Bammidipati et al. so that we could compare our model prediction to experimental data [12]. The data were collected at pressures between 10 and 40 Torr (5.3310 3 Pa) and temperatures between 1273 and 1373 K. An isothermal zone was maintained around the substrate. Deposition rates on the substrate were measured at different pressures, temperatures, flow rates (50–400 SCCM), and diluent concentration (0–87.5%). Diluent gases used in the study included hydrogen, nitrogen, and argon. At these conditions the dominant hydrocarbon species are acetylene, benzene, and ethylene; hydrogen (H 2 ) makes up nearly 83% of the reactant mixture. The reduced reaction sets are considered
677
reasonable as long as they predict the concentrations of the three main hydrocarbon species to within 10% of the general mechanism. Sensitivity, rate, and Damkohler analyses selected the same reactions for the bulk gas phase adjacent to the growing surface. Table 1 shows the dominant reactions in decreasing order of importance. Fig. 1 shows the error of the reduced mechanism as a function of the number of reactions within the mechanism. The concentration of benzene in the reactor becomes significant only at high residence times. At lower residence times, acetylene is the dominant species contributing to the growth of carbon. At conditions where acetylene is the dominant growth species, the reactions involving benzene and its precursors are not important. Reactions 1–10 in Table 1 are required to accurately predict the C 2 species concentrations. The errors lie well below 10% and can be used to predict the concentrations of the dominant C 2 species. Further reduction from these ten reactions leads to errors much greater
Table 1 Reduced mechanism deduced from sensitivity, rate, and Damkohler analysis Arrhenius expression (mol, cm, s) Reactions to predict C1 –C2 1. H1CH 3 (1M)⇔CH 4 (1M) 2. H1CH 4 ⇔CH 3 1H 2 3. H1C 2 H 3 ⇔C 2 H 3 (1M) 4. H1C 2 H 3 ⇔H 2 1C 2 H 2 5. H1C 2 H 4 1(M)⇔C 2 H 5 (1M) 6. H1C 2 H 4 ⇔C 2 H 3 1H 2 7. H1 C 2 H 6 ⇔C 2 H 5 1H 2 8. 2CH 3 (1M)⇔C 2 H 6 (1M) 9. 2CH 3 ⇔H1C 2 H 5 10. CH 3 1C 2 H 4 ⇔C 2 H 3 1CH 4
1.27310 16 T 20.6 exp(192.7 / T) 6.6310 8 T 1.6 exp(5454.9 / T) 5.6310 12 exp(1207.7 / T) 3310 13 1.08310 12 T 0.5 exp(915.9 / T) 1.32310 6 T 2.5 exp(6159.4 / T) 1.15310 8 T 1.9 exp(3789.3 / T) 2.12310 16 T 1.0 exp(312 / T) 4.99310 12 T 20.1 exp(5334.1 / T) 2.27310 5 T 2.0 exp(4629.6 / T)
Reactions to predict C2 –C6 species 11. C 2 H1H 2 ⇔H1C 2 H 2 12. H1C 2 H 4 (1M)⇔C 2 H 2 (1M) 13. C 4 H 3 1C 2 H 3 ⇔C 6 H 6 14. C 3 H 3 1C 3 H 3 ⇔C 6 H 6 15. C 3 H 4 ⇔C 3 H 3 1H 16. CH 3 1C 2 H 2 ⇔C 3 H 4 1H 17. C 2 H1C 2 H 2 ⇔C 4 H 3 18. C 2 H 2 1C 2 H 2 ⇔C 4 H 4 19. C 2 H 3 1C 2 H 2 ⇔C 4 H 4 1H 20. C 2 H1C 2 H 4 ⇔C 4 H 4 1H 21. CH 3 1C 2 H 6 ⇔C 2 H 5 1CH 4 22. C 3 H 4 1C 3 H 3 ⇔C 6 H 6 1H 23. C 2 H 4 (1M)⇔H 2 1C 2 H 2 (1M) 24. C 4 H 4 1C 2 H 2 ⇔C 6 H 6 25. C 3 H 4 1C 2 H⇔C 3 H 3 1C 2 H 2 26. H1C 2 H 5 (1M)⇔C 2 H 6 (1M) 27. H1C 2 H 3 (1M)⇔C 2 H 4 (1M) 28. H1C 2 H 5 ⇔H 2 1C 2 H 4
4.07310 5 T 2.4 exp(100.6 / T) 1.00310 17 T 21.0 2.87310 14 exp(412.6 / T) 3.0310 11 1.00310 17 exp(35225.4 / T) 6.74310 19 T 22.1 exp(15896.7 / T) 1.00310 13 2.45310 14 exp(23349.4 / T) 2.0310 12 exp(25 16.1 / T) 1.21310 13 6.14310 6 T 1.7 exp(5258.7 / T) 2.2310 11 exp(1006.4 / T) 8.0310 12 T 0.4 exp(44670.9 / T) 4.47310 11 exp(15141.91T) 1.00310 13 5.21310 17 T 21.0 exp(795.1 / T) 6.08310 12 T 0.3 exp(140.9 / T) 2.0310 12
Reactions are listed in order of decreasing importance. Reactions with a third body (1M) have rate constants that are pressure dependent. The pressure dependency can be expressed in either the Lindemann or Troe form. Rate parameters were taken from the Gas Research Institute database and verified using the NIST database [19,20].
678
N. Birakayala, E. A. Evans / Carbon 40 (2002) 675 – 683
Fig. 1. Error in the species distribution predicted by the mechanism shown in Table 1. Error is based on comparison to the full mechanism (47 reactions) prediction.
than 10%; the predictions of the reduced models proposed by Bammidipati et al., for example, differed from our general mechanism prediction by up to 100% [12]. To predict the benzene concentrations accurately, reactions 11–28 are also required. If an error of nearly10% is acceptable, reactions 1–16 can be used to predict the concentrations of species up to and including benzene. Fig. 2 shows a reaction network based on our reduced mechanism for the pyrolysis of methane at pressures between 10 and 40 Torr (5.3310 3 Pa) and temperatures between 1273 and 1373 K; the larger arrows indicate the dominant path between methane and the products for this set of conditions. Sensitivity and rate analyses were used to validate the use of Damkohler numbers for selecting reactions in the gas phase. Sensitivity and rate analyses are based on the kinetics alone while the Damkohler number incorporates the influence of transport. To select reactions within the
Fig. 2. Network showing the reaction pathways during the pyrolysis of methane at CVD/ CVI conditions. The larger arrows indicate the dominant reaction path. The network is based on reactions shown in Table 1.
porous substrate, therefore, we used the Damkohler number as a guide for selecting the dominant reactions [14]. The dominant gas phase reaction set inside the porous region was the same as the set shown in Table 1.
5. Surface kinetics In order to predict the change in porosity of the C–C composite, surface reactions are included; parameters for surface reactions are not available at conditions for the manufacture of C–C composites. To estimate the kinetic parameters for the surface reactions we have used the deposition rate data of Bammidipati et al. in order to fit the surface reaction parameters; the reported deposition rates were measured on non-porous substrates for different experimental conditions [12]. We tested the model by comparing predicted and measured deposition rates; we did not compare gas phase compositions. For a detailed mechanism, the potential species, which react on the surface, should first be identified. To identify the possible growth species, we applied the gas phase kinetic model to simulate the gas phase conditions during the carbon deposition experiments of Bammidipati et al. [12]. These authors predict significant velocity gradients in the reactor in both the radial and axial directions. Allowing for such velocity gradients, our simulations were performed at an average temperature of 1395 K and a pressure of 40 Torr (5.3310 3 Pa). The
N. Birakayala, E. A. Evans / Carbon 40 (2002) 675 – 683 Table 2 Surface mechanism for carbon deposition during CVD/ CVI 1. C 2 H 2 12C(S)⇒2C(S)12C(D)1H 2 2. C 6 H 6 16C(S)⇒6C(S)16C(D)13H 2 3. C 2 H 4 12C(S)⇒2C(S)12C(D)12H 2 4. H1CH(S)⇔C(S)1H 2 5. H1C(S)⇔CH(S)
679
sticking probabilities for the adsorbing reactant. The sticking probabilities in our model were obtained by matching the simulation results to the measured deposition rates of Bammidipati et al. [12]. Sticking probability, which is a function of coverage (u ), is defined as Number of molecules that stick or react g (u ) 5 ]]]]]]]]]]]] Number of molecules that impinge on a surface (1)
model gives the distribution of various hydrocarbons along the length of the reactor at the given conditions of temperature, pressure, flow rate, and initial composition [22]. The species distribution was obtained at various flow rates of methane and the amount of carbon available for deposition from respective species was calculated as a function of flow rate. Based on this analysis, the major contribution for carbon deposition is from methane, acetylene, ethylene, benzene, and ethane. To be classified as a growth species, apart from the gas phase abundance, the species should also have a high reaction probability (Ks ). In recent work, Becker and Huttinger suggested acetylene, ethylene, and benzene as the growth species based on their respective reaction probabilities [5]. Based on these two analyses, we have selected C 2 H 2 , C 2 H 4 , and C 6 H 6 as the major growth species participating in the surface mechanism. The resulting surface mechanism is shown in Table 2. The rate constants for the surface reactions are given as
Fig. 3 shows the comparison of the model predictions with the experimental growth rates at 40 Torr (5.3310 3 Pa) and 1373 K. The simulation results are qualitatively consistent with the observed growth rates. The sticking probabilities used were:
gC 2 H 2 55310 25 gC 2 H 4 55310 27 gC 6 H 6 54310 24 (2) The sticking probability is a function of surface coverage. Various models have been proposed to relate the variation in the sticking probability with coverage and / or pressure. As the pressure increases the number of open active sites decreases as a result of increased coverage by any adsorbed species. With fewer open sites for adsorption the surface reaction rates decrease. The decrease in the surface reaction rates with pressure can be accommodated through the sticking probability; an alternative approach would be
Fig. 3. Comparison of the model to experimental data. P540 torr (5.3310 3 Pa), T51373 K, CH 4 : 100%. Experimental data with error bars are taken from Bammidipati et al. [12].
N. Birakayala, E. A. Evans / Carbon 40 (2002) 675 – 683
680
Table 3 Sticking coefficients used to fit the model predictions to the experimental rates Pressure (Torr)
C2H2 C2H4 C6H6
10 (1.3310 3 Pa)
40 (5.3310 3 Pa)
5310 24 4310 24 4310 23
0.5310 24 4310 27 4310 24
to use additional surface reactions to account for the blocking of active surface sites. The simple Langmuir model predicts a decrease of the sticking coefficient with increasing coverage.
g (u ) ]] 5 (1 2 u )n g (0)
(3)
Table 4 Effective binary diffusivities (EBD) for species in methane Species
EBD (m 2 / s)
H H2 CH 4 C2H2 C2H6
6.9310 22 4.5310 22 1.72310 22 1.38310 22 1.25310 22
the concentration profile of the species in the pores are fast compared to the changes in pore geometry, exist in the substrate [1]. This assumption requires gas phase species reaching equilibrium at relatively short time scales compared to the densification time. This is true in ICVI processes, where the densification times typically are 100– 300 h. From the pseudo steady-state hypothesis, the continuity equation reduces to
OR K
where g (0) is the sticking coefficient at zero coverage and n is the number of sites that is needed to hold the adsorbate. For n 5 1, Eq. (3) shows a linear decrease and for n . 1 it shows a nonlinear decrease. Bammidipati et al. [12] have reported a decrease in the sticking coefficients with increasing pressure to account for the increased coverage. For fitting the data at pressures lower than 40 torr, the sticking coefficients are increased. The sticking coefficients used in the current work are given in Table 3.
6. Pore model A cylindrical pore model has been chosen to predict the changes in porosity as deposition occurs. The pore model can be used to predict the concentration and deposition profiles within the pore. The cylindrical pore model is sufficient for the current work, since the objective of the work is to understand qualitatively the difference in using the various kinetic models. We have followed the development of Moene et al. to set-up the pore diffusion model [13]. The first assumption is that the diffusive flux is dominant in the substrate. This assumption is valid for the ICVI processes, where the thermal and flow gradients are absent. The flux term (Ni ) can be written as Ni 5 2 De,i (=Ci )
(4)
The mixture averaged diffusivity, Dei , is equal to the effective binary diffusivity of the species with respect to methane, De,i,CH 4 . The effective binary diffusivities used in this model for some representative species are given in Table 4. The effective diffusivity is calculated as the harmonic mean of binary and Knudsen diffusivities [13,14]. Pseudo steady-state conditions, at which the changes in
=? (ANi ) 5 (VTu )
OR S
Gi
1 (VT a)
k 51
si
(5)
k 51
where A is the cross-sectional area, Ni is the flux of i, VT is the total volume, u is the porosity, a is the internal surface area per unit volume, and R Gi and R Si are the net rate of decomposition of i in the gas phase and on the surface, respectively. Heat effects due to the gas phase reactions and carbon deposition are assumed to be small and a uniform temperature profile exists throughout the pore. Since very low deposition rates are maintained in the CVI reactors, the low conversion of the source gases within the pores does not affect the temperatures. For the case of cylindrical pores representing the pores in the substrate, the flow is modeled as a one-dimensional problem. The aspect ratios of the pore, defined as the ratio of the length over the diameter of the pore, vary from 50 to 250. Hence, the diffusion in the radial direction can be neglected and a uniform concentration profile over the cross-section can be assumed. Incorporating the above assumptions and combining Eqs. (4) and (5), the mass balance on a component i, within a cylindrical pore of radius r, over a small length (DX) gives ≠C S D 1 D (pr )S] D ≠x
≠Ci 2 Dei (p r 2 ) ] ≠x
2
x
ei
i
x 1Dx
5 R Gi (p r 2 Dx) 1 R Si (2p r Dx)
(6)
For an infinitesimal element, DX →0, Eq. (6) becomes, ≠(2Dei r 2 ≠Ci / ≠x) 2 ]]]]] 2 R Gi (r ) 2 R Si (2r) 5 0 ≠x
(7)
The equation is further simplified by assuming that the diffusivity and the pore radius do not change considerably in a very thin section. If the pore is divided into a fine grid with several thin sections, the diffusivity and radius remain constant within the section. However, the diffusivity and
N. Birakayala, E. A. Evans / Carbon 40 (2002) 675 – 683
the radius differ from section to section. The equation is solved at different instances of time independently of the equation describing the growth rate. With these simplifications, d 2 Ci 2 Dei r 2 ]] 2 R Gi (r 2 ) 2 R Si (2r) 5 0 dx 2
(8)
The radius of the pore as a function of the deposition rate is deduced from Eq. (9). The radius as a function of the pore length and time is determined by solving Eqs. (8) and (9) simultaneously.
O
M I ≠r ] 5 2 ]s R ≠t rs i 51 si
(9)
The boundary conditions for the problem are as follows: Ci (0,t) 5 Ci,bulk , the concentration of the species i in the bulk; Ci (L,t) 5 Ci,bulk and r(x,0) 5 r0, initial radius of the pore. The equations are solved with the boundary conditions to obtain the variation of radius as a function of the pore length and time.
7. Discussion For the comparison of different kinetic mechanisms, a pore aspect ratio of 100 has been used to distinguish between the mechanisms. The differences between the model predictions for two different kinetic mechanisms are
681
shown in Fig. 4; the densification profile is shown for the same time step in each case. The overall kinetic equation results are determined using the direct deposition of methane while the detailed kinetics combine the mechanisms given in Table 1 (reactions 1–16) and Table 2. The detailed kinetic model predicts a non-uniform deposition profile. A more uniform profile is predicted with the single reaction when methane is depositing carbon directly without the formation of any intermediates. A higher diffusivity of methane compared to that of acetylene (|19%) at the same conditions may have contributed to the uniformity. Also, because of the higher concentration of methane at the surface predicted by the detailed mechanism, the deposition within the pore occurs at a faster rate leading to pore closure at much shorter times. Thus, a completely different densification profile is obtained when the kinetics are changed, as expected. The deposition reaction results in a concentration gradient within the pore decreasing the concentrations of the growth species away from the opening (boundary) of the pore. The matrix forms rapidly at the surface leading to residual porosity. A uniform densification of the composite is possible, however, by operating at conditions of low reaction rates and high transport rates. Typically, low temperatures are employed in CVI reactors to reduce the density gradients. The results of Middleman [9] suggest that inside-out densification can be achieved by making use of the methane pyrolysis kinetics. If the growth species, which deposit carbon, were produced within the pore at a rate higher than its dissociation on the surface, a
Fig. 4. Final pore radii for different kinetic models: a single deposition reaction and the kinetics shown in Tables 1 and 2. The initial pore radius was 2.24310 26 m with an aspect ratio5100. The dimensionless pore lengths 0 and 1 are the two ends of the cylindrical pore.
682
N. Birakayala, E. A. Evans / Carbon 40 (2002) 675 – 683
Fig. 5. Initial evolution of pore structure resulting from different acetylene concentrations (Cg [5]kg-mol m 3 ) near the pore boundary. The nominal residence time within the reactor was 16 h. The initial pore radius was 2.24310 26 m with an aspect ratio550. The dimensionless pore lengths 0 and 1 are the two ends of the cylindrical pore.
uniform densification would be obtained. For example, the simulation is carried out at low acetylene concentrations at the entrance to the pore. Fig. 5 shows the deposition profiles for low acetylene concentrations at very short times. At very low or no acetylene at the pore opening, an inside-out densification is achieved. As the acetylene concentration is increased, the profile flattens and finally residual porosity sets in; the concentration gradients within the pore increase as the concentration of the growth precursors increase in the external gas phase. This result suggests that inside-out densification may be achieved when the growth precursor (acetylene) is produced within the pore. The low acetylene concentrations in the reactor are possible at high flow rates or with high diluent (H 2 , N 2 , Ar) concentrations; however, nitrogen and argon have been shown to have a negative influence on the deposition process [14,17]. Simulations with a single gas phase reaction would not result in similar predictions as there is no production of any intermediates depositing carbon and the source gas is assumed to deposit carbon directly. At atmospheric pressures and infinitely long times, Benzinger and Huttinger have demonstrated that the degree of pore filling reaches a maximum with respect to the partial pressure of the growth precursor (methane). As the partial pressure of methane is increased above 150 Torr there is a rapid decrease in the degree of pore filling [16]; the trend above 150 Torr agrees qualitatively with the results presented here. Our future work will investigate the
possibility of such a maximum at the conditions used in this work.
8. Conclusions The importance of accurately representing the kinetics during CVD/ CVI of C–C composites is clear. We have demonstrated that the dominant kinetics can be identified and used to predict deposition profiles within the pores. Furthermore, conditions for uniform or inside-out densification can be suggested based on the kinetics. A more rigorous pore model is required to estimate the fabrication times and substrate properties as a function of fabrication time. Likewise, more experimental data are required to validate the model predictions. The same approach can be used to model other composite materials.
References [1] Vaidyaraman S, Lackey WJ, Agrawal PK, Starr TL. Carbon 1996;34:1123. [2] McAllister P, Wolf EE. AIChE J 1993;39:1196. [3] Becker A, Huttinger KJ. Carbon 1998;36:177. [4] Becker A, Huttinger KJ. Carbon 1998;36:201. [5] Becker A, Huttinger KJ. Carbon 1998;36:213. [6] Becker A, Huttinger KJ. Carbon 1998;36:225.
N. Birakayala, E. A. Evans / Carbon 40 (2002) 675 – 683 [7] Bruggert M, Hu Z, Huttinger KJ. Carbon 1999;37:2021. [8] Antes J, Hu Z, Zhang W, Huttinger KJ. Carbon 1999;37:2031. [9] Middleman S. J Mater Res 1989;4:1513. [10] Sheldon W. J Mater Res 1990;5:2729. [11] Lackey WJ, Vaidyaraman S, Beckloff BN, Moss TS, Lewis JS. J Mater Res 1998;13:2251. [12] Bammidipati S, Stewart GD, Elliott Jr. JR, Gokoglu SA, Purdy MJ. AIChE J 1996;42:3123. [13] Moene R, Dekker JP, Makkee M, Schoonman J, Mouljin JA. J Electrochem Soc 1994;141:282. [14] Benzinger W, Huttinger KJ. Carbon 1996;34:1465. [15] Benzinger W, Huttinger KJ. Carbon 1998;36:1033. [16] Benzinger W, Huttinger KJ. Carbon 1999;37:181. [17] Benzinger W, Huttinger KJ. Carbon 1999;37:931. [18] Benzinger W, Huttinger KJ. Carbon 1999;37:841. [19] Gas Research Institute Database, version 2.11. http: / / www.gri.org
683
[20] Westley F, Frizzell DH, Herron JT, Hampson RF, Mallard WG. NIST Chemical Kinetics Database Version 5.0, NIST Standard Reference Data, Gaithersburg, MD 20899, USA, 1992. [21] Murphy DB, Carroll RW, Klonowski JE. Carbon 1997;35:1819. [22] Kee RJ, Rupley FM, Miller JA. The CHEMKIN Thermodynamic Database, Sandia National Laboratories document SAND87-8215B, 1990. [23] Rota R, Bonini F, Servida A, Morbidelli M, Carra S. Chem Eng Sci 1994;49:4211. [24] Brock EE, Savage PE, Barker JR. Chem Eng Sci 1998;53:857. [25] Wolden C, Gleason KK, Howard JB. Combustion Flame 1994;96:75. [26] Frenklach M. Modeling. In: Gardiner WC, editor, Combustion Chemistry, Berlin: Springer, 1984, Chapter 7.
A reduced reaction model for carbon CVD/ CVI processes Narayana Birakayala, Edward A. Evans* Chemical Engineering Department, University of Akron, Akron, OH 44325 -3906, USA Received 26 June 2000; accepted 3 June 2001
Abstract Carbon–carbon composites are produced by chemical vapor deposition / chemical vapor infiltration (CVD/ CVI) processes. Models of carbon–carbon composite production processes will help reduce production costs. Reliable process models must, however, include details of the gas phase kinetics in order to identify optimal conditions. We have combined detailed gas phase kinetics, surface kinetics, and a pore closure model to predict pore geometry changes with respect to time. To determine the dominant gas phase kinetics, we reduced a large set of reactions to a minimal set using a sensitivity, rate, and dimensional analysis approach. These robust and relatively fast techniques can be used under a variety of conditions, including those within the pores of the composite. The process model shows that the deposition profile depends on the kinetic model chosen. Using the more realistic reaction model, conditions for uniform, or inside-out, densification can be suggested. 2002 Elsevier Science Ltd. All rights reserved. Keywords: A. Pyrolytic carbon; B. Chemical vapor deposition; Chemical vapor infiltration
1. Introduction Carbon–carbon (C–C) composites are a class of synthetic, pure carbon materials consisting of a carbon matrix reinforced by carbon fibers. C–C composites combine the unique properties of bulk carbon (e.g., low-density, high thermal conductivity, and low thermal expansion) with the excellent mechanical properties of carbon fibers [1]. C–C composites are used in a variety of high temperature structural applications. The high thermal and mechanical shock resistance and high stiffness of these composites contribute to their use as a material for advanced aerospace applications; applications include rocket nose cones and nozzles, heat shields for re-entry vehicles, and heat sinks and radiators for satellite applications. C–C composite materials are also used in jet engine rotors and stators, where higher operating temperatures increase the efficiency of the heat engines. C–C composites can be used to significantly reduce engine weight and size which results in lower fuel consumption [2]. The excellent friction and wear characteristics of these composites combined with the high thermal conductivity make them particularly useful for aircraft brake discs. The high thermal conductivity and thermal capacity of C–C *Corresponding author. Fax: 11-330-972-5856. E-mail address: [email protected] (E.A. Evans).
composites helps in absorbing and conducting away large quantities of heat. Production costs of C–C composites are high which limits their use in some markets. Our objective, therefore, is to develop models of the C–C composite production process that will help reduce production costs. We have focused our efforts on understanding the details of the kinetics involved during C–C composite manufacturing; identification of the dominant gas phase kinetics for C–C composite fabrication is required for improved modeling of the production process [1,3–8].
2. Background Carbon–carbon composites are produced by chemical vapor deposition / chemical vapor infiltration (CVD/ CVI) processes. In the CVD/ CVI process the densification of a porous carbon preform is achieved through carbon deposition. Model development for improving CVD/ CVI technology has received much attention in the literature over the last decade. Most of the published research on the modeling of C–C composites has emphasized the modeling of flow through porous media [2,10,11]. These models, which use sophisticated pore and flow models, typically represent the gas phase kinetics by a single overall reaction; more generally, the reactions have been modeled as reaching equilibrium. Simple reaction kinetics can be
0008-6223 / 02 / $ – see front matter 2002 Elsevier Science Ltd. All rights reserved. PII: S0008-6223( 01 )00184-1
676
N. Birakayala, E. A. Evans / Carbon 40 (2002) 675 – 683
used when the time evolution of the species distribution is not considered. The equilibrium assumption is valid for long residence times and a well-mixed reactor. To obtain ideal conditions for densification, however, the reaction times should be longer than the transport times within the system — a reaction-limited regime. Using a single gas phase reaction and surface reaction for a generic CVD/ CVI process, Middleman has shown that inside-out densification is possible under a reaction-limited regime [9]. In a reaction-limited regime, the growth species have time to diffuse before reacting which results in a more uniform deposition profile. In a reaction-limited regime the details of the kinetics must be known. Most CVD/ CVI processes are carried out under conditions closer to plug flow rather than well mixed [12] (Peclet number |1). Under these conditions, the species distribution should be predicted by the kinetics rather than by an equilibrium point. The CVD/ CVI modeling of titanium composites by Moene et al. [13] showed that the reaction rate profile within the porous region depended strongly on the kinetic model chosen. Accurate kinetics are therefore required to predict the species distribution within the porous structure of the carbon preform. To model the composite production process, a three-dimensional numerical model may be required. There are limitations to the extent of the reaction mechanism that can be used when three-dimensional numerical models are required. Most three-dimensional models can only handle a small set of reactions [20]. This set of reactions should be chosen in a methodical, robust fashion. Furthermore, multiple techniques for selecting the reactions should be used to verify the importance of those reactions. Finally, techniques for selecting reactions should be based upon both the characteristic reaction and transport times within the system. A model based on a reliable kinetic model can be used to improve process efficiency and predict the composite properties with more confidence. Huttinger and co-workers have recently published several articles which detail the kinetics involved during C–C composite production [3–8,14–18]. It is difficult, however, to understand the reaction pathways and the dominant kinetics from a large set of reactions. We have used sensitivity and rate analyses to identify the dominant reactions and then used them to model the deposition of carbon into a C–C composite. The analysis is based on the kinetic parameters in the literature and not on parameters estimated or measured by the authors.
3. Model development Three components of the model that have been developed are the gas phase kinetics, the surface kinetics, and the pore closure model. All three components are required to predict the pore geometry changes with respect to time. Our focus is on developing the gas phase kinetics; our
results can be applied to other flow and deposition models. Future work will focus on improving the surface kinetic model and the pore closure model.
4. Gas phase kinetics To determine the dominant gas phase kinetics, we reduced a large set of reactions to a minimal set using a sensitivity and rate analysis approach. We started with a set of gas phase reactions and their associated rate parameters from the Gas Research Institute (GRI) database; the NIST database was used to check the rate parameters [19,20]. We included all of the species in the mechanisms used by Becker and Huttinger [3] but used the reactions and rate parameters found in the GRI and NIST databases. Experimental results in the literature and our own thermodynamic calculations indicate the importance of hydrocarbon species containing up to six carbon atoms [21]; other aromatic species and higher hydrocarbons were left out. To include all of these species, the complete mechanism is made up of 47 reversible reactions with 19 species. The rate constant is expressed in a modified Arrhenius form in which the pre-exponential term includes a temperature dependence. Pressure dependency is included for some of the reactions and Troe expressions are used for the calculation of the pressure-dependent rate coefficients [19]. Thermodynamic properties are taken from the CHEMKIN thermodynamic database [22]. The predictions of the general gas phase mechanism at long times were checked against equilibrium predictions and the predictions at short times were assumed to be correct. Reaction sets of this size (greater than 20) cannot be efficiently used with current finite difference modeling packages. Furthermore, it is difficult to visualize the reaction pathways with such a large number of reactions. The dominant reactions must be identified. To identify the dominant reactions we have applied sensitivity, rate, and dimensional analyses; sensitivity and rate analysis techniques have been used previously to identify the dominant reactions in a large set [23–26]. Briefly, a sensitivity analysis approach selects reactions based on the response of the model solution to changes in the kinetic parameters (e.g., activation energy); rate analysis selects reactions based on their contribution to the overall production rate of the important species. We have also used a dimensional analysis approach based on the Damkohler number to determine the important reactions within the porous region. The Damkohler number is the ratio of the reaction rate constant to the diffusion coefficient of the reactants. Reactions with the largest Damkohler numbers are selected because they are typically the fastest reactions in a parallel set and would be the dominant terms in a species balance equation. To reduce the large reaction mechanism to a smaller set, we established criteria for a good comparison of the
N. Birakayala, E. A. Evans / Carbon 40 (2002) 675 – 683
reduced mechanism to the large mechanism. The chosen reactions represent the general mechanism only when the deviations of the predictions of the reduced mechanism from the complete mechanism are below 10% for the dominant species. The dominant species depend on the reaction conditions chosen. We have chosen the same conditions as those used in the work of Bammidipati et al. so that we could compare our model prediction to experimental data [12]. The data were collected at pressures between 10 and 40 Torr (5.3310 3 Pa) and temperatures between 1273 and 1373 K. An isothermal zone was maintained around the substrate. Deposition rates on the substrate were measured at different pressures, temperatures, flow rates (50–400 SCCM), and diluent concentration (0–87.5%). Diluent gases used in the study included hydrogen, nitrogen, and argon. At these conditions the dominant hydrocarbon species are acetylene, benzene, and ethylene; hydrogen (H 2 ) makes up nearly 83% of the reactant mixture. The reduced reaction sets are considered
677
reasonable as long as they predict the concentrations of the three main hydrocarbon species to within 10% of the general mechanism. Sensitivity, rate, and Damkohler analyses selected the same reactions for the bulk gas phase adjacent to the growing surface. Table 1 shows the dominant reactions in decreasing order of importance. Fig. 1 shows the error of the reduced mechanism as a function of the number of reactions within the mechanism. The concentration of benzene in the reactor becomes significant only at high residence times. At lower residence times, acetylene is the dominant species contributing to the growth of carbon. At conditions where acetylene is the dominant growth species, the reactions involving benzene and its precursors are not important. Reactions 1–10 in Table 1 are required to accurately predict the C 2 species concentrations. The errors lie well below 10% and can be used to predict the concentrations of the dominant C 2 species. Further reduction from these ten reactions leads to errors much greater
Table 1 Reduced mechanism deduced from sensitivity, rate, and Damkohler analysis Arrhenius expression (mol, cm, s) Reactions to predict C1 –C2 1. H1CH 3 (1M)⇔CH 4 (1M) 2. H1CH 4 ⇔CH 3 1H 2 3. H1C 2 H 3 ⇔C 2 H 3 (1M) 4. H1C 2 H 3 ⇔H 2 1C 2 H 2 5. H1C 2 H 4 1(M)⇔C 2 H 5 (1M) 6. H1C 2 H 4 ⇔C 2 H 3 1H 2 7. H1 C 2 H 6 ⇔C 2 H 5 1H 2 8. 2CH 3 (1M)⇔C 2 H 6 (1M) 9. 2CH 3 ⇔H1C 2 H 5 10. CH 3 1C 2 H 4 ⇔C 2 H 3 1CH 4
1.27310 16 T 20.6 exp(192.7 / T) 6.6310 8 T 1.6 exp(5454.9 / T) 5.6310 12 exp(1207.7 / T) 3310 13 1.08310 12 T 0.5 exp(915.9 / T) 1.32310 6 T 2.5 exp(6159.4 / T) 1.15310 8 T 1.9 exp(3789.3 / T) 2.12310 16 T 1.0 exp(312 / T) 4.99310 12 T 20.1 exp(5334.1 / T) 2.27310 5 T 2.0 exp(4629.6 / T)
Reactions to predict C2 –C6 species 11. C 2 H1H 2 ⇔H1C 2 H 2 12. H1C 2 H 4 (1M)⇔C 2 H 2 (1M) 13. C 4 H 3 1C 2 H 3 ⇔C 6 H 6 14. C 3 H 3 1C 3 H 3 ⇔C 6 H 6 15. C 3 H 4 ⇔C 3 H 3 1H 16. CH 3 1C 2 H 2 ⇔C 3 H 4 1H 17. C 2 H1C 2 H 2 ⇔C 4 H 3 18. C 2 H 2 1C 2 H 2 ⇔C 4 H 4 19. C 2 H 3 1C 2 H 2 ⇔C 4 H 4 1H 20. C 2 H1C 2 H 4 ⇔C 4 H 4 1H 21. CH 3 1C 2 H 6 ⇔C 2 H 5 1CH 4 22. C 3 H 4 1C 3 H 3 ⇔C 6 H 6 1H 23. C 2 H 4 (1M)⇔H 2 1C 2 H 2 (1M) 24. C 4 H 4 1C 2 H 2 ⇔C 6 H 6 25. C 3 H 4 1C 2 H⇔C 3 H 3 1C 2 H 2 26. H1C 2 H 5 (1M)⇔C 2 H 6 (1M) 27. H1C 2 H 3 (1M)⇔C 2 H 4 (1M) 28. H1C 2 H 5 ⇔H 2 1C 2 H 4
4.07310 5 T 2.4 exp(100.6 / T) 1.00310 17 T 21.0 2.87310 14 exp(412.6 / T) 3.0310 11 1.00310 17 exp(35225.4 / T) 6.74310 19 T 22.1 exp(15896.7 / T) 1.00310 13 2.45310 14 exp(23349.4 / T) 2.0310 12 exp(25 16.1 / T) 1.21310 13 6.14310 6 T 1.7 exp(5258.7 / T) 2.2310 11 exp(1006.4 / T) 8.0310 12 T 0.4 exp(44670.9 / T) 4.47310 11 exp(15141.91T) 1.00310 13 5.21310 17 T 21.0 exp(795.1 / T) 6.08310 12 T 0.3 exp(140.9 / T) 2.0310 12
Reactions are listed in order of decreasing importance. Reactions with a third body (1M) have rate constants that are pressure dependent. The pressure dependency can be expressed in either the Lindemann or Troe form. Rate parameters were taken from the Gas Research Institute database and verified using the NIST database [19,20].
678
N. Birakayala, E. A. Evans / Carbon 40 (2002) 675 – 683
Fig. 1. Error in the species distribution predicted by the mechanism shown in Table 1. Error is based on comparison to the full mechanism (47 reactions) prediction.
than 10%; the predictions of the reduced models proposed by Bammidipati et al., for example, differed from our general mechanism prediction by up to 100% [12]. To predict the benzene concentrations accurately, reactions 11–28 are also required. If an error of nearly10% is acceptable, reactions 1–16 can be used to predict the concentrations of species up to and including benzene. Fig. 2 shows a reaction network based on our reduced mechanism for the pyrolysis of methane at pressures between 10 and 40 Torr (5.3310 3 Pa) and temperatures between 1273 and 1373 K; the larger arrows indicate the dominant path between methane and the products for this set of conditions. Sensitivity and rate analyses were used to validate the use of Damkohler numbers for selecting reactions in the gas phase. Sensitivity and rate analyses are based on the kinetics alone while the Damkohler number incorporates the influence of transport. To select reactions within the
Fig. 2. Network showing the reaction pathways during the pyrolysis of methane at CVD/ CVI conditions. The larger arrows indicate the dominant reaction path. The network is based on reactions shown in Table 1.
porous substrate, therefore, we used the Damkohler number as a guide for selecting the dominant reactions [14]. The dominant gas phase reaction set inside the porous region was the same as the set shown in Table 1.
5. Surface kinetics In order to predict the change in porosity of the C–C composite, surface reactions are included; parameters for surface reactions are not available at conditions for the manufacture of C–C composites. To estimate the kinetic parameters for the surface reactions we have used the deposition rate data of Bammidipati et al. in order to fit the surface reaction parameters; the reported deposition rates were measured on non-porous substrates for different experimental conditions [12]. We tested the model by comparing predicted and measured deposition rates; we did not compare gas phase compositions. For a detailed mechanism, the potential species, which react on the surface, should first be identified. To identify the possible growth species, we applied the gas phase kinetic model to simulate the gas phase conditions during the carbon deposition experiments of Bammidipati et al. [12]. These authors predict significant velocity gradients in the reactor in both the radial and axial directions. Allowing for such velocity gradients, our simulations were performed at an average temperature of 1395 K and a pressure of 40 Torr (5.3310 3 Pa). The
N. Birakayala, E. A. Evans / Carbon 40 (2002) 675 – 683 Table 2 Surface mechanism for carbon deposition during CVD/ CVI 1. C 2 H 2 12C(S)⇒2C(S)12C(D)1H 2 2. C 6 H 6 16C(S)⇒6C(S)16C(D)13H 2 3. C 2 H 4 12C(S)⇒2C(S)12C(D)12H 2 4. H1CH(S)⇔C(S)1H 2 5. H1C(S)⇔CH(S)
679
sticking probabilities for the adsorbing reactant. The sticking probabilities in our model were obtained by matching the simulation results to the measured deposition rates of Bammidipati et al. [12]. Sticking probability, which is a function of coverage (u ), is defined as Number of molecules that stick or react g (u ) 5 ]]]]]]]]]]]] Number of molecules that impinge on a surface (1)
model gives the distribution of various hydrocarbons along the length of the reactor at the given conditions of temperature, pressure, flow rate, and initial composition [22]. The species distribution was obtained at various flow rates of methane and the amount of carbon available for deposition from respective species was calculated as a function of flow rate. Based on this analysis, the major contribution for carbon deposition is from methane, acetylene, ethylene, benzene, and ethane. To be classified as a growth species, apart from the gas phase abundance, the species should also have a high reaction probability (Ks ). In recent work, Becker and Huttinger suggested acetylene, ethylene, and benzene as the growth species based on their respective reaction probabilities [5]. Based on these two analyses, we have selected C 2 H 2 , C 2 H 4 , and C 6 H 6 as the major growth species participating in the surface mechanism. The resulting surface mechanism is shown in Table 2. The rate constants for the surface reactions are given as
Fig. 3 shows the comparison of the model predictions with the experimental growth rates at 40 Torr (5.3310 3 Pa) and 1373 K. The simulation results are qualitatively consistent with the observed growth rates. The sticking probabilities used were:
gC 2 H 2 55310 25 gC 2 H 4 55310 27 gC 6 H 6 54310 24 (2) The sticking probability is a function of surface coverage. Various models have been proposed to relate the variation in the sticking probability with coverage and / or pressure. As the pressure increases the number of open active sites decreases as a result of increased coverage by any adsorbed species. With fewer open sites for adsorption the surface reaction rates decrease. The decrease in the surface reaction rates with pressure can be accommodated through the sticking probability; an alternative approach would be
Fig. 3. Comparison of the model to experimental data. P540 torr (5.3310 3 Pa), T51373 K, CH 4 : 100%. Experimental data with error bars are taken from Bammidipati et al. [12].
N. Birakayala, E. A. Evans / Carbon 40 (2002) 675 – 683
680
Table 3 Sticking coefficients used to fit the model predictions to the experimental rates Pressure (Torr)
C2H2 C2H4 C6H6
10 (1.3310 3 Pa)
40 (5.3310 3 Pa)
5310 24 4310 24 4310 23
0.5310 24 4310 27 4310 24
to use additional surface reactions to account for the blocking of active surface sites. The simple Langmuir model predicts a decrease of the sticking coefficient with increasing coverage.
g (u ) ]] 5 (1 2 u )n g (0)
(3)
Table 4 Effective binary diffusivities (EBD) for species in methane Species
EBD (m 2 / s)
H H2 CH 4 C2H2 C2H6
6.9310 22 4.5310 22 1.72310 22 1.38310 22 1.25310 22
the concentration profile of the species in the pores are fast compared to the changes in pore geometry, exist in the substrate [1]. This assumption requires gas phase species reaching equilibrium at relatively short time scales compared to the densification time. This is true in ICVI processes, where the densification times typically are 100– 300 h. From the pseudo steady-state hypothesis, the continuity equation reduces to
OR K
where g (0) is the sticking coefficient at zero coverage and n is the number of sites that is needed to hold the adsorbate. For n 5 1, Eq. (3) shows a linear decrease and for n . 1 it shows a nonlinear decrease. Bammidipati et al. [12] have reported a decrease in the sticking coefficients with increasing pressure to account for the increased coverage. For fitting the data at pressures lower than 40 torr, the sticking coefficients are increased. The sticking coefficients used in the current work are given in Table 3.
6. Pore model A cylindrical pore model has been chosen to predict the changes in porosity as deposition occurs. The pore model can be used to predict the concentration and deposition profiles within the pore. The cylindrical pore model is sufficient for the current work, since the objective of the work is to understand qualitatively the difference in using the various kinetic models. We have followed the development of Moene et al. to set-up the pore diffusion model [13]. The first assumption is that the diffusive flux is dominant in the substrate. This assumption is valid for the ICVI processes, where the thermal and flow gradients are absent. The flux term (Ni ) can be written as Ni 5 2 De,i (=Ci )
(4)
The mixture averaged diffusivity, Dei , is equal to the effective binary diffusivity of the species with respect to methane, De,i,CH 4 . The effective binary diffusivities used in this model for some representative species are given in Table 4. The effective diffusivity is calculated as the harmonic mean of binary and Knudsen diffusivities [13,14]. Pseudo steady-state conditions, at which the changes in
=? (ANi ) 5 (VTu )
OR S
Gi
1 (VT a)
k 51
si
(5)
k 51
where A is the cross-sectional area, Ni is the flux of i, VT is the total volume, u is the porosity, a is the internal surface area per unit volume, and R Gi and R Si are the net rate of decomposition of i in the gas phase and on the surface, respectively. Heat effects due to the gas phase reactions and carbon deposition are assumed to be small and a uniform temperature profile exists throughout the pore. Since very low deposition rates are maintained in the CVI reactors, the low conversion of the source gases within the pores does not affect the temperatures. For the case of cylindrical pores representing the pores in the substrate, the flow is modeled as a one-dimensional problem. The aspect ratios of the pore, defined as the ratio of the length over the diameter of the pore, vary from 50 to 250. Hence, the diffusion in the radial direction can be neglected and a uniform concentration profile over the cross-section can be assumed. Incorporating the above assumptions and combining Eqs. (4) and (5), the mass balance on a component i, within a cylindrical pore of radius r, over a small length (DX) gives ≠C S D 1 D (pr )S] D ≠x
≠Ci 2 Dei (p r 2 ) ] ≠x
2
x
ei
i
x 1Dx
5 R Gi (p r 2 Dx) 1 R Si (2p r Dx)
(6)
For an infinitesimal element, DX →0, Eq. (6) becomes, ≠(2Dei r 2 ≠Ci / ≠x) 2 ]]]]] 2 R Gi (r ) 2 R Si (2r) 5 0 ≠x
(7)
The equation is further simplified by assuming that the diffusivity and the pore radius do not change considerably in a very thin section. If the pore is divided into a fine grid with several thin sections, the diffusivity and radius remain constant within the section. However, the diffusivity and
N. Birakayala, E. A. Evans / Carbon 40 (2002) 675 – 683
the radius differ from section to section. The equation is solved at different instances of time independently of the equation describing the growth rate. With these simplifications, d 2 Ci 2 Dei r 2 ]] 2 R Gi (r 2 ) 2 R Si (2r) 5 0 dx 2
(8)
The radius of the pore as a function of the deposition rate is deduced from Eq. (9). The radius as a function of the pore length and time is determined by solving Eqs. (8) and (9) simultaneously.
O
M I ≠r ] 5 2 ]s R ≠t rs i 51 si
(9)
The boundary conditions for the problem are as follows: Ci (0,t) 5 Ci,bulk , the concentration of the species i in the bulk; Ci (L,t) 5 Ci,bulk and r(x,0) 5 r0, initial radius of the pore. The equations are solved with the boundary conditions to obtain the variation of radius as a function of the pore length and time.
7. Discussion For the comparison of different kinetic mechanisms, a pore aspect ratio of 100 has been used to distinguish between the mechanisms. The differences between the model predictions for two different kinetic mechanisms are
681
shown in Fig. 4; the densification profile is shown for the same time step in each case. The overall kinetic equation results are determined using the direct deposition of methane while the detailed kinetics combine the mechanisms given in Table 1 (reactions 1–16) and Table 2. The detailed kinetic model predicts a non-uniform deposition profile. A more uniform profile is predicted with the single reaction when methane is depositing carbon directly without the formation of any intermediates. A higher diffusivity of methane compared to that of acetylene (|19%) at the same conditions may have contributed to the uniformity. Also, because of the higher concentration of methane at the surface predicted by the detailed mechanism, the deposition within the pore occurs at a faster rate leading to pore closure at much shorter times. Thus, a completely different densification profile is obtained when the kinetics are changed, as expected. The deposition reaction results in a concentration gradient within the pore decreasing the concentrations of the growth species away from the opening (boundary) of the pore. The matrix forms rapidly at the surface leading to residual porosity. A uniform densification of the composite is possible, however, by operating at conditions of low reaction rates and high transport rates. Typically, low temperatures are employed in CVI reactors to reduce the density gradients. The results of Middleman [9] suggest that inside-out densification can be achieved by making use of the methane pyrolysis kinetics. If the growth species, which deposit carbon, were produced within the pore at a rate higher than its dissociation on the surface, a
Fig. 4. Final pore radii for different kinetic models: a single deposition reaction and the kinetics shown in Tables 1 and 2. The initial pore radius was 2.24310 26 m with an aspect ratio5100. The dimensionless pore lengths 0 and 1 are the two ends of the cylindrical pore.
682
N. Birakayala, E. A. Evans / Carbon 40 (2002) 675 – 683
Fig. 5. Initial evolution of pore structure resulting from different acetylene concentrations (Cg [5]kg-mol m 3 ) near the pore boundary. The nominal residence time within the reactor was 16 h. The initial pore radius was 2.24310 26 m with an aspect ratio550. The dimensionless pore lengths 0 and 1 are the two ends of the cylindrical pore.
uniform densification would be obtained. For example, the simulation is carried out at low acetylene concentrations at the entrance to the pore. Fig. 5 shows the deposition profiles for low acetylene concentrations at very short times. At very low or no acetylene at the pore opening, an inside-out densification is achieved. As the acetylene concentration is increased, the profile flattens and finally residual porosity sets in; the concentration gradients within the pore increase as the concentration of the growth precursors increase in the external gas phase. This result suggests that inside-out densification may be achieved when the growth precursor (acetylene) is produced within the pore. The low acetylene concentrations in the reactor are possible at high flow rates or with high diluent (H 2 , N 2 , Ar) concentrations; however, nitrogen and argon have been shown to have a negative influence on the deposition process [14,17]. Simulations with a single gas phase reaction would not result in similar predictions as there is no production of any intermediates depositing carbon and the source gas is assumed to deposit carbon directly. At atmospheric pressures and infinitely long times, Benzinger and Huttinger have demonstrated that the degree of pore filling reaches a maximum with respect to the partial pressure of the growth precursor (methane). As the partial pressure of methane is increased above 150 Torr there is a rapid decrease in the degree of pore filling [16]; the trend above 150 Torr agrees qualitatively with the results presented here. Our future work will investigate the
possibility of such a maximum at the conditions used in this work.
8. Conclusions The importance of accurately representing the kinetics during CVD/ CVI of C–C composites is clear. We have demonstrated that the dominant kinetics can be identified and used to predict deposition profiles within the pores. Furthermore, conditions for uniform or inside-out densification can be suggested based on the kinetics. A more rigorous pore model is required to estimate the fabrication times and substrate properties as a function of fabrication time. Likewise, more experimental data are required to validate the model predictions. The same approach can be used to model other composite materials.
References [1] Vaidyaraman S, Lackey WJ, Agrawal PK, Starr TL. Carbon 1996;34:1123. [2] McAllister P, Wolf EE. AIChE J 1993;39:1196. [3] Becker A, Huttinger KJ. Carbon 1998;36:177. [4] Becker A, Huttinger KJ. Carbon 1998;36:201. [5] Becker A, Huttinger KJ. Carbon 1998;36:213. [6] Becker A, Huttinger KJ. Carbon 1998;36:225.
N. Birakayala, E. A. Evans / Carbon 40 (2002) 675 – 683 [7] Bruggert M, Hu Z, Huttinger KJ. Carbon 1999;37:2021. [8] Antes J, Hu Z, Zhang W, Huttinger KJ. Carbon 1999;37:2031. [9] Middleman S. J Mater Res 1989;4:1513. [10] Sheldon W. J Mater Res 1990;5:2729. [11] Lackey WJ, Vaidyaraman S, Beckloff BN, Moss TS, Lewis JS. J Mater Res 1998;13:2251. [12] Bammidipati S, Stewart GD, Elliott Jr. JR, Gokoglu SA, Purdy MJ. AIChE J 1996;42:3123. [13] Moene R, Dekker JP, Makkee M, Schoonman J, Mouljin JA. J Electrochem Soc 1994;141:282. [14] Benzinger W, Huttinger KJ. Carbon 1996;34:1465. [15] Benzinger W, Huttinger KJ. Carbon 1998;36:1033. [16] Benzinger W, Huttinger KJ. Carbon 1999;37:181. [17] Benzinger W, Huttinger KJ. Carbon 1999;37:931. [18] Benzinger W, Huttinger KJ. Carbon 1999;37:841. [19] Gas Research Institute Database, version 2.11. http: / / www.gri.org
683
[20] Westley F, Frizzell DH, Herron JT, Hampson RF, Mallard WG. NIST Chemical Kinetics Database Version 5.0, NIST Standard Reference Data, Gaithersburg, MD 20899, USA, 1992. [21] Murphy DB, Carroll RW, Klonowski JE. Carbon 1997;35:1819. [22] Kee RJ, Rupley FM, Miller JA. The CHEMKIN Thermodynamic Database, Sandia National Laboratories document SAND87-8215B, 1990. [23] Rota R, Bonini F, Servida A, Morbidelli M, Carra S. Chem Eng Sci 1994;49:4211. [24] Brock EE, Savage PE, Barker JR. Chem Eng Sci 1998;53:857. [25] Wolden C, Gleason KK, Howard JB. Combustion Flame 1994;96:75. [26] Frenklach M. Modeling. In: Gardiner WC, editor, Combustion Chemistry, Berlin: Springer, 1984, Chapter 7.