Influence of temperature and pressure on carbon black size distribution during allothermal cracking of methane

abstract

Allothermal cracking of methane is a suitable and eco-friendly way to simultaneously produce hydrogen and carbon black. The economic viability of the process relies on the ability to produce carbon black having well-defined characteristics, particularly concerning the particle size. A model for the study of the carbon particle size distribution during thermal cracking of methane has been developed. The model takes into account: heat transfer by conduction, convection, particle and gas radiation, homogenous and heterogenous reactions of methane dissociation, nucleation, and growth of solid carbon particle formed. The model alleges nanoparticles are in thermal equilibrium and do not impact the flow. A parametric study is made on operating pressure and temperature. As a result, the increase of the pressure and temperature increases the yield of thermal methane cracking. Moreover, results show a particle size distribution becoming narrower with increasing temperature and/or pressure. In these conditions, the particles population tends to be monodispersed.

Copyright © 2016 American Association for Aerosol Research

EDITOR:

• Mark T. Swihart

1. Introduction

Total world production of carbon black is about 10 million tons per year and the demand is constantly growing. Ninety percent of this production is nowadays produced by the furnace process, Lockwood and Vanniekerk (1995). The furnace process is based on an incomplete combustion of hydrocarbons and is among the most polluting industrial process since a significant part of the feedstock is burnt for the only reason to reach the necessary high thermal conditions for carbon black production (usually in the range of 1600 to 2200 K) (Fabry et al. 2001; Fulcheri et al. 2002).

Fulcheri and Schwob (1995) suggested an eco-friendly concept for carbon black production using plasma methane cracking. Following the latter publication and during about 20 years, researches on carbon nanomaterial plasma processing have been carried out at MINES Paris Tech: Fabry et al. (2001), Fulcheri et al. (2002), Fulcheri and Schwob (1995), and Moreno-Couranjou et al. (2009). The plasma process for methane cracking consists in replacing the flame process usually made by incomplete combustion (furnace process), by a thermal plasma generated by an electric discharge. The injection of hydrocarbons in such a high temperature medium induces thermal cracking reactions leading to the production of carbon black particles and hydrogen. Researches in this domain yielded to the development of an original three-phase AC plasma technology (Fabry et al. 2001; Fulcheri et al. 2002; Moreno-Couranjou et al. 2009). In parallel, an industrial attempt of plasma pathway for carbon black production using a DC technology has been performed during the 1990s by the KVAERNER company (Gaudernack and Lynum 1996). High temperature oxygen-free allothermal cracking of methane is a promising route for a large-scale co-production of hydrogen and carbon black with little emissions in contrary to the usual furnace process (Fulcheri and Schwob 1995; Steinfeld 2005; Maag et al. 2009; Rodat et al. 2010). Nevertheless, the economic viability of this technique critically relies on obtaining the desired properties of the carbon black produced (Abanades and Flamant 2007), and therefore on the control of the formation and growth of carbon particles.

