Modelling the flame synthesis of silica nanoparticles from tetraethoxysilane

Chemical Engineering Science 70 (2012) 54–66

Contents lists available at ScienceDirect

Chemical Engineering Science journal homepage: www.elsevier.com/locate/ces

Modelling the flame synthesis of silica nanoparticles from tetraethoxysilane Shraddha Shekar, Markus Sander, Rebecca C. Riehl, Alastair J. Smith, Andreas Braumann, Markus Kraft  Department of Chemical Engineering and Biotechnology, University of Cambridge, New Museums Site, Pembroke Street Cambridge, CB2 3RA, UK

a r t i c l e i n f o

a b s t r a c t

Article history: Received 21 December 2010 Received in revised form 21 March 2011 Accepted 8 June 2011 Available online 30 June 2011

This work proposes a kinetic model and an inception pathway for the flame synthesis of silica nanoparticles from tetraethoxysilane (TEOS). The kinetic model for the decomposition of TEOS is developed by generating reactions involving species that were reported in high concentrations at equilibrium. Flux and sensitivity analyses are then performed to identify the main reaction pathways. The parameters for these reactions are systematically fitted to experimental data using low discrepancy (LD) sequences and response surfaces. The main product of TEOS decomposition is deduced to be silicic acid (Si(OH)4). To increase computational efficiency, the kinetic model has been reduced by determining the level of importance (LOI) of each species and retaining only the important ones. This reduced kinetic model is then coupled to a detailed population balance model using an operator splitting technique. New particle inception and surface growth steps have been incorporated into the particle model in which particles form and grow by the interaction of Si(OH)4 monomers. Coagulation and sintering of particles are also included in the model and the material dependent sintering parameters have been determined by fitting the model to experimental values of collision and primary particle diameters using LD sequences. The particle size distributions and computer-generated TEMstyle images have been generated and good agreement with experiments is observed. The gas-phase reactor composition and the temporal evolution of particle size at different temperatures are also presented. & 2011 Elsevier Ltd. All rights reserved.

Keywords: Silica nanoparticle Flame synthesis TEOS decomposition Kinetics Population balance Parameter estimation

1. Introduction Silica nanoparticles find widespread applications in many fields such as catalysis, chromatography, ceramics, and as stabilising and binding agents (Iler, 1979). Recent studies have also reported their use in biotechnology for bio-sensing and drug delivery (Slowing et al., 2007). It has been demonstrated that the rates of gas phase reactions which lead to particle formation are important in determining the final product properties (Pratsinis, 1988), and thus a detailed understanding of these processes is essential. The formation of silica nanoparticles from tetraethoxysilane (TEOS) using thermal decomposition methods such as flame spray pyrolysis (FSP) is one of the most popular routes for the synthesis of these particles (Stark and Pratsinis, 2002). Despite its industrial importance, the current knowledge of gas-phase decomposition of TEOS to eventually form silica particles remains incomplete. A comprehensive understanding of this process requires the successful coupling of two detailed models: one to determine

 Corresponding author. Tel.: þ 44 1223 762784; fax: þ44 1223 334796.

E-mail address: [email protected] (M. Kraft). URL: http://como.cheng.cam.ac.uk (M. Kraft). 0009-2509/$ - see front matter & 2011 Elsevier Ltd. All rights reserved. doi:10.1016/j.ces.2011.06.010

the kinetics of the gas-phase and the other to determine the particle dynamics. Previous kinetic studies conducted in this field suggested that the main reaction pathway for the disintegration of TEOS leads to the formation of silica molecules in the gas-phase. Chu et al. (1991) studied the low temperature (700–820 K) decomposition of TEOS and suggested a six-centre decomposition, which formed diethoxy silicate ððC2 H5 OÞ2 SiQOÞ and equal amounts of ethanol and ethylene. They carried out experiments in the presence of toluene (a chain inhibitor), and these results indicated the possible contributions from a radical-induced mechanism. Takeuchi et al. (1994) studied this system in a low pressure chemical vapour deposition reactor at 950 K and 13 Pa and detected the presence of hexaethoxysiloxane ((C2H5O)3Si–O– Si(OC2H5)3). They inferred that TEOS reacted with silanol to form hexaethoxysiloxane and ethanol. Herzler et al. (1997) studied the decomposition of a few alkoxy silanes (TEOS, tetramethoxysilane (TMOS) and tetrapropoxysilane) using a heated single pulse shock-tube and analysed the reactants and final stable products using gas chromatography. The main products were shown to be ethanol and ethylene. Based on this experimental data, Herzler et al. (1997) proposed a reaction pathway involving the decomposition of TEOS to form silicates and silyl acids. This model predicted the activation barrier for the formation of ethanol and

S. Shekar et al. / Chemical Engineering Science 70 (2012) 54–66

silicate from silanol to be less than 200 kJ. However, this is in contrast to the thermochemistry calculated by Ho and Melius (1995) using BAC-MP4 methods, where high activation barriers for the formation of ethanol and silicate from silanol were predicted. Their results were consistent with the thermochemistry calculations reported by Phadungsukanan et al. (2009) for the TEOS system, where they found the equilibrium concentration of gas-phase silica at atmospheric pressure to be insignificant. The study of particle dynamics requires an understanding of the following processes: inception, coagulation, surface growth and sintering (Morgan et al., 2005). Jang (1999), Smolı´k and Moravec (1995) and Seto et al. (1997) have previously investigated the formation and sintering of silica nanoparticles from TEOS in the vapour phase. Jang (1999) synthesised and characterised silica nanoparticles from TEOS vapour in a lab-scale diffusion flame reactor. Seto et al. (1997) studied the formation of silica particles in a heated pipe and determined the particle structure and size distribution. Agglomerates were generated at 900 1C and the change in shape of the particles was evaluated by heating the particles to different final temperatures. Due to the lack of a detailed kinetic model in previous studies, the flame synthesis of silica nanoparticles from TEOS was assumed to proceed via the sticking of gaseous SiO2 monomers. The aim of this work is twofold: (i) to develop a new kinetic model which addresses the thermodynamic inconsistencies of the previous models; (ii) to develop a new particle model which incorporates new inception and surface growth steps. The paper is structured as follows. The kinetic model development for TEOS decomposition is first described in Section 2. Section 2.1 describes the generation of reaction mechanism, Section 2.2 outlines the optimisation techniques for obtaining the reaction rate parameters and Section 2.3 discusses the reduction of this kinetic model. The particle model for the formation of silica nanoparticles is described in Section 3. In Section 3.1 and Section 3.2, the internal variables of the particle model and various particle processes are described, followed by the optimisation of free parameters in the particle model in Section 3.3. The results of the kinetic and particle models are then presented and compared to experimental observations in Section 3.4. The paper concludes with a discussion of the model’s performance and suggests potential improvements to it, along with suitable steps for future research.

2. Kinetic model This section outlines the steps involved in the development of a kinetic model for the TEOS system. These steps consist of the generation of a gas phase reaction mechanism, obtaining the optimised rate parameters for these reactions and reducing this comprehensive mechanism to a succinct form suitable to be coupled to population balance equations in Section 3. 2.1. Gas-phase reaction mechanism The first step towards generating a reaction mechanism requires the identification of key intermediate species and formulating the possible reactions between them. Each reaction is described by both a forward and a backward rate expressions. The rate constant for the forward reaction is usually expressed in modified Arrhenius form: kf ¼ AT b expðEA =RTÞ,

ð1Þ

where A is the temperature independent pre-exponential factor, b is the temperature coefficient, EA is the activation energy, T is the temperature and R is the universal gas constant. For reversible

55