Computational fluid dynamic (CFD) simulations have already been performed in order to have better insights of the thermal methane cracking and the particle formation for different reactor configurations. Kogan et al. (2007) performed a two-dimensional (2D) axis-symmetric simulation of the confined tornado flow reactor. Kogan et al. (2007) were thereby able to identify certain basic characteristics of the confined tornado flow that were observed in the laboratory, such as the slight deposition or not of solid particles on the reactor windows according to the carrier gas’ viscosity. Deme (2002) and Gonzalez-Aguilar et al. (2003) made a three-dimensional (3D) modeling of methane dissociation in a tubular reactor. Gonzalez and Deme showed how the radiative heat transfer varies in the presence of carbon particles. The theory selected to model the absorption coefficient induced by carbon particles was the Mie theory associated with the Rayleigh limit. Similarly Ravary et al. (2003) also used this assumption to simulate a plasma reactor for the production of nano-sized carbon materials. Patrianakos et al. (2011), by the use of an only one-dimensional simulation, added characteristics of a polydisperse particle population. As a result, they were able to model size dependent heterogenous reactions on the surface of the particle and also on the wall of the reactor. Patrianakos et al. were also able to dispose of the Rayleigh limit assumption to model the absorption induced by carbon particle by knowing the particle size distribution (PSD) and directly using the Mie theory. The PSD in this model was calculated while solving the population balance equation (PBE), Friedlander (2000), using a sectional method, so called the class method first used by Gelbard et al. (1980). The PBE describes the dynamic of particle formation and growth. The PBE is an integro-differential equation depending on spatial position and on an internal variable, which is usually the particle size (Ramkrishna 2000). It is computationally expensive to solve such an equation. Patrianakos’s model has been validated with experimental data of the degree of methane conversion. By using the same model, Patrianakos et al. (2012) have been able to emphasize the effect of seeding on the carbon production. Lastly, a two-dimensional model of methane cracking presented by Caliot et al. (2012) has been used for the simulation of a tubular reactor (Rodat et al. 2010). The two-dimensional axis-symmetric simulation of Caliot et al. describes fluid flow, conduction, convection, radiation heat transfer, gas-phase kinetics, particle formation, growth, and coagulation. Just like Patrianakos, Caliot used a class method to solve the population balance equation of carbon particles. Caliot et al. (2012) also added a detailed multi-gray radiative model with a particle size population dependency over the absorption coefficient. They observed the influence of radiation absorption by solid particles and how the reaction zone is consequently shifted.

Improvements in theses latter models can still be done; especially in the way carbon particles initially form in the gas phase, i.e., nucleation (Lahaye 1992; Berezkin 2001; Frenklach 2002; Fincke et al. 2002; Richter et al. 2005; Rodat et al. 2009). But going any further in the complexity of the model would be fruitless without having a better understanding on the impact of variables such as thermal boundary conditions and operating pressure within the reactor. Indeed, all of those previous simulations used different pressure values and they all agreed on the fact that temperature plays a major role on results without explicitly showing the role of varying operating conditions on the obtained particle size distribution.

This study models the methane cracking in order to see the influence of temperature and pressure on the carbon particle formation and growth. Two-dimensional axis-symmetric laminar fluid flow simulations with conduction, convection, and radiation heat transfer have been performed. The PBE is solved using the class method taking into account nucleation, growth, and coagulation of solid carbon nanoparticles. After having compare model results with results published in the literature for a tubular geometry (Caliot et al. 2012), a parametric study is presented varying the temperature and the operating pressure. The particle size distribution is chosen as the main comparative output result in order to obtain the direct characterization of solid carbon particles produced according to the varied parameters. Simulation results are interpreted considering the model assumptions. This study constitutes a preliminary step in the aim to model methane cracking by plasma process while controlling particle size distribution. Further model improvements are discussed in the conclusion.

2. Formulation of the model

2.1. Governing equations

The reaction of methane cracking is modeled using the single overall kinetic reaction given by Equation (1[1] ):[1]

This study is limited to the modeling of laminar flow where methane cracking happens. Assuming the very small volume fraction of the solid phase produced by Equation (1[1] ), the flow is assessed not to be disturbed by the nanoparticle population and embedded solid particles have hypothetically the same velocity that the gas phase. In the same way, the small size of these solid particles leads to a high specific surface area (several tens ), which indicates that thermal inertia of these particles can be neglected. In consequence, the temperature of the particle is equal to the gas phase temperature. Given the previous assumptions, the biphasic flow is treated as a single-phase fluid flow. The system of local steady state governing equations describing the flow is written as[2] [3] [4] [5]

with :[6] [7] [8] [9]

where is the total energy, is the total enthalpy, is the enthalpy of the species , and is the diffusion coefficient of the gas species . Chemical source terms in the energy and species mass conservation equations take into account homogenous and heterogenous reactions of methane cracking as follow:[10] [11] [12]

Parameters of the total volumetric heterogenous reaction rate, denoted , are described in Section 2.2.2. The heterogenous reaction represents the methane conversion on the surface of solid carbon particles while the homogenous reaction refers to the direct production of carbon nanoparticle by methane dissociation. Although carbon particles are not treated as an independent phase, their contribution in the density, global enthalpy, and conductivity in the only phase has to be established and will depend of the considered case. At this state, the phase can be interpreted as a gas-particle mixture. As mentioned previously in this section, the flow is assumed not to be disturbed by the particles and therefore the presence of particles does not influence the single-phase fluid’s viscosity.

2.2. Population balance equation

The theory of polydispersed particle dynamics refers to the study of entities such as solid particles, droplets or bubbles, characterized by one or more properties. The common given name of this research field is Aerosol Dynamics. Aerosol Dynamics has been studied since the beginning of the past century: Fuchs (1964), Friedlander (1982), Warren and Seinfeld (1984), Phanse and Pratsinis (1989), Pratsinis (1988), Pratsinis (1993), Friedlander (2000), Seinfeld et al. (2003). A Population Balance Equation (PBE) (Ramkrishna 2000) describes particle dynamics under the influence of various physical and chemical phenomena, i.e., convection, diffusion, nucleation, surface growth, coagulation. A well-known formulation is given by Friedlander (2000), known as the General Dynamic Equation for Aerosol, and can be seen in Equation (13[13] ).[13]

where represents the number density function, is the particle volume, is the critical particle volume, is the gas velocity, is the particle diffusion coefficient. is the nucleation rate, is the Dirac delta function, is the heterogenous condensation rate, and is the collision frequency function for coagulation. The class method, also called discrete method, is a sectional method that was developed by Gelbard et al. (1980) in order to return a numerical solution of the PBE. It is based on representing the continuous Particle Size Distribution (PSD) in terms of a set of discrete size classes or bins. By discretizing Equation (13[13] ) and considering the steady state, the PBE, Equation (13[13] ) becomes Equation (14[14] ).[14]

The advantages of this method are numerical robustness and the ability to return the PSD directly without any profile assumption. The drawbacks are that bins must be defined a priori and that a large number of bins may be required. The discrete method is therefore computationally expensive but requires, in theory, less assumptions. In this study, the particle size distribution function is discretized with a geometric progression on the volume. Equation (15[15] ) shows this discretization. M different particle sizes are considered:[15]

The volume fraction and the mass fraction of particle in the bin , i.e., size , can be written as follow:[16] [17] [18]

where represents the concentration of the particle size [], and stands for the number density function.

The next subsections provide the expressions and details of the terms involved in Equation (14[14] ).

2.2.1. Diffusion coefficient

[19] [20] [21] [22] [23]

The formulation for the diffusion coefficient is known as the Stokes–Einstein equation. The latter equation is corrected by a Slip Correction Factor () accounting for non-continuum effects of the alleged particle phase. Scalar values for the are empirical and specific to carbon particles (Caliot et al. 2012; Colombo et al. 2012).

2.2.2. Heterogenous reaction rate

[24] [25] [26] [27] [28] [29] [30]

The above heterogenous reaction formulation is based on a two-step mechanism: diffusion transport of methane to the surface of particles and first-order heterogenous reaction on their surface. Thereby, the mass transfer diffusion coefficient, denoted , and the kinetic reaction rate, denoted , influence the growth rate and can be seen in Equation (24[24] ). The Sherwood number, noted by , is considered equal to 2 for a stationary sphere in the continuum regime. The above heterogenous reaction formulation is described in greater detail by Patrianakos et al. (2011).

2.2.3. Coagulation source term

The coagulation source term in the PBE is expressed by[31] [32]

where the coagulation kernel is[33] [34] [35] [36] [37] [38] [39]

The above formulation of the coagulation kernel, denoted , is known as the Fuchs formulation and is described by Fuchs (1964). This is an extension of the Smoluchowski theory in order to account for the transition from the free molecule regime to the continuum range. As a result, it can be applied to particle sizes smaller as well as larger than the mean free path of the carrier gas. It is a general formulation for the coagulation kernel and can be applied to the complete range of Knudsen particle numbers. This theory assumes neutral particles and instantaneous fusion after all collisions.

2.2.4. Nucleation source term