reactions, the backward rate is calculated from thermochemistry. Thermochemical information is therefore required for each species, and is usually expressed as a polynomial in temperature for specific heat, enthalpy, and entropy (Gordon and McBride, 1976). The thermochemistry of over 160 species in the TEOS system is reported by Phadungsukanan et al. (2009). Generating a reaction mechanism involving such a large number of species is computationally intractable. It is therefore necessary to make simplifying assumptions as outlined by West et al. (2007) in the case of TiO2 systems. As a first approximation, thermodynamic equilibrium calculations reported by Phadungsukanan et al. (2009) were used to provide clues as to which species are most significant, and thus likely to play a prominent role in the reaction mechanism. From the chemical equilibrium calculation, the most stable intermediate species in the temperature range 1100–1500 K and at atmospheric pressure were found to be ethoxy and methoxy silanes (Si(OH)x(OCH3)y(OC2H5)z, such that x þ y þz ¼ 4Þ. The most stable Si species at equilibrium for a temperature range of 800–2200 K are Si(OH)4, which has an equilibrium mole fraction orders of magnitude higher than any other monomer species. Gaseous SiO2 is highly unstable at equilibrium, and is thus likely to exist only in very low concentrations. The dimerisation and dehydration of Si(OH)4 are therefore a more favourable route for the formation of silica nanoparticles than the dimerisation of gaseous SiO2 molecules. The equilibrium calculations were used in this way to identify likely intermediate species and key steps in the reaction mechanism. A preliminary study of bond energies (Ho and Melius, 1995) was also performed to indicate the following order of bond strengths: (Si–O) 4 (O–C) 4 (C–C). Based on the results from the equilibrium calculation and the dissociation energies for different bonds in the TEOS molecule, a TEOS decomposition mechanism has been systematically developed. These ‘‘lumped’’ non-elementary reaction pathways can be grouped into key categories: (i) De-methylation: Removal of methylene radicals (:CH2) from adjacent ethoxy groups attached to the central Si atom to form methoxyethoxysilanes (Si(OH)x(OCH3)y(OC2H5)4  x  y). The TEOS molecule undergoes step-wise de-methylation to subsequently form Si(OH)4 (reactions R1–R4 in Table 2). (ii) De-ethylation: Removal of terminal C2H4 from each ethoxy branch to form ethoxysilanes (Si(OH)x(OC2H5)4  x) and ultimately form Si(OH)4 (reactions R15–R18 in Table 2). (iii) De-methylation and disproportionation: Removal of terminal methyl radicals (.CH3) to form methoxy silanol radicals Si(.OCH2)x(OC2H5)4  x. These radicals could either undergo a disproportionation reaction to form species with a vinyl group ð2CHQCH2 Þ attached to the central Si atom (reactions R10–R14 in Table 2). Alternatively, these radicals could undergo further de-methylation to form Si(OH)4 (reactions R5–R9 in Table 2). (iv) De-ethylation and dehydration: Removal of water or ethanol molecules from ethoxy silanes Si(OH)x(OC2H5)4  x to form the corresponding silicates ð2SiQOÞ. The reaction mechanism of Herzler et al. (1997) proposes this pathway to form OQSiQO in gas-phase. The reaction rules listed above have been used to develop a heuristic mechanism to achieve the step-wise decomposition of TEOS to eventually form a silicic acid molecule Si(OH)4. The main elementary reaction pathway proposed by Herzler et al. (1997) is also included in the mechanism as an alternative route for TEOS decomposition. The full mechanism can be roughly divided into two different sections: the TEOS decomposition with silicon intermediates described above and C2 carbon chemistry which is based on the ethanol mechanism as proposed by Marinov (1999).

56

S. Shekar et al. / Chemical Engineering Science 70 (2012) 54–66

All reactions are considered to be reversible to ensure thermodynamic consistency. As a starting point to perform flux and sensitivity analyses in Section 2.1.1, the rate parameters of the TEOS decomposition reactions are first assumed to have preexponential factor, A, as 1013 s  1, and both the temperature coefficient, b, and the activation energy, EA are assumed to be zero. The rate parameters are then refined by fitting them to experimental data from Herzler et al. (1997) in Section 2.2. The system of ordinary differential equations describing this mechanism is solved using an in-house open source gas-phase chemistry solver written in Cþþ called Sprog (Celnik, 2007). 2.1.1. Flux and sensitivity analyses The analysis of reaction fluxes enables us to draw important conclusions about the underlying chemical mechanism (Zse´ly et al., 2005). The instantaneous flux of element Z from species j to species k through reaction i at each reaction time step can be defined as: nZ,j nZ,k ri , NZ,i

Si,k ðtÞ ¼ ðAi =Xkmax Þ

ð2Þ

where nZ,j and nZ,k are the number of atoms of Z in a molecule of species j and k respectively. NZ,i is the sum of atoms of Z in all species on either side of reaction i and ri is the rate of reaction i. The instantaneous flux when integrated over the residence time of the reactor gives the time integrated flux. The time integrated flux of element Si for different species involved in TEOS decomposition is calculated for an initial mixture of 420 ppm TEOS in Ar at constant temperature and pressure (1200 1C, 2 bar). The resulting flux diagram is illustrated in Fig. 1. The thickness of each arrow in this diagram is proportional to the flux associated with the corresponding reaction. Fig. 1 shows that the net elemental flux of Si in the direction of Si(OH)4 formation is appreciably higher than in the direction of SiO2 formation. Sensitivity analysis is a quantitative tool to identify the importance of each reaction in a mechanism resulting in the formation and/or destruction of a given species. In kinetic model development, sensitivities are usually expressed as the logarithmic derivative of species concentration (Xk) with respect to each

@Xk ðtÞ , @Ai

ð3Þ

where Si,k ðtÞ is the sensitivity parameter of the kth species with respect to the ith reaction at time t, Ai is the pre-exponential factor of the ith reaction and Xk(t) is the mole fraction of the kth species at a time t. Sensitivity analysis has been implemented in Cþþ by solving a system of NR ðNS þ1Þ ODEs using the CVODES solver where NR is the number of reaction species and NS is the number of parameters with respect to which the sensitivities are to be calculated. Examination of the magnitude and sign of the sensitivity coefficients Si,k ðtÞ provides a useful quantitative tool for assessing the importance of the competing reactions i on the concentration of species k. Fig. 2 depicts the maximum normalised sensitivities of

0.06 Normalized Sensitivity Coefficient

Zi,j,k ¼

reaction coefficient Ai. The normalised sensitivity coefficient defined by Miller et al. (1985) is calculated as

Si(OH4) C2 H 4

0.05 0.04 0.03 0.02 0.01 0 −0.01 −0.02

High Low Sensitivity Sensitivity (R11−R14) (R15−R18)

High Sensitivity (R1−R10)

−0.03 −0.04 0

5

10

15

Si(OH)(OC2H5)3

Si(.OCH2)(OC2H5)3

Si(O.)(OH) (OC2H5)2 Si(OCH3)4

Si(OH)(.OCH2)(OC2H5)3

Si(OCH3)(OCH=CH2)(OC2H5)3 Si(O.)(OH)2(OC2H5) Si(.OCH2)(OCH3)(OCH=CH2)(OC2H5)3 Si(O.)(OH)3

Si(O.)(OH)2

Si(OCH3)(OCH=CH2)(OC2H5)3

Si(OH)3(OC2H5)

Si(OH)2(OCH3)2

Si(OH)(OCH3)2(OC2H5)3 Si(OH)2(.OCH2)(OC2H5)

SiO2

25

30

Fig. 2. Maximum normalised sensitivity coefficients for ethylene C2H4 and silicic acid Si(OH)4.

Si(OCH3)2(OC2H5)2

Si(OH)2(OC2H5)2

20

Reaction Number (R1−R30)

Si(OC2H5)4

Si(O.)(OC2H5)3

Low Sensitivity (R19−R30)

Si(OH)4 Fig. 1. Integrated flux diagram depicting main reaction pathways.

Si(OH)3(.OCH2)

S. Shekar et al. / Chemical Engineering Science 70 (2012) 54–66

the two main decomposition products, C2H4 and Si(OH)4, for each of the reactions for an initial mixture of 420 ppm TEOS in Ar at constant temperature and pressure (1200 1C, 2 bar). From sensitivity analysis calculations performed in this work, it is apparent that the maximum value of the normalised sensitivity parameter for Si(OH)4 and C2H4 to reactions R1–R10 and R15–R18 (see Table 2) is 10 orders of magnitude higher than that of elementary reactions proposed by Herzler et al. (1997) (R19–R30). Sensitivity to disproportionation reactions of methoxy silanol radicals to form species with vinyl groups (R11–R14) are also found to be low. Thus, flux and sensitivity analyses have enabled us to identify the main reaction pathways (R1–R18) for which the rate parameters are optimised as described in Section 2.2. Note: Herzler et al. (1997) analyse the gas phase products resulting from the shock-tube decomposition of TEOS and report the main product to be C2H4 followed by ethanol (C2H5OH). The reaction mechanism proposed in this study concludes that C2H4 is the main hydrocarbon product, and the formation of C2H5OH is not a direct consequence of TEOS decomposition. However, we hypothesise that the formation of ethanol can be explained by the presence of H2O in the system. We envisage that the particle processes release H2O back into the gas-phase which shifts the following equilibrium towards the direction of ethanol formation: C2 H5 OH þM"C2 H4 þ H2 Oþ M:

ðR31Þ

57

2.2.1. Low discrepancy series The method of LD series enables us to search for the best values of some objective function in a given domain, by doing repetitive evaluations of an objective function at several quasirandom data points in the domain. Data points in the 24dimensional parameter space of the current system are generated using the Halton (1960) low discrepancy series. The model is then evaluated at these points and the objective function

F1 ðxÞ ¼

Nexp X

ð/Yiexp S/Yisim SðxÞÞ2

ð5Þ

i¼1

computed, where Nexp is the number of experimental points from Herzler et al. (1997). /Yiexp S and /Yisim S are the experimental and simulated yield of C2 H4 respectively. The bounds for the free parameters in the current work were estimated based on physical considerations. For unimolecular reactions these are given by: 1011 r Ai r 1018 , and 0 rEAj =DHrxn,j r10, where DHrxn,j is the maximum enthalpy of reaction for the jth bond. Due to the large disparity between the upper and lower bounds, logarithmic transformations as described by Braumann et al. (in press) have been used, such that  1hk x xk ¼ xk,up k,low : ð6Þ xk,up

The current model is limited in its functionality concerning the treatment of ethanol. The validity of the current kinetic model can be justified by observing that the concentration of Si(OH)4 in the product is independent of ethanol generation.

hk A ½0,1 is the kth term given by the Halton sequence with kA f1, . . . ,Kg and K is the number of dimensions. xk,low and xk,up are respectively the lower and upper bounds for the kth parameter xk. The set of parameters x that minimises this objective function is then determined:

2.2. Estimation of rate parameters

xn ¼ argminfF1 ðxÞg:

To draw meaningful comparisons between model calculations and experimental measurements of species concentration profiles, it is necessary to have an accurate estimate of the reaction parameter values appearing in the modified Arrhenius expression, Eq. (1). A parameter fitting procedure is therefore adopted in order to reliably estimate these values, namely, the pre-exponential factor A, and the activation energy EA. The temperature coefficient b has been kept at a constant value of zero. The experimental results that are used for parameter estimation are reported by Herzler et al. (1997) for an initial composition of 420 ppm TEOS in Ar and a temperature range of 1210–1290 K. Due to a lack of experimental kinetics data from literature, model optimisation in this work has been performed using only one set of experimental results mentioned above. The availability of more literature data would enable the optimisation of multiple objective functions, which would result in better parameter estimation and applicability of the model to a wider range of process conditions. The parameter space of the fitting procedure used in this work has 24 dimensions, and contains elements x of the form

For this work, a low discrepancy series of 2000 data points in a 24-dimensional space is used to evaluate the objective function F1 ðxÞ described by Eq. (5) and the set of parameters xn that minimises F1 ðxÞ is determined.

x ¼ ðAi ,EAj Þ,

ð7Þ

x

ð4Þ

where i A f1, . . . ,18g parameterises the pre-exponential factors, Ai, for each of the 18 reactions and j A f1, . . . ,6g parameterises the activation energies, EAj , for each of the six bond types. The parameter fitting consists of the following steps: (i) Performing a pre-scan of the parameter space using a low discrepancy sequence. (ii) Construction of a quadratic response surface around the best points from LD sequence using a factorial design to produce a surrogate model. (iii) Optimisation of the surrogate model to minimise an objective function.

2.2.2. Response surface methodology The next step is the construction of a surrogate model around the best point obtained from the LD series. To facilitate this construction, it is desirable to reduce the dimension of the parameter space, which is realised by performing a sensitivity analysis (as described in Section 2.1.1). It is found that the model response (concentration of C2H4) is most sensitive to the preexponential factors (A) and activation energies (EA) of R1 and R15 defined in Table 2. The normalised sensitivity coefficients of the model response to these reactions are tabulated in Table 1. A surrogate model is then constructed using quadratic response surface methodology in a 4-dimensional space. This methodology has been proposed by Sheen et al. (2009) and used by Braumann et al. (2010) to optimise a granulation model. Prior to the construction of the surrogate model, it is convenient to perform a logarithmic transformation and normalise the para~ where x~ k A f1,0,1g and k A f1, . . . ,4g. The meter space x to x, ~ is bounds for x~ are assumed to be 7 1%. A surrogate model (ZðxÞ Table 1 Normalised sensitivity coefficients of C2H4 concentration with respect to reactions R1 and R15. Parameter

Normalised sensitivity coefficient

A1 EA1 A15 EA15

0.6563 1.2237 1.3056 3.1692

58

S. Shekar et al. / Chemical Engineering Science 70 (2012) 54–66

then constructed and the simulated C2H4 yield, /Yisim SðxÞ is approximated by a second order response surface as: ~ ¼ b0 þ /Yisim SðxÞ  ZðxÞ

4 X

bk x~ k þ

k¼1

4 X 4 X

bkl x~ k x~ l :

ð8Þ

Nexp X

ð½/Yiexp Smðx~ 0 ,cÞ2 þ½sexp si ðx~ 0 ,cÞ2 Þ: i

ð11Þ

i¼1

x~ ¼ x~ 0 þ cx,

ð9Þ

where x is Gaussian distributed with mean 0 and variance 1. The Gaussian distribution of the parameters x~ also influences the model ~ The mean m can be calculated (Braumann et al., 2010) response ZðxÞ. using: 4 X

F2 ðx~ 0 ,cÞ ¼

k ¼ 1 lZk

The coefficients b0 , bk , and bkl have been fitted by evaluating the model at 81 data points (3 data points for each of the 4 dimensions ¼ 34 data points). It is assumed that the parameters x~ are Gaussian distributed and have a mean x~ 0 and a standard deviation c:

mðx~ 0 ,cÞ ¼ Zðx~ 0 Þ þ

objective function used in the response surface optimisation also optimises the variances:

bkk ck2 :

ð10Þ

k¼1

It is also possible to express the variance in terms of the coefficients

b0 , bk , and bkl (Braumann et al., 2010). In contrast to Eq. (5), the

Minimising this objective function leads to the optimal set of parameters: n

ðx~ 0 ,cn Þ ¼ argminfF2 ðx~ 0 ,cÞg:

ð12Þ

x~ 0 ,c

Fig. 3 shows the product yield of ethylene as reported experimentally, and that obtained after the response surface optimisation. The bounds of uncertainties for the model response after optimisation are also presented. A relative error of 10% has been assumed for each experimental measurement. The rate parameters of the TEOS decomposition kinetic model obtained after optimisation are tabulated in Table 2 alongside their uncertainty bounds. The uncertainties in parameters for Reactions R1 and R15 are calculated using a response surface whilst uncertainties for other reaction parameters are assumed to be the upper and lower bounds of a priori values used in RSM.

2.6

[C2H4]/[TEOS consumed]

2.4

2.3. Gas phase mechanism reduction

Experimental values from Herzlera with error bounds Model response from optimised RSM Uncertainities in model response for optimised RSM

1.6

Prior to coupling the optimised kinetic model to the stochastic particle model, it is desirable to reduce the chemical mechanism to a concise form to minimise computational times and thus improve efficiency. The mechanism was reduced using the Level ˚ 2009) of Importance (LOI) technique as introduced by (Løvas, which combines time scale analysis with a sensitivity analysis with respect to a target parameter. Each species is assigned an LOI parameter, which is defined as

1.4

LOIi j ¼ Si,j  ti

1.2

where LOI ji is the level of importance of species i with respect to parameter j. Si,j is the sensitivity of species i to sensitivity target j and ti is the chemical lifetime of species i. The sensitivity of species i to sensitivity target j is calculated as     X @X   NR @Xi @Ak  Si,j ¼  i    ,  @Xj @Ak @Xj  k¼1

2.2 2 1.8

1 1210

1220

1230

1240 1250 1260 Temperature (K)

1270

1280

1290

Fig. 3. Optimised model response from response surface methodology and associated uncertainties. (a) Experimental values are from Herzler et al. (1997).

ð13Þ

Table 2 Kinetic model for TEOS decomposition proposed in the current work. The uncertainties in the rate parameters correspond to a relative error of 10% in experimental measurement. No.

Reaction pathway

log(A) (s  1)

b

EA =R (K  1)

R1 R2 R3 R4

SiðOC2 H5 Þ4 "SiðOCH3 Þ2 ðOC2 H5 Þ2 þ C2 H4 SiðOCH3 Þ2 ðOC2 H5 Þ2 "SiðOCH3 Þ4 þ C2 H4 SiðOCH3 Þ4 "SiðOHÞ2 ðOCH3 Þ2 þ C2 H4 SiðOHÞ2 ðOCH3 Þ2 "SiðOHÞ4 þ C2 H4

12.9 7 0.75% 12.2 7 1% 12.1 7 1% 11.9 7 1%

0 0 0 0

266517 2.87% 266517 10% 266517 10% 266517 10%

R5 R6 R7 R8 R9

SiðOC2 H5 Þ4 "Sið:OCH2 ÞðOC2 H5 Þ3 þ :CH3 Sið:OCH2 ÞðOC2 H5 Þ3 "SiðOHÞð:OCH2 ÞðOC2 H5 Þ2 þ C2 H4 SiðOHÞð:OCH2 ÞðOC2 H5 Þ2 "SiðOHÞ2 ð:OCH2 ÞðOC2 H5 Þ þ C2 H4 SiðOHÞ2 ð:OCH2 ÞðOC2 H5 Þ"SiðOHÞ3 ð:OCH2 Þ þ C2 H4 SiðOHÞ3 ð:OCH2 Þ"SiðOHÞ4 þ C2 H4

13.4 7 1% 15.3 7 1% 15.0 7 1% 14.8 7 1% 11.6 7 1%

0 0 0 0 0

499757 10% 27103 7 10% 27103 7 10% 27103 7 10% 442297 10%

R10 R11 R12 R13 R14

Sið:OCH2 ÞðOC2 H5 Þ3 "SiðOCH3 ÞðOCHQCH2 ÞðOC2 H5 Þ2 þ :H SiðOCH3 ÞðOCHQCH2 ÞðOC2 H5 Þ2 "Sið:OCH2 ÞðOCH3 ÞðOCHQCH2 ÞðOC2 H5 Þþ :CH3 Sið:OCH2 ÞðOCH3 ÞðOCHQCH2 ÞðOC2 H5 Þ"SiðOCH3 Þ2 ðOCHQCH2 Þ2 þ :H SiðOCH3 Þ2 ðOCHQCH2 Þ2 "SiðOHÞðOCH3 Þ2 ðOCHQCH2 Þþ C2 H2 SiðOHÞðOCH3 Þ2 ðOCHQCH2 Þ"SiðOHÞ2 ðOCH3 Þ2 þ C2 H2

14.9 7 1% 12.2 7 1% 14.6 7 1% 15.6 7 1% 15.3 7 1%

0 0 0 0 0

222397 10% 499757 10% 222397 10% 469177 10% 469177 10%

R15 R16 R17 R18

SiðOC2 H5 Þ4 "SiðOHÞðOC2 H5 Þ3 þ C2 H4 SiðOHÞðOC2 H5 Þ3 "SiðOHÞ2 ðOC2 H5 Þ2 þ C2 H4 SiðOHÞ2 ðOC2 H5 Þ2 "SiðOHÞ3 ðOC2 H5 Þ þ C2 H4 SiðOHÞ3 ðOC2 H5 Þ"SiðOHÞ4 þ C2 H4

14.3 7 0.99% 14.0 7 1% 13.9 7 1% 13.7 7 1%

0 0 0 0

312987 2.88% 312987 10% 312987 10% 312987 10%

S. Shekar et al. / Chemical Engineering Science 70 (2012) 54–66

where Xi is the mole fraction of species i and Xj is the mole fraction of the target species. For this case, Si(OH)4 was chosen as the target, because this is the product which we seek to accurately model. Ak is the pre-exponential reaction rate parameter of the kth reaction and NR are total number of reaction. The chemical lifetime of a species ti is approximated using the Jacobian matrix (Tura´nyi et al., 1993) as

ti ¼ 

1 Xi ¼ PN , 0 n00 Þn0 o Jii ð n k ¼ 1 ik ik ik k

Reactor Concentration (mol/cm3)

where Jii is the diagonal of the Jacobian and n0ik n00ik is the net stoichiometric coefficient of species i with respect to parameter k. Rigourous computational singular perturbation (CSP) analysis shows that the eigenvalues of the Jacobian matrix (as described above) give a more exact evaluation of the species’ lifetimes (Lam and Goussis, 1994). However, the eigenvalues become very similar to the inverse of Jacobian matrix diagonals as the lifetimes get shorter. This allows the two quantities to be used interchangeably with negligible error. A threshold value for the LOI of each species is defined; if a species has an LOI lower than the threshold, it is excluded from the mechanism. Every reaction associated with that species is also deemed unimportant and thus excluded. The full mechanism prior to reduction, contained 79 species and 401 reactions which can be roughly divided into two different sections: the TEOS decomposition with silicon intermediates and C2 carbon chemistry. The carbon chemistry was based on the ethanol mechanism as proposed by Marinov (1999). A skeletal mechanism was developed for process conditions identical to those in experiments by Seto et al. (1997) where they synthesised silica nanoparticles in a heated wall reactor at 900 1C with an initial TEOS concentration of 250 ppm and N2 as the diluent gas. The evolution of important species, such as TEOS, TMOS, C2H4, and Si(OH)4, all showed negligible error. This skeletal mechanism containing 27 species and 58 reactions was created using both LOI reduction and elemental flux analysis. The species that were removed using flux analysis were silicon product species whose evolution was orders of magnitude less than Si(OH)4, the main product species. Fig. 4 depicts the concentration profiles of key species with both full and reduced mechanisms and complete agreement is observed between the two. The kinetic model for TEOS decomposition obtained after LOI reduction consisting of a silicon submechanism and a hydrocarbon submechanism is provided in Appendix A.

59

3. Particle model The population balance model which describes the formation of silica nanoparticles is solved using a stochastic particle algorithm developed by Eibeck and Wagner (2000) and improved by Goodson and Kraft (2002) and Patterson and Kraft (2007). This model incorporates submodels that describe particle formation, growth and sintering. The chemical reactions of the precursor (TEOS) lead to the formation of Si(OH)4 monomers (as described in Section 2) which form particles by inception. These particles grow through surface reactions with monomer, in addition to condensation, coagulation and sintering (Ulrich, 1971; Pratsinis and Vemury, 1996; Morgan et al., 2005; Sander et al., 2009) resulting in the formation of web-like silica nanoparticles. An in-house developed open source solver written in Cþþ called MOPS (Celnik, 2007) is used to couple the kinetic model to the particle model. The details of the particle model used in the current work are outlined in this section. 3.1. Type space The collection of variables that completely describe the properties of a particle constitutes the type space of the entity. This is represented by a finite dimensional vector and is defined for the current system by: Pi ¼ Pi ðp1 , . . . ,pn ,C,I,SÞ,

ð14Þ

where iA f1, . . . ,Ng, and N is the total number of particles in the system. In this model it is assumed that each particle Pi consists of n primary particles pj, where j A f1, . . . ,ng. Both N and n vary with time. Each primary pj is described by its volume: pj ¼ pj ðvj Þ:

ð15Þ

The sum of the volumes of all the primaries is equal to the volume V of the particle: V¼

n X

vj :

ð16Þ

j¼1

C, I and S are square matrices describing the sintering between the primary particles (Sander et al., 2009). ½Ci,j contains the common surface of any two primaries pi and pj, ½Ii,j contains the sum of the initial surfaces of primary particles pi and pj before sintering and ½Si,j contains the surface of a sphere with the same volume as primary particles pi and pj. The dimension of these matrices is the number of primaries belonging to the particle. 3.2. Particle processes

10-10

The particle can undergo the following processes: 10-12

 Inception: We propose a new particle inception step that 10

-14

10-16 10-18 10-20 0.0001

Si(OH)4 Full mechanism Si(OH)4 Reduced Mechanism C2H4 Full mechanism C2H4 Reduced Mechanism TEOS Full Mechanism TEOS Reduced Mechanism

0.001

0.01

0.1

Time (s) Fig. 4. Species concentration profiles for full mechanism (401 reactions) and reduced mechanism (58 reactions) at experimental conditions of Seto et al. (1997).

proceeds by the interaction of two silicic acid molecules resulting in the loss of a water molecule and the formation of a dimer. This dimer is considered to be the first particle. In the current system, collision of two molecules in the gas phase leads to a dehydration reaction. However, in this model, inception has been treated as a coagulation event and a spherical particle consisting of one primary containing the volume of the two molecules is inserted into the particle ensemble. An inception event increases the total number of particles N by 1: moleculeþ molecule-PN þ 1 ðp1 ,C,I,SÞ:

ð17Þ

This can be represented by a chemical reaction as shown in Fig. 5.

60

S. Shekar et al. / Chemical Engineering Science 70 (2012) 54–66

H H

H H O

+

Si O

O

H

O

O

O O

O

H

Si O

H

H

O

H

Si

O H

O

H

H

H Si

H O

O

O

Si

-H2O

H

O H

H O O

H

O

H Si

O

O

O

H

H

O H

[primary particle]

[monomers] Fig. 5. Inception of a primary particle from gas-phase monomers.

Si

Si O

Si O

OH O

O

Si

Si O

H

H

O

HO OH

Si

-H2O

OH

Si

O

Si

HO

Si O

O

OH

Si Fig. 6. Surface reaction between a particle and a gas-phase molecule.

to the number of primaries of particle Pi plus the number of primaries of particle Pj. The matrices C, I and S are updated according to the method described by Sander et al. (2009) to represent the connectivity of particles Pi and Pj. The free molecular collision kernel (Eq. (23)) is used to calculate the coagulation rate between two particles Pi and Pj.

For an inception event, the sintering matrices C, I and S have dimension 1  1 and the only entry is 0. Due to the small size of the colliding gas-phase monomers, the free molecular kernel has been used to calculate the rate of inception as follows: Rinception ¼



1 fm 2 2 K NA C , 2

ð18Þ

where K fm is the free molecular regime coagulation kernel constant, NA is Avogadro’s number, C is the gas-phase concentration of the incepting species. Surface reaction: A surface reaction is the reaction of a gasphase molecule (Si(OH)4) on an active site on a particle Pi as shown in Fig. 6. A molecule from the gas-phase collides with a particle, sticks to its surface increasing its mass. In the current model, surface reactions are implemented so that the state space of a particle Pj is altered as follows: Pj þ molecule-Pj ðp1 , . . . ,pn ,C,I,SÞ:

ð19Þ

This is modelled by increasing the volume of a selected primary pi of particle Pj as pi ðvÞ-pi ðv þ dvÞ,

ð20Þ

where dv is the volume of the added gas-phase molecule. The rate of this process is calculated using the free molecular collision kernel for the particle and a single molecule from the gas phase. The rate takes the form rffiffiffiffiffiffiffiffiffiffiffi pkB T ðd þ dc Þ2 , Rsurfrxn ¼ 2:2aC ð21Þ 2m



where C, m and d are the gas-phase concentration, mass and collision diameter of the condensing species respectively. dc is the collision diameter of the particle and a the efficiency of the collision. Coagulation: The coagulation of particles is the collision and coalescence of two particles Pi and Pj to form a new particle Pk such that the primary particles of the new particle Pk are the same as those of Pi and Pj. This is implemented in the model as: Pi þ Pj -Pk ðp1 , . . . ,pnðPi Þ ,pðnðPi Þ þ 1Þ , . . . ,pnðPk Þ ,C,I,SÞ,

ð22Þ

where nðPi Þ is the number of primaries belonging to particle Pi and nðPk Þ ¼ nðPi Þ þ nðPj Þ the number of primaries belonging to particle Pk. Therefore the number of primary particles l is equal

K fm ði,jÞp



1 1 þ mðiÞ mðjÞ

1=2

ðdc ðiÞ þ dc ðjÞÞ2 ,

ð23Þ

where m is the mass and dc the collision diameter. dc of a particle P is calculated as dc ¼



6V ð1=df Þ n , A

ð24Þ

where the fractal dimension df is assumed to be 1.8 (Schaefer and Hurd, 1990). Sintering: The sintering process is modelled by calculating sintering levels for neighbouring primary particles (Sander et al., 2009). It is assumed that the excess agglomerate surface area over that of a spherical particle with the same mass decays exponentially (Xiong and Pratsinis, 1993). A sintering step reduces the surface area as well as the sphericity of a particle. The sintering level sði,jÞ between two primaries pi and pj is defined using the equation proposed by Sander et al. (2009) as sði,jÞ ¼

½Sij ½Cij

½S

 ½I ij

½S 1 ½I ij ij

ij

:

ð25Þ

The sintering level is updated by calculating the value of matrix C for each pair of neighbouring primaries (each non-zero entry of matrix CÞ: d½Ci,j 1 ¼  ð½Ci,j ½Si,j Þ: dt ts

ð26Þ

The characteristic sintering time ts is calculated using the formula by Tsantilis et al. (2001) as

ts ¼ AS  di,j  exp

   dp,min ES 1 , T di,j

ð27Þ

S. Shekar et al. / Chemical Engineering Science 70 (2012) 54–66

3.3. Estimation of sintering parameters The material dependant sintering parameters AS , ES and dp,min in Eq. (27) have been determined by fitting the model to the experiments of Seto et al. (1997) who generate silica particles at 900 1C and observe the change in primary and collision diameters by heating the particles to six different final temperatures. The kinetic parameters estimated in Section 2.2 for the experimental conditions of Herzler et al. (1997) are used for modelling the experimental conditions of Seto et al. (1997) because these process conditions are observed to be similar to each other. The simulation used to reproduce this experiment involves: (i) decomposition of TEOS to generate monomers and formation of particles in the first reactor at 900 1C from an initial gas mixture of TEOS and N2; (ii) heating the particles within 0.1 s to different final temperatures Tf ; (iii) allowing the particles to remain at temperature Tf for a certain period of time. The sintering parameters are optimised using a low discrepancy series as described in Section 2.2. The 3-dimensional parameter space is given by: y ¼ ðAS ,ES ,dp,min Þ:

ð28Þ

The objective function is given by:

FðyÞ ¼

Nexp X

sim 2 exp sim 2 ð½/dexp ci S/dci SðyÞ þ ½/dpi S/dpi SðyÞ Þ,

ð29Þ

i¼1

where Nexp is the number of experimental points from Seto et al. exp sim sim (1997). /dexp ci S, /dpi S, /dci S and /dpi S are the experimental and simulated collision and primary particle diameters respectively. The set of parameters that minimises Eq. (29) is determined as n

y ¼ argminfFðyÞg:

ð30Þ

Parameter

Value (current model)

Value (Sander et al., 2009)

AS (m  1 s) ES (K) dp,min (m)

1.16  10  16 12.7  104 4.4  10  9

6.5  10  15 10.2  104 10.2  10  9

120

Average particle diameter (nm)

The particles’ processes described above are implemented for an ensemble of stochastic particles using a Linear Process Deferment Algorithm (Patterson et al., 2006). The kinetic model described in Section 2 and the particle model described in this section are coupled using an operator splitting technique developed by Celnik et al. (2007, 2009) to generate the overall model.

Table 3 Values of sintering parameters in the particle model compared to values calculated by Sander et al. (2009).

dp from experiments(a) dc from experiments(a) dp from model (Full mechanism) dc from model (Full mechanism) dp from model (Skeletal mechanism) dc from model (Skeletal mechanism)

100 80 60 40 20 0 800

1000

1200 1400 Temperature (°C)

1600

1800

Fig. 7. Optimised model response, i.e., average primary diameter (dp Þ and collision diameter (dc Þ at different final temperatures. (a) Experimental values are from Seto et al. (1997).

10−9 10−10 Concentration (mol/cm3)

where di,j is the minimum diameter of the two neighbouring primaries pi and pj. dp,min is the critical diameter below which the primaries are assumed to be liquid like and they sinter instantaneously. The parameters AS , ES and dp,min are free parameters in the model. The sintering time is determined by the smallest of two neighbouring primaries (Yadha and Helble, 2004). Primaries pi and pj are completely sintered if the sintering level sði,jÞ (Eq. (25)) is greater than 0.95. This is included in the model by replacing primary pi and pj by a new primary pk. The volume of the new primary pk is the sum of the volumes of the primaries pi and pj and the neighbouring primaries of the new primary pk are the neighbours of the primaries pi and pj.

61

10−11 10−12 10−13 10−14 10−15 0

0.01 0.02 0.03 0.04 0.05 0.06 0.07 0.08 0.09 0.1 Time (s)

y

In this work, 2000 data points are generated using Halton low discrepancy series (Halton, 1960) in the 3-dimensional parameter space given by y. The objective function FðyÞ is evaluated at these data points and the set of parameters yn that minimises FðyÞ is determined. The resulting values of yn are shown in Table 3 and are similar to those reported by Sander et al. (2009). The underprediction in the value of AS is envisaged to be due to the inclusion of the chemistry of TEOS decomposition in the particle model, which leads to an increase in the residence time of the reactor. The optimised model response, i.e., the average collision diameter and particle diameter at different final temperatures is

Fig. 8. Gas phase reactor composition versus time. Initial composition is 250 ppm TEOS in N2, T¼900 1C and P ¼1.01 bar.

shown in Fig. 7. The model response using full and reduced mechanisms (Section 2.3) are also compared and perfect agreement is observed between the two. 3.4. Analysis of results The values of the free parameters in Section 3.3 have been estimated by fitting the model to the experimental values of average primary and collision diameters from Seto et al. (1997). To establish the consistency of the model with experiments, in

62

S. Shekar et al. / Chemical Engineering Science 70 (2012) 54–66

this subsection we analyse various other model results and compare them with experimental observations. 3.4.1. Gas phase reactor composition The time evolution of the prominent gaseous silicon species in the reactor for the first 100 ms is shown in Fig. 8. The first reactor is at constant temperature and pressure (900 1C, 1.01 bar) with an initial mixture of 250 ppm TEOS in N2. The residence time tr in the first reactor is set to 1 s in the following simulations. Fig. 8 suggests that reaction proceeds via the step-wise decomposition of TEOS to form ethoxy or methoxy silanols as intermediates which undergo further decomposition to ultimately form Si(OH)4. These Si(OH)4 monomers act as seeds for particle formation and are rapidly consumed from the gas phase to particulate phase. 80 Reactor 1

Reactor 2

Average particle diameter (nm)

70

3.4.2. Temporal evolution of particle size Chemical reactions in gas-phase occur in the first reactor and Si(OH)4 monomers are formed as depicted in Fig. 8. In the first reactor, particle formation occurs as a result of collision of two silicic acid molecules to form one primary particle with the volume of the formed dimer. When introduced into the second reactor at high temperature, the particles rapidly start to sinter. The time evolution of particles at different final temperatures (Tf) is shown in Fig. 9. In the first reactor at a temperature of 900 1C, the primary particles, once formed, remain constant in size. However, as the particles are introduced in the second reactor, the primary particle diameter rapidly increases indicating greater levels of sintering at higher temperatures. The collision diameter of the particle increases in the first reactor due to surface reaction and coagulation. As soon as the particles are introduced in the second reactor, the sintering process predominates, reducing the overall size of the particle. At a temperature of about 2000 1C, the collision diameter and primary diameter become identical and particles become spherical, thus no further sintering is possible.

Tf = 900°C

60

Tf = 1500°C

50

Collision diameter (s) Tf = 2000°C

40 30 Primary particle diameter (s)

20

Tf = 1500°C Tf = 900°C

10

3.4.3. Particle size distributions The size distributions with respect to primary particle diameter (dp Þ and collision diameter (dc Þ at three different final temperatures are given in Fig. 10. At 900 1C, the primary particle size distribution has a small standard deviation and mean. As the furnace temperature increases, the mean and standard deviation also increase. The corresponding collision diameter distribution gets narrower with increasing temperature due to the sintering of the particles. In the comparison it is assumed that the mobility diameter (experimental values from Seto et al. (1997)) is equal to the collision diameter.

0 0.4

0.6

0.8

1 1.2 Time (s)

1.4

1.6

1.8

2

2.2

Fig. 9. Time evolution of average particle size at different final temperatures.

T = 900°C

T = 1500°C

100

10−1

100

102

100

100

10−1

10−1

dp (nm)

102

100

dp (nm) 101

101

100

100

100

10−1 10−2 10−3 100

10−1

104

10−3 100

10−1 10−2

10−2

102 dc (nm)

102 dp (nm)

101

dN/dlog (dc)

dN/dlog (dc)

T = 1750°C 101

101

dN/dlog (dp)

dN/dlog (dp)

101

100

3.4.4. Computer generated TEM images The computer TEM images are generated by assuming that each particle is created by sticking together its constituting primary particles in a random fashion. These particles are projected onto a

dN/dlog (dp)

0.2

dN/dlog (dc)

0

102 dc (nm)

104

10−3 100

102 dc (nm)

104

Fig. 10. Particle size distribution at different furnace temperatures. The circles are experimental points from Seto et al. (1997) and the solid lines are from simulation.

S. Shekar et al. / Chemical Engineering Science 70 (2012) 54–66

63

Fig. 11. TEM style images at different final temperatures from simulation and experiment. (a) Experimental images are from Seto et al. (1997).

plane to simulate a collection of particles as generated using a TEM grid in an experiment. Computer generated TEM style images are compared to experimental values from Seto et al. (1997) as demonstrated in Fig. 11. From Fig. 11, it is apparent that the particles reduce in size with increasing temperature and become rounder. The particles are nearly spherical at 1750 1C. The particle size distributions and the TEM images calculated by the model are found to be in close agreement with the experimental data.

4. Conclusion A comprehensive model has been generated for the flame synthesis of silica nanoparticles from TEOS. Using recently published thermochemical data for TEOS intermediates, the present work outlines a thermodynamically consistent kinetic model for decomposition of TEOS generated from first principles. Initially, a rudimentary kinetic model was constructed which was found to be at least qualitatively correct, and suggested that silicic acid Si(OH)4 was the main product of the decomposition. Flux and sensitivity analyses were used to identify the dominant reaction pathways, which led to the conclusion that the decomposition of TEOS to form silicates is insignificant under the process conditions of current work. The model predictions, i.e., the final product concentrations, were sensitive to the assumed pre-exponential rate parameters, and therefore, in order to improve the performance of the model, rate parameter optimisation was carried out. Low discrepancy series and response surface methodology were applied in order to estimate the rate parameters of the model and to determine the uncertainties in model parameters given the

experimental error. A detailed transition state search analysis is an essential next step in order to determine exact values of the rate parameters. The kinetic model developed was then coupled to a stochastic particle model that includes the processes through which Si(OH)4 monomers ultimately form silica aggregates. The material dependent sintering parameters in the particle model were estimated by fitting the model to experimental values of collision and primary particle diameters using low discrepancy series. The particle size distributions and computer-generated TEM images at different temperatures were compared with experimental values and good agreement was found. Work is currently underway to develop a detailed population balance model, where the particle jump processes are considered in conjunction with chemical reactions and the chemical units (i.e., numbers of Si, O and H units) within each particle in the population are tracked. By coupling a new kinetic model, generated from quantum chemical calculations, to a new population balance model with primary particle tracking, the current work demonstrates the feasibility of using first principles modelling to understand complex particle synthesis processes.

Nomenclature Greek characters

f

Z b

objective function surrogate model parameter of response surface

64

m s n a t

S. Shekar et al. / Chemical Engineering Science 70 (2012) 54–66

mean (model prediction) standard deviation (uncertainty in model response) net stoichiometric coefficient collision efficiency chemical life time (s)

Roman characters kf T A b EA R nZ r S X NR x Nexp Y x~ x~ 0 c LOI J P p N n V v A ½C, ½I, ½S Rinception Kfm NA C Rsurfrxn kB m d dc dp df

rate constant of forward chemical reaction temperature (K) pre-exponential factor (s  1) temperature coefficient (–) activation energy (J mol  1] gas constant (J mol  1 K  1) number of atoms of element Z (–) rate of reaction (mol m  3 s  1) normalised sensitivity coefficient (–) mole fraction (–) number of reactions parameter space for kinetic model optimisation number of experimental points model response normalised parameter space for response surface mean of normalised parameter space standard deviation of normalised parameter space level of importance parameter (s) Jacobian matrix type space of a particle type space of a primary particle number of particles in a system number of primary particles in a particle volume of a particle (m3) volume of a primary particle (m3) surface area of a particle (m2) sintering matrices (m2) rate of inception (m  3 s  1) free molecular collision kernel (m3 s  1) Avogadro constant (mol  1) concentration of gas phase (mol m  3) rate of surface reaction (m  3 s  1) Boltzmann constant (J K  1) mass (kg) diameter (m) collision diameter (m) primary particle diameter (m) fractal dimension (–)

s AS ES dp,min y

sintering level (–) sintering rate parameter (m  1 s) sintering temperature parameter (K) minimum diameter below which sintering is instantaneous (m) parameter space for sintering parameters optimisation

Superscripts and subscripts exp sim i, j, k low up

experimental values simulated values Counting variables lower limit upper limit

Abbreviations TEOS FSP TMOS BAC MP4 LOI LD RSM PSD TEM

tetraethoxysilane flame spray pyrolysis tetramethoxysilane bond additivity correction Møller–Plesset theory, 4th order level of importance low discrepancy response surface methodology particle size distribution transmission electron microscopy

Acknowledgments The authors acknowledge financial support provided by HeiQ materials, the Cambridge Commonwealth Trust and Murray Edwards BP Centenary bursary. MK is grateful for the support of the University of Duisburg-Essen. The authors also thank members of the Computational Modelling Group for their guidance and support.

Appendix A Kinetic model for TEOS decomposition Table A.1 depicts the skeletal mechanism for the gas phase decomposition of TEOS with N2 as a diluent. This mechanism consists of 58 reactions and 27 species.

Table A.1 Reaction mechanism for TEOS decomposition. The references marked [*] are reactions proposed in this study. No.

Reaction pathway

log(A)

Si submechanism 1 SiðOC2 H5 Þ4 "SiðOCH3 Þ2 ðOC2 H5 Þ2 þ C2 H4 2 SiðOCH3 Þ2 ðOC2 H5 Þ2 "SiðOCH3 Þ4 þ C2 H4 3 SiðOCH3 Þ4 "SiðOHÞ2 ðOCH3 Þ2 þ C2 H4 4 SiðOHÞ2 ðOCH3 Þ2 "SiðOHÞ4 þ C2 H4 5 SiðOC2 H5 Þ4 "SiðOHÞðOC2 H5 Þ3 þ C2 H4 6 SiðOHÞðOC2 H5 Þ3 "SiðOHÞ2 ðOC2 H5 Þ2 þ C2 H4 7 SiðOHÞ2 ðOC2 H5 Þ2 "SiðOHÞ3 ðOC2 H5 Þþ C2 H4 8 SiðOHÞ3 ðOC2 H5 Þ"SiðOHÞ4 þ C2 H4

12.5 12.2 12.1 11.9 14.3 14.0 13.9 13.7

Hydrocarbon submechanism 9 H þ H þ M"H2 þ M 10 H þ H þ H2 "H2 þ H2 11 CH3 þ H þ M"CH4 þ M Low pressure limit: Troe parameters: a¼ 0, T nnn ¼ 1:00E15,T n ¼ 1:00E15,T nn ¼ 40 12 CH4 þ H"CH3 þ H2

12 11 15 31 4

EA =R (K  1)

Reference

0 0 0 0 0 0 0 0

26651 26651 26651 26651 31298 31298 31298 31298

[n] [n] [n] [n] [n] [n] [n] [n]

1 1 0 4 3

0 0 0 1061 4406

b

Marinov et al. (1996) Marinov et al. (1996) Marinov et al. (1996) Marinov et al. (1996)

S. Shekar et al. / Chemical Engineering Science 70 (2012) 54–66

65

Table A.1 (continued ) No. 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58

Reaction pathway

log(A)

b

EA =R (K  1)

Reference

CH3 þ H"CH2 þ H2 CH3 þ M"CH þ H2 þ M CH3 þ M"CH2 þ H þ M CH2 þ H"CH þ H2 CH2 þ CH3 "C2 H4 þ H CH2 þ CH2 "C2 H2 þ H CH2 þ C2 H2 "H2 CCCH þ H CH2 ðSÞþ M"CH2 þ M CH2 ðSÞþ CH4 "CH3 þ CH3 CH2 ðSÞþ H2 "CH3 þ H CH2 ðSÞþ C2 H2 "H2 CCCH þ H CH2 ðSÞþ C2 H4 "AC3 H5 þ H CH2 ðSÞþ H"CH þ H2 CH2 ðSÞþ CH3 "C2 H4 þ H CH þ C2 H2 "C3 H2 þ H CH þ CH2 "C2 H2 þ H CH þ CH3 "C2 H3 þ H CH þ CH4 "C2 H4 þ H C2 H5 þ H"C2 H4 þ H2 C2 H5 þ H"CH3 þ CH3 C2 H4 þ CH3 "C2 H3 þ CH4 C2 H4 þ H"C2 H3 þ H2 C2 H4 þ H þ M"C2 H5 þ M Low pressure limit: Troe parameters: a¼ 1, T nnn ¼ 1:00E15,T n ¼ 95,T nn ¼ 200 C2 H4 "C2 H2 þ H2 C2 H3 þ H þ M"C2 H4 þ M Low pressure limit: Troe parameters: a¼ 0.782, T nnn ¼ 208,T n ¼ 2663,T nn ¼ 6095 C2 H3 þ H"C2 H2 þ H2 C2 H3 þ C2 H"C2 H2 þ C2 H2 C2 H3 þ CH"CH2 þ C2 H2 C2 H3 þ CH3 "AC3 H5 þ H C2 H3 þ CH3 "C3 H6 C2 H3 þ CH3 "C2 H2 þ CH4 C2 H2 þ CH3 "C2 H þ CH4 C2 H2 þ M"C2 H þ H þ M C2 H2 þ H þ M"C2 H3 þ M Low pressure limit: Troe parameters: a¼ 1, T nnn ¼ 1:00E15,T n ¼ 675,T nn ¼ 1:00E þ 15 C2 H þ H2 "C2 H2 þ H C3 H6 "C2 H2 þ CH4 PC3 H5 þ H"C3 H6 C3 H6 þ H"C2 H4 þ CH3 C3 H6 þ H"AC3 H5 þ H2 C3 H6 þ H"PC3 H5 þ H2 C3 H6 þ CH3 "AC3 H5 þ CH4 C3 H6 þ CH3 "PC3 H5 þ CH4 AC3 H5 þ H"C3 H6 PC3 H5 þ H"AC3 H5 þ H H2 CCCH þ H"C3 H2 þ H2 H2 CCCH þ CH3 "C3 H2 þ CH4

14 15 16 18 14 14 13 13 14 14 14 14 13 13 14 14 13 14 14 13 1 6 12 34 5 13 30 14 13 14 3 57 13 11 17 11 40 6 12 14 13 5 6 0 0 26 14 14 13

0 0 0 2 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 4 6 0 5 0 0 4 0 0 0 4  13 0 0 0 1 7 2 0 0 0 3 3 4 4 4 0 0 0

7604 41528 46031 0 0 0 3323 0 0 0 0 0 0 0 0 0 0 0 4028 0 4784 852 917 2240 0 141 1672 0 0 0 2859 6982 0 8706 53880 1304 3312 435 35249 0 656 1255 6186 2858 6470 2753 0 504 0

Marinov Marinov Marinov Marinov Marinov Marinov Marinov Marinov Marinov Marinov Marinov Marinov Marinov Marinov Marinov Marinov Marinov Marinov Marinov Marinov Marinov Marinov Marinov

References Braumann, A., Man, P.L.W., Kraft, M., 2010. Statistical approximation of the inverse problem in multivariate population balance modelling. Ind. & Eng. Chem. Res. 49, 428–438. Braumann, A., Man, P.L.W., Kraft, M. The inverse problem in granulation modelling – two different statistical approaches. AIChE Journal, in press. doi:10.1002/aic. 12526. Celnik, M.S., 2007. On the Numerical Modelling of Soot and Carbon Nanotube Formation. Ph.D. Thesis. University of Cambridge. Celnik, M.S., Patterson, R.I.A., Kraft, M., Wagner, W., 2009. A predictor-corrector algorithm for the coupling of stiff odes to a particle population balance. J. Comp. Phys. 228, 2758–2769. Celnik, M.S., Patterson, R.I.A., Kraft, M., Wagner, W., 2007. Coupling a stochastic soot population balance to gas-phase chemistry using operator splitting. Combust. Flame 148, 158–176. Chu, J.C., Breslin, J., Wang, N., Lin, M., 1991. Relative stabilities of tetramethyl orthosilicate and tetraethyl orthosilicate in the gas phase. Mater. Lett. 12, 179–184. Eibeck, A., Wagner, W., 2000. An efficient stochastic algorithm for studying coagulation dynamics and gelation phenomena. SIAM J. Sci. Comput. 22, 802–821. Goodson, M., Kraft, M., 2002. An efficient stochastic algorithm for simulating nanoparticle dynamics. J. Comput. Phys. 183, 210–232.

et et et et et et et et et et et et et et et et et et et et et et et

al. al. al. al. al. al. al. al. al. al. al. al. al. al. al. al. al. al. al. al. al. al. al.

(1996) (1996) (1996) (1996) (1996) (1996) (1996) (1996) (1996) (1996) (1996) (1996) (1996) (1996) (1996) (1996) (1996) (1996) (1996) (1996) (1996) (1996) (1996)

Marinov et al. (1996) Marinov et al. (1996) Marinov Marinov Marinov Marinov Marinov Marinov Marinov Marinov Marinov

et et et et et et et et et

al. al. al. al. al. al. al. al. al.

(1996) (1996) (1996) (1996) (1996) (1996) (1996) (1996) (1996)

Marinov et al. (1996) Marinov et al. (1996) Marinov et al. (1996) Marinov et al. (1996) Marinov et al. (1996) Marinov et al. (1996) Marinov et al. (1996) Marinov et al. (1996) Marinov et al. (1996) Marinov et al. (1996) Pauwels et al. (1995) Marinov et al. (1997)

Gordon, S., McBride, B.J., 1976. Computer Program for Calculation of Complex Chemical Equilibrium Composition, Rocket Performance, Incident and Reflected Shocks and Chapman–Jouget Detonations. Technical Report. NASA Washington DC. Halton, J.H., 1960. On the efficiency of certain quasi-random sequences of points in evaluating multi-dimensional integrals. Numer. Math. 2, 84–90. Herzler, J., Manion, J.A., Tsang, W., 1997. Single-pulse shock tube study of the decomposition of tetraethoxysilane and related compounds. J. Phys. Chem. A 101, 5500–5508. Ho, P., Melius, C.F., 1995. Theoretical study of the thermochemistry of molecules in the Si–O–H–C system. J. Phys. Chem. 99, 2166–2176. Iler, R.K., 1979. The Chemistry of Silica: Solubility, Polymerization, Colloid and Surface Properties and Biochemistry of Silica. WILEY, Toronto, Canada. Jang, H.D., 1999. Generation of silica nanoparticles from tetraethylorthosilicate (TEOS) vapor in a diffusion flame. Aerosol Sci. Tech. 30, 477–488. Lam, S.H., Goussis, D.A., 1994. The CSP method for simplifying kinetics. Int. J. Chem. Kin. 26, 461–486. ˚ T., 2009. Automatic generation of skeletal mechanisms for ignition combusLøvas, tion based on level of importance analysis. Combust. Flame 156, 1348–1358. Marinov, N.M., 1999. A detailed chemical kinetic model for high temperature ethanol oxidation. Int. J. Chem. Kin. 31, 183–220. Marinov, N.M., Castaldi, M.J., Melius, C.F., Tsang, W., 1997. Aromatic and polycyclic aromatic hydrocarbon formation in a premixed propane flame. Comb. Sci. Tech. 128, 295–342.

66

S. Shekar et al. / Chemical Engineering Science 70 (2012) 54–66

Marinov, N.M., Pitz, W.J., Westbrook, C.K., Castaldi, M.J., Senkan, S.M., 1996. Modeling of aromatic and polycyclic aromatic hydrocarbon formation in premixed methane and ethane flames. Comb. Sci. Tech. 116, 211–287. Miller, J.A., Branch, M.C., McLean, W.J., Chandler, D.W., Smooke, M.D., Kee, R.J., 1985. The conversion of HCN to NO and N2 in H2–O2–HCN–Ar flames at low pressure. Int. Symp. Comb. 20, 673–684. Morgan, N., Wells, C., Kraft, M., Wagner, W., 2005. Modelling nanoparticle dynamics: coagulation, sintering, particle inception and surface growth. Combust. Theor. Model. 9, 449–461. Patterson, R.I.A., Kraft, M., 2007. Models for the aggregate structure of soot particles. Combus. Flame 151, 160–172. Patterson, R.I.A., Singh, J., Balthasar, M., Kraft, M., Norris, J., 2006. The linear process deferment algorithm: a new technique for solving population balance equations. SIAM J. Sci. Comput. 28, 303–320. Pauwels, J.F., Volponi, J.V., Miller, J.A., 1995. The oxidation of allene in a lowpressure H2/O2/Ar–C3H4 flame. Comb. Sci. Tech. 110, 249–276. Phadungsukanan, W., Shekar, S., Shirley, R., Sander, M., West, R.H., Kraft, M., 2009. First-principles thermochemistry for silicon species in the decomposition of tetraethoxysilane. J. Phys. Chem. A 113, 9041–9049. Pratsinis, S.E., 1988. Simultaneous nucleation, condensation, and coagulation in aerosol reactors. J. Colloid Interface Sci. 124, 416–428. Pratsinis, S.E., Vemury, S., 1996. Particle formation in gases: a review. Powder Technol. 88, 267–273. Sander, M., West, R.H., Celnik, M.S., Kraft, M., 2009. A detailed model for the sintering of polydispersed nanoparticle agglomerates. Aerosol Sci. Tech. 43, 978–989. Schaefer, D.W., Hurd, A.J., 1990. Growth and structure of combustion aerosols: fumed silica. Aerosol Sci. Tech. 12, 876–890. Seto, T., Hirota, A., Fujimoto, T., Shimada, M., Okuyama, K., 1997. Sintering of polydisperse nanometer-sized agglomerates. Aerosol Sci. Tech. 27, 422–438.

˚ T., 2009. Spectral uncertainty quantification, Sheen, D.A., You, X., Wang, H., Løvas, propagation and optimization of a detailed kinetic model for ethylene combustion. Proc. Combust. Inst. 32, 535–542. Slowing, I.I., Trewyn, B.G., Giri, S., Lin, V.S.Y., 2007. Mesoporous silica nanoparticles for drug delivery and biosensing applications. Adv. Funct. Mater. 17, 1225–1236. Smolı´k, J., Moravec, P., 1995. Gas phase synthesis of fine silica particles by oxidation of tetraethylorthosilicate vapour. J. Mater. Sci. Lett. 14, 387–488. Stark, W.J., Pratsinis, S.E., 2002. Aerosol flame reactors for manufacture of nanoparticles. Chem. Eng. Prog. 126, 103–108. Takeuchi, H., Izami, H., Kawasaki, A., 1994. Decomposition study on TEOS. Mater. Res. Soc. Symp. Proc. 334, 45–55. Tsantilis, S., Briesen, H., Pratsinis, S.E., 2001. Sintering time for silica particle growth. Aerosol Sci. Tech. 34, 237–246. Tura´nyi, T., Tomlin, A.S., Pilling, M.J., 1993. On the error of the quasi-steady-state approximation. J. Phys. Chem. 97, 163–172. Ulrich, G.D., 1971. Theory of particle formation and growth in oxide synthesis flames. Comb. Sci. Tech. 4 (1), 47–57. West, R.H., Celnik, M.S., Inderwildi, O.R., Kraft, M., Beran, G.J.O., Green, W.H., 2007. Towards a comprehensive model of the synthesis of TiO2 particles from TiCl4. Ind. Eng. Chem. Res. 46, 6147–6156. Xiong, Y., Pratsinis, S.E., 1993. Formation of agglomerate particles by coagulation and sinterings – part I. A two-dimensional solution of the population balance equation. J. Aerosol Sci. 24, 283–300. Yadha, V., Helble, J., 2004. Modeling the coalescence of heterogenous amorphous particles. J. Aerosol Sci. 35, 665–681. Zse´ly, I.G., Vira´g, I., Tura´nyi, T., 2005. Investigation of a Methane Oxidation Mechanism via the Visualization of Element Fluxes. Technical Report. Depart¨ os ¨ University (ELTE), Budapest, Hungary. ment of Physical Chemistry, Eotv