The nucleation source term is presented in Equation (40[40] ) and is related to the single overall kinetic reaction assumption. The nucleation source term formulation assumes that every methane dissociation event produces one carbon nanoparticle of the smallest size considered (i.e., first bin). Equation (40[40] ) is based on an Arrhenius law and can have different sets of kinetic parameters (Gonzalez-Aguilar et al. 2004; Trommer et al. 2004; Patrianakos et al. 2011). This formulation weakly describes the initial carbon formation process but has the advantage to be simple and close to experimental observations. The kinetic monomer production is in reality much more complex than a single one-step chemical reaction (Frenklach 2002; Fincke et al. 2002; Richter et al. 2005). In addition, nucleation theory usually involves clusters formation and determines a critical volume of such structure using a supersaturation ratio and operating conditions (Doring and Becker 1935; Katz and Wiedersich 1971; Girshick and Chiu 1989). The present model constitutes only a preliminary step to the modeling of thermal methane cracking. A more relevant and precise nucleation model involving polycyclic aromatic hydrocarbons formation will be developed in future investigations.[40]

2.3. Radiation transport model

Methane cracking processes for carbon black and hydrogen production work at a temperature range lower than 2500 K (Fulcheri et al. 2002; Rodat et al. 2011). Consequently, radiations are mainly in the infrared spectrum (more than 90% of the black body emitted energy at 2200 K has wavelengths larger than 1 , 84% at 2500 K). Experimental results of Fulcheri et al. (2002) show that carbon particle sizes are less than 1 . Therefore, the particle size parameter is less than one. In this case, light scattering is then negligible in comparison with light absorption for carbon black particles, Gonzalez-Aguilar et al. (2004). The gas-particle mixture is here considered as a nongray gas that does not scatter the radiation.

The radiative heat transfer equation is solved by using the “spectral-line-weighted sum-of-gray-gases” (SLWSGG) in a nongas medium described by Modest (2013). The Radiative Transport Equation (RTE) for the th gray gas is expressed as[41]

where and are, respectively, the emission weighting factor and absorption coefficient of the mixture for the th gray component. The absorption coefficient for carbon particle, denoted , is calculated using the Mie theory. The subroutine BHMIE, available in the appendix of Bohren and Huffman (1983), is computed in order to return the wavelength absorption efficiency for each particle size, denoted . Figure 1 shows variation over the temperature for different particle size. The global absorption coefficient is then computed in the mixture property using the expression:[42]

with :[43]

Figure 1. Wavelength average absorption efficiency for different particle diameters.

2.4. Thermodynamic properties and transport coefficients

Material properties in Equations (2[2] ) to (14[14] ) are determined using the software T&Twinner (2009). These calculations are performed for the specific heat, the thermal conductivity and the viscosity. Eventually the temperature dependence of all these properties is implemented using polynomial functions in the CFD software. Default FLUENT data of solid carbon are kept for the thermal properties of solid particles.

The ideal gas law is used to determine the density of the gas-carbon particle mixture. The laminar viscosity is calculated using the semi-empirical rule of Wilke Bird et al. (2007) and considering gas species only. The effective thermal conductivity, denoted in the governing Equation (4[4] ), is calculated with a parallel model between gas species conductivities and carbon particle conductivity. This parallel model can be seen in Equation (44[44] ). The gas species conductivity is calculated using kinetic theory, Bird et al. (2007), while carbon particle conductivity is equal to the carbon pure conductivity. The mass diffusivity of gas species is calculated using the Chapman–Enskog procedure (Bird et al. 2007). The diffusive flux of carbon particle in the energy governing equation, Equation (4[4] ), is neglected:[44]

3. Simulation

3.1. Case study

Figure 2. Reactor scheme of the case study with wall temperature boundary condition.

In order to validate the model, a case study from Caliot et al. (2012) is computed as a benchmark. The reactor, shown in Figure 2, is a cylinder made of graphite with a length L of 0.6 m and a diameter D of 15 mm. The temperature profile of the wall is imposed. The flow is assumed axisymmetric.

3.2. Boundary conditions

The gas inlet is composed of 50 of argon and 50 of methane in mole fraction. Argon is the carrier gas and remains the highest concentrated species within the entire reactor. The feed gas rates at standard conditions (273.15 K and 101,325 Pa) are 3.727 slpm of methane and 3.727 slpm of argon. The wall temperature profile is separated into the two following zones;[45] [46]

In the first part, the temperature increases linearly up to . In the second part the temperature stays at . For this reference case, is equal to 1800 K and the operating pressure is set to 41 kPa. The wall emissivity is set to corresponding to the emissivity of graphite at high temperature. The inlet and outlet emissivity are null. Whatever the particle size, the particle concentration at the inlet and along the wall is fixed to zero whereas at the outlet only the gradient is set to zero. Kinetic parameters are , for the homogenous reaction, Caliot et al. (2012), and , for the heterogenous reaction of methane dissociation, Caliot et al. (2012). The minimal considered volume of particle is determined using Equation (47[47] ) in order to be close to the carbon atom size. This minimal volume, denoted , is equal to 8.79E-30 . The geometric ratio for the volume discretization, , is set to 2 in order to reach the biggest size of carbon particle with a computationally acceptable number of bins. By setting up M equal to 46, the maximum diameter is then 8.39 and the expected particle size range is adequately covered (Caliot et al. 2012). Among gas species, only methane is assumed to participate in the absorption of radiations. The emission weighting factor and absorption coefficient of methane are computed from the absorption distribution function model (ADF) by Caliot et al. (2010), which is derived from the narrow band radiative properties of methane at high temperature provided by Perrin and Soufiani Perrin and Soufiani (2007).[47]

3.3. Numerical methods

The 14.5 version of ANSYS FLUENT commercial finite volume software has been used for the implementation of the model. Equations (2[2] ) to (5[5] ) are solved using the FLUENT energy and species transport models. A Coupled pressure-based solver is used, using coupled pressure-velocity scheme. Gradients are determined using the least-square cell-based method. The standard scheme is applied for the pressure interpolation. First-order upwind spatial discretization is used for all the governing equations, except for the energy equation, since the use of higher level schemes does not change simulation results. The RTE, Equation (41[41] ), is solved using the Discrete Ordinates (DO) FLUENT radiative model. The angular discretization takes place in only four octants of the space (due to the 2D axis-symmetric geometry). A 4 4 discretization per octant is made, which results in 64 different directions. User Defined Scalar (UDS) equations have been set up in FLUENT for every different particle size considered. UDS equations describe the sectional method of the Friedlander equation. The Population Balance add-on model available in FLUENT (from version 6.3 to date) could not have been used in this modeling due to the too small volume fraction of carbon particles. This critical limitation has been found during investigations of this present study and could not have been overcome due to the hard coded structure of FLUENT. User Defined Functions (UDFs) have been coded in order to compute the thermal conductivity, the absorption coefficient, the extinction coefficient, and also the mass diffusivity of the solid particles and gas mixture. UDS values are used and coupled with the single phase, by using also UDFs, in order to take into account for the chemical reactions, the particle absorption and the particle mass diffusion. Simulations have been run on a Linux high-performance computing cluster, architecture x8664. Each simulation uses four processors at 2.0 GHz. The computation time depends on the parametric case considered and is about 10–30 h for the coarsest mesh used, and more than 200 h for the finest one.

3.4. Computational mesh

The geometry in Figure 2 is oriented such that the x-direction lies along the length of the reactor and represents also the axis of symmetry. The y-direction lies along the radius of the reactor. Three different meshes have been used in this study. The first mesh is denoted as follow: 600X8, which means 600 cells along the x-direction and 8 along the y. Using the same nomenclature, the second mesh is denoted by 1200X15 and the third one by 1200X120. Both, second and third meshes, have a gradually refinement close to the wall of the reactor, i.e., the cell length along the y-direction are smaller when cells are closer to the wall. The smallest thickness along y, for the finest mesh, have a length of 3.24E-5 . In the x-direction, constant grid spacing has been used for all different meshes.

Figure 3. Temperature profiles for the case study reactor with an operating pressure of 41 kPa and T_w,max equals to 1800 K.

The coarse mesh (600X8) is used for the parametric study after being validated by comparing characteristic results between this mesh and a finer one 1200X15. Figures 3a and b show agreement between the temperature profiles for the two different meshes. Likewise, Figures 4a and b convey very similar results of the particle size distribution at the end of the reactor, which is the key result of the present study. The 1200X120 mesh is used in order to check the results of the present model with those reported by Caliot et al. (2012) at a very close position to the wall, particularly the absorption coefficient. This small difference between the two coarser meshes and the finest one close to the wall has no impact on the concerned output results of this parametric study. Additionally, the absorption coefficient value at the end of the reactor for the gas mixture is mainly dependent on the particle population. The absorption coefficient for methane is less than 4 whereas, at the end of the reactor, the absorption coefficient for the solid particle population reaches 50 . Similar absorption coefficient profiles can be seen in Caliot’s case. Figure 5 presents a comparison between the temperature profiles from Figure 3a and the temperature profiles from the case study reference Caliot et al. (2012). It can be seen that the temperature gradient is slightly less important from this study than from the reference Caliot et al. (2012). This discrepancy can be explained by the use of different thermodynamic properties and transport coefficients database sources; HSC software in (Caliot et al. 2012) and TTwinner in this study. Nevertheless a good agreement is found and confirms the relative accuracy of the present model.

Figure 4. Particle size distribution at the end of the case study reactor and for different radial positions. The operating pressure is 41 kPa and equals to 1800 K.

Figure 5. Comparison of temperature profiles with the temperature profiles of the reference (Caliot et al. 2012) denoted by the notation [R].

4. Results and discussion

4.1. Nucleation and growth study

Results of this subsection refer to the base case where is equal to 1800 K and the operating pressure is set to 41 kPa. Figures 6a–f present the mass concentration of particle for different specific particle sizes. As a result, it can be observed the particle formation and the growth along the reactor from the inlet to the outlet. The nucleation of particles can be seen in Figure 6a. Smallest particles are mostly generated close to the wall at around 0.2 m from the inlet because this is the first place where the temperature becomes high enough to decompose methane. Nucleation of particles happens at the center of the reactor as well but at a later distance downstream and in less proportion due to the lower temperature and the higher gas velocity at the center. The growth of medium size particles, in the particle range considered (Figures 6b and c) happens in the middle of the reactor. This higher concentration of medium size particles at the center of the reactor can be explained by the diffusion of small particles toward the center. The diffusion process tries to homogenize the small particle concentration. Small particles move toward the symmetric axis of the reactor while growing. Once they arrive at the center, they fuse with the other new nucleus, which just have been created due to the shifted nucleation happening at the center. This coagulation process taking place mostly between small and big particle is the consequence of the Fuch formulation of the coagulation kernel, Fuchs (1964). In fact, the highest probability of coagulation is found between particles having the highest difference in size. In other words, small particles have significantly greater chance to collide and merge with the near biggest particles, resulting in bigger particles. Eventually, the biggest particles of the reactor can be found close to the wall (Figures 6e and f) they represent the past part of small particles that stayed close to the wall during the diffusion process and because of the lower velocity at this position, they had much more time to grow and combine with each other.

Figure 6. Mass concentrations for different particle sizes [kg.m] inside the case study reactor, with an operating pressure of 41 kPa and using equals 1800 K.

Figure 7 shows the particle size distribution at different axial positions along the reactor. The different curves per graphic stand for different radial positions. The radial position is measured from the center to the wall. Figure 7 illustrates the formation and growth of particle from the beginning to the end of the reactor and how it occurs earlier close to the wall of the reactor. At the very end of the process, Figure 7f, a local maximum slightly appears. The particle size seems to tend to a narrower particle size distribution center around 100–300 nm of diameter depending on the axial position.

Figure 7. Particle size distributions at different axial positions of the case study reactor, with an operating pressure of 41 kPa and using equals 1800 K.

4.2. Parametric study results

In this section, simulations have been performed considering the same reference case study but with different operating pressures and maximal wall temperatures. The operating pressure, , within the reactor is the first parameter of this parametric study. In order to keep the same normal feed condition, the inlet velocity value is changed for each different operating pressure used. For example, the inlet velocity at absolute pressure equal to 41,000 Pa is 1.724 ; whereas for an absolute pressure of 101,325 Pa the inlet velocity becomes 0.7030 . is the second variable of the parametric study. Thermodynamic properties and transport coefficients have been recalculated using the same database for every case. This database is detailed in Section 2.4. Table 1 gives the total conversion rate of the process at the end of the reactor for different thermal and pressure conditions. This rate has been calculated according to two different ways, either by determining the methane dissociation rate (ratio between methane outlet and inlet mass flux) and by calculating the carbon particle formation rate (ratio between mass carbon flux from solid particle at the outlet and mass carbon flux from methane at the inlet). As expected, those two different methods give dissociation and conversion rate in agreement with each other, in other words, the carbon mass is conserved in the process. Results of the parametric study are displayed in Table 1. Increasing the temperature or/and the pressure leads to greater methane decomposition and, considering the model assumptions, greater carbon particle formation. Table 1 also gives two different mean particle diameters calculated at the outlet of the reactor. The first one, denoted is the length mean diameter. The length mean diameter increases when temperature and pressure go up. The second mean diameter is the Sauter’s diameter. The Sauter’s diameter is defined as the diameter of a sphere having the same volume/surface ratio as the particle, it is given by Equation (48[48] ), which is the ratio between third- and second-order moments. The Sauter’s diameter can be imagined as the diameter of particle powder with the same mean specific surface as the particle population of interest. The Sauter’s diameter is widely used in the field of reacting gas-solid flows, since the specific surface area is of most interest in reaction rate or flow interaction calculations. In this study, it has been supposed that all particles are spherical and so the specific surface only depends on the particle size distribution. Table 1 shows that the Sauter’s diameter is approximately constant with the temperature but increases with the pressure:[48]

Table 1. Characteristic values of carbon particle formation during allothermal cracking of methane inside the case study reactor.

[K]	Pressure [kPa]	Conversion rate	[]	[]
1600
1600	101
1600	161
1800
1800	101
1800	161	100
1900
1900	101	100
1900	161	100
2000		100
2000	101	100
2000	161	100

Figure 8 presents the particle size distribution at the end of the reactor with different operating pressures and different maximum wall temperatures. The figure also shows the particle size distribution for different radial positions. The first significant observation is that the temperature has more influence on the particle size distribution than the pressure in this range of temperature and pressure. For high temperature and pressure conditions, the coagulation rate is more important compare to the nucleation rate in a way there is no existing monomer at the end of the reactor. Every new emerging solid particle instantly fuses by coagulation with an already existing particle. Consequence of the coagulation process, an interesting observation is the narrowing of the particle range at the end of the reactor for the highest temperature and pressure conditions. It can be expected that for even higher temperature and pressure, the size population of carbon particle tends to be even more narrow, which would be advantageous in the desire of high value carbon black production.

Figure 8. Influence of the maximum wall temperature boundary condition and operating pressure on the particle size distribution at the outlet of the case study reactor.

5. Conclusion

A model for nucleation and growth of carbon nanoparticles in a gas phase has been developed and successfully coupled to a CFD software. The model considers energy and mass conservation, species transport, homogenous chemical reaction for monomer production, heterogenous chemical reaction for particle growth, particle coagulation, radiative heat transfer. This study brings numerical results of the formation and growth of carbon nanoparticles inside the thermal methane cracking reactor under different temperature and pressure conditions. Subsequently, the pressure and more importantly the temperature are observed to be significant parameters in the methane thermal cracking process. As expected, the higher the walls temperature the higher the methane conversion. Increasing pressure also increases the cracking yield. More importantly, the present study gives insights of the effect of varying operating pressure and temperature on the particle size distribution. Particularly, in this case, an increase of pressure and/or temperature provide a narrower particle size population what is definitively a key point regarding production of high value carbon black.

The objective to develop a tool, which gives consistent results in comparison with the literature for methane cracking simulation, has been achieved. In addition, important insights along with temperature and pressure variations have been found. However, the results have to be dealt with caution due to the model limitations. The simplicity of the chemical kinetic used leads to a possible bias. First, the simple one-step reaction expression, used to predict the monomer production, cannot describe truthfully the complexity of the initial particle formation. Although the kinetic parameters come from experimental data fit, literature agrees on the involvement of a much more complex chemical scheme (Frenklach 2002; Fincke et al. 2002; Richter et al. 2005). Second, other sets of different kinetic parameters could have been chosen (Gonzalez-Aguilar et al. 2003; Trommer et al. 2004), and this choice have a great influence on the output results. Knowing these limitations, the present model constitutes only a preliminary step toward the modeling of thermal methane cracking. Improvements are expecting for the future, particularly concerning the nucleation part. The next step will be to include a more advanced chemical kinetic model involving polycyclic aromatic hydrocarbons formation (Fincke et al. 2002); coupled with physical nucleation model (Girshick and Chiu 1989). Other interesting improvements in the presented model are looked forward such as heterogenous reactions on the reactor wall, Patrianakos et al. (2011), and plasma conditions including: (i) methane–hydrogen mixing under highly turbulent conditions, (ii) much higher reacting temperatures, (iii) much smaller residence times. The authors are currently working on this next step in order to properly apply this model to a plasma torch geometry developed in PERSEE laboratory Mines ParisTech. Main difficulties will rely on maintaining an affordable computational cost while describing all the complexity of the process. A validation will eventually follow by comparing experimental data obtained from the PERSEE plasma torch with simulation results.

Nomenclature

Latin symbols

	=	th gray body emission weighting factor []
	=	heat capacity species []
DCH4	=	diffusion coefficient for methane []
	=	diameter particle size []
	=	diffusion coefficient gas species []
	=	binari diffusivity for system []
	=	diffusion coefficient particle []
	=	diffusion coefficient particle size []
	=	total energy []
EA,het	=	heterogenous reaction activation energy []
EA,hom	=	activation energy homogenous reaction []
	=	diameter of Fuchs’ absorbing sphere []
	=	heterogenous condensation rate []
	=	mass transfer diffusion coefficient methane []
	=	growth rate by condensation for particle []
	=	enthalpy of formation species []
	=	nucleation rate in PBE []
	=	radiance []
	=	mass diffusif flux species []
	=	boltzmann constant []
	=	carbon thermal conductivity []
keff	=	effective thermal conductivity []
kgas	=	gas phase thermal conductivity []
khet,0	=	heterogenous reaction kinetic parameter []
khom,0	=	kinetic parameter homogenous reaction []
	=	knudsen number for particle size []
	=	number of bins []
	=	molecular weight carbon []
MCH4	=	methane molecular weight []
	=	molecular weight species []
	=	number density function []
	=	avogadro number []
	=	average number density function []
nspec	=	number of species []
	=	pressure []
Qabs,i	=	wavelength average absorption efficiency []
	=	ideal gas constant []
	=	geometric ratio for discretization []
	=	radiation unit direction vector []
Schem	=	chemical reactions energy source term []
	=	surface of particle size []
Si,coag	=	mass source term particle size []
Si,chem	=	chemical reactions mass source term []
Snuc	=	nucleation source term []
Srad	=	radiation energy source term []
SCFi	=	slip correction factor for particle size []
	=	sherwood number for spherical particle []
	=	mixture temperature []
	=	time []
Twall	=	temperature at the wall []
Tw,max	=	maximum temperature x at the wall []
	=	mixture velocity []
	=	mean thermal velocity of particle size i []
	=	volume particle []
	=	volume particle size []
	=	molar fraction species []
	=	mass fraction gas species []

Greek symbols

	=	volume fraction particle size []
	=	collision frequency function of coagulation []
het	=	enthalpy of heterogenous reaction []
	=	dirac delta function []
	=	absorption coefficient due to gas species []
	=	absorption coefficient for the th gray component []
	=	absorption coefficient for carbon particles []
λCH4	=	methane mean free path []
	=	mean free path for all species in gas phase []
	=	gas phase molecular viscosity []
	=	stoichiometric coefficient species []
	=	carbon density []
	=	gas phase density []
	=	Stefan–Boltzmann constant []