Synthesis of silicon nanoparticles with a narrow size distribution: A theoretical study

Journal of Aerosol Science 44 (2012) 46–61

Contents lists available at SciVerse ScienceDirect

Journal of Aerosol Science journal homepage: www.elsevier.com/locate/jaerosci

Synthesis of silicon nanoparticles with a narrow size distribution: A theoretical study William J. Menz a, Shraddha Shekar a, George P.E. Brownbridge a, Sebastian Mosbach a, b ¨ , Wolfgang Peukert b, Markus Kraft a,n Richard Kormer a b

Department of Chemical Engineering and Biotechnology, University of Cambridge, New Museums Site, Pembroke Street, Cambridge CB2 3RA, United Kingdom Institute of Particle Technology, Friedrich-Alexander-Universit¨ at Erlangen-N¨ urnberg, Cauerstraße 4, Erlangen 91058, Germany

a r t i c l e in f o

abstract

Article history: Received 3 May 2011 Received in revised form 28 September 2011 Accepted 14 October 2011 Available online 25 October 2011

This work presents a study of the processes involved in synthesis of narrowly distributed silicon nanoparticles from the thermal decomposition of silane. Two models are proposed, one which simultaneously solves the kinetic mechanism of Swihart & Girshick (1999, Journal of Physical Chemistry B 103, 64–76) while adjusting the sintering parameters; and another which adjusts the kinetic and surface growth mechanisms while neglecting coagulation and sintering. The models are applied to simulate the ¨ centreline of the hot-wall reactor and process conditions of Kormer et al. (2010, Journal of Aerosol Science 41, 998–1007). Both models are shown to give good agreement with experimental PSDs at a range of process conditions. However, it is reported that an unphysical sintering process is obtained when attempting to use Swihart & Girshick’s kinetic mechanism, while solving for the sintering parameters. The model with adjusted gas-phase and surface growth processes gives better quantitative and qualitative agreement with experimental results. It is therefore recommended that further study into the kinetic and heterogeneous growth mechanisms be conducted in order to better understand the fundamental processes occurring in this hot-wall reactor. & 2011 Elsevier Ltd. All rights reserved.

Keywords: Silicon Nucleation Parameter Estimation nanoparticles Size distribution Particle model

1. Introduction When synthesised to a controlled diameter, silicon nanoparticles can be used as photoluminescent standards or as ¨ ‘building blocks’ for use in hierarchically ordered systems (Kormer et al., 2010; Li et al., 2003). Understanding the mechanism by which they are formed is therefore important in developing improved techniques of manufacture. The mechanism is thought to consist of complex interconnected processes including gas phase reactions, nucleation, ¨ coagulation, condensation and surface reaction (Dang & Swihart, 2009; Kormer et al., 2010; Nijhawan et al., 2003; Swihart & Girshick, 1999). The gas-phase mechanism for the thermal decomposition of silane has attracted considerable attention over the past 80 years. The mechanism was first studied by Hogness, Wilson & Johnson in 1936 who reported that silane decomposed by a simple single-step reaction into only pure silicon and hydrogen (Hogness et al., 1936). This view was subsequently modified by Purnell and Walsh, who proposed that it was a multi-body reaction which could also form di- and tri-silanes

n

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).

0021-8502/$ - see front matter & 2011 Elsevier Ltd. All rights reserved. doi:10.1016/j.jaerosci.2011.10.005

W.J. Menz et al. / Journal of Aerosol Science 44 (2012) 46–61

47

(Purnell & Walsh, 1966). Since then, the mechanism has been continuously expanded with larger and more complex processes involving the formation of silenes (doubly bound silane species, e.g. SiH2SiH2, suffixed ‘A’), silylenes (radical species, e.g. SiH3SiH, suffixed ‘B’), silyl anions and three-dimensional polycyclic silicon hydride species (Bhandarkar et al., 2000; Giunta et al., 1990; Swihart & Girshick, 1999; Vepˇrek et al., 1993). Ho et al. (1994) used laser-induced fluorescence measurements of silicon atoms during the chemical vapour deposition (CVD) of silicon from silane at atmospheric pressure to determine kinetic parameters for the formation of silicon hydride species with up to Si3H8. This work formed the basis for the mechanism of Swihart & Girshick, who used the parameters of Ho et al. and reaction rules to generate reactions for silane, silene and silylene species up to 20 silicon atoms. The result was a total of 2615 reactions with 221 silicon hydride species (Swihart & Girshick, 1999). More recently, automatic mechanism generation was used to generate a mechanism with 8076 reactions among 1398 species for clusters up to eight silicon atoms (Wong et al., 2004) and quantum chemical calculations were undertaken to investigate the molecular transition steps involved in the decomposition of silane (Adamczyk et al., 2009, 2010). Modelling the formation of silicon nanoparticles also requires treatment of the particle growth processes as well as the gas-phase chemistry. Swihart et al. used their chemical mechanism (Swihart & Girshick, 1999) as the basis for aerosol dynamics models including particle growth by surface reactions, coagulation and diffusion (Nijhawan et al., 2003; Talukdar & Swihart, 2004). These models used the method of moments and the sectional method to provide solutions to the aerosol general dynamic equation (GDE) (Friedlander, 2000; Frenklach & Harris, 1987). Dang & Swihart (2009) expanded upon these by preparing a two-dimensional model for the decomposition of silane driven by laser heating in a tubular reactor. This model included the chemical kinetics and particle processes of the previous model as well as terms for laser heating, fluid dynamics and thermophoresis. When particle models are capable of describing particles with more detail than as simple spheres, a term for sintering is typically included. This provides kinetics governing the gradual coalescence of neighbouring particles into a sphere. Kruis et al. (1993) adapted the expression developed by Kobata et al. (1991) for the grain-boundary diffusion of titania to the sintering of silicon nanoparticles. However, this expression was never validated in a population balance model. Zachariah & Carrier (1999) conducted molecular dynamics (MD) calculations on the sintering of small (up to 2.6 nm) silicon particles. It was reported that sintering was best described by a viscous flow model for temperatures above 830 1C. This model, however, required scaling to reproduce the MD predictions. There is limited work available on numerical modelling of the synthesis of nanoparticles with a narrow particle size distribution (PSD). Tsantilis & Pratsinis (2004) used a moving sectional model to solve the aerosol GDE for the formation of ¨ titania, reporting that accounting for surface reaction leads to narrower PSDs. Kormer et al. presented experimental ¨ ¨ (Kormer et al., 2010a) and theoretical (Kormer et al., 2010b) studies on the synthesis of silicon nanoparticles with a narrow size distribution: this is discussed in more detail in Section 2 as it is central to the model development and analysis of results presented in this paper. ¨ The present work aims to build on the experimental and theoretical results of the studies by Kormer et al. (2010a, 2010b) by further investigating the mechanism of narrowly distributed nanoparticle synthesis and assessing which of the available literature information is most reliable. ¨ The structure of the paper is as follows: Section 2 outlines the previous model of Kormer et al. which is central to the present work. Section 3 describes the philosophy behind multi-scale model development and the mathematical formulation of the models employed. Section 4 presents two models for silicon nanoparticle synthesis and explains the parameter estimation technique used. The results of these models are discussed in Section 5, compared in Section 5.3 and critically analysed in Section 5.4. The paper is concluded by assessing the significance of the models as well as future avenues in which work could be pursued. 2. Previous work ¨ In the experimental work of Kormer et al. (2010a) silicon nanoparticles were produced in the range of 20–40 nm and with a geometric standard deviation (GSD) of 1.06–1.30. It was found that these narrow distributions could only be obtained at high temperatures (900–1100 1C) and low pressures (1–3 mbar silane with 25 mbar total pressure). A brief study of the variation of the process conditions was also conducted in this work. Six of these cases are used in the present work and are summarised in Table 1. The base-case studied has a residence time of 84 ms, a plateau Table 1 ¨ Cases of Kormer et al. (2010a) used in the present work, the base-case highlighted in bold.

t (ms)

T (1C)

pSiH4 (mbar)

/dpri S (nm)

GSD (–)

84 84 84 420 420 420

1100 1100 1100 1100 1000 900

1.0 2.0 3.2 1.0 1.0 1.0

26 27 32 32 29 22

1.07 1.25 1.30 1.07 1.08 1.08

48

W.J. Menz et al. / Journal of Aerosol Science 44 (2012) 46–61

temperature of 1100 1C and a partial pressure of silane of 1.0 mbar. The other cases selected for use here incorporate variations in pressure, temperature and residence time. ¨ The theoretical paper (Kormer et al., 2010b) used a population balance model to study the particle processes leading to the narrowly distributed particles obtained experimentally. The model treated the gas-phase with a single-step homogeneous reaction according to Eq. (1), where a modified form of the kinetic parameters of Petersen & Crofton (2003) was applied SiH4 -Si þ 2H2

ð1Þ

Nucleation of particles from the homogeneous vapour was assumed to occur through a homogeneous nucleation process. This assumed that particles were formed through the collision of monomers or subcritical clusters forming clusters with size larger than a critical cluster diameter dcrit. Particles could subsequently grow through surface reaction (heterogeneous reaction of the precursor monomer with particle surface) or condensation (collision and coalescence of subcritical clusters with a particle). The collision, adhesion and subsequent coalescence of particles is labelled as coagulation. With appropriate scaling of the gas-phase decomposition and surface growth rates, the model showed excellent agreement with experimental data. It was reported that after a short burst of particle nucleation, growth occurs primarily through surface reaction and condensation processes. A key role was attributed to condensation, as numerical calculations could only match experimental observations if condensation was included in the model. However, a one-dimensional particle model was applied, asserting that all particles could be modelled as spheres. The use of a spherical particle model assumes that a collision of two particles creates a new spherical particle with volume equal to the sum of the previous particles. This effectively asserts that sintering occurs instantaneously; an assumption which is invalid where non-spherical particles are obtained, such as at high-pressures. Finally, the model only included a single-step gas-phase reaction, despite evidence of a much more complicated gasphase decomposition process (Ho et al., 1994; Houf et al., 1993; Nijhawan et al., 2003; Petersen & Crofton, 2003; Swihart & Girshick, 1999). The present work assesses these issues by developing two models: one which investigates the sintering particles and another which addresses the gas-phase decomposition and surface growth process. This is achieved by using a state-of-the-art particle model and a robust parameter optimisation methodology. 3. Modelling methodology The multi-scale model described here is characterised by three main components: an application model (AM), an instrumental model (IM) and a data model (DM) (Kraft & Mosbach, 2010). The DM describes the experimental apparatus and settings as well as the raw data obtained. The raw data are then processed into useful information—this comprises the IM. The AM is typically formulated from mathematical equations and aims to reproduce the underlying phenomenon measured by the DM and interpreted by the IM. The inter-relationship of the present work’s AM and the IM ¨ of Kormer et al. (2010a) is displayed in Fig. 1. This paper is primarily concerned with the development of an AM to accurately predict PSDs for the synthesis of silicon nanoparticles with a narrow size distribution. This is achieved through the use of an operator-splitting code (MOPS) which solves the gas-phase chemistry and particle population balance while periodically exchanging information between them. MOPS has already been applied to model the synthesis of silica (Grosschmidt et al., 2002; Shekar et al., in press) and titania (Lavvas et al., 2011; West et al., 2009) nanoparticles.

Fig. 1. Schematic of the model used in the present work.

W.J. Menz et al. / Journal of Aerosol Science 44 (2012) 46–61

49

The kinetic mechanism is solved by use of an in-house-developed software Sprog to solve the system of ordinary differential equations (ODEs) describing the gas-phase chemistry (Celnik, 2007). It calculates the chemical reaction rates and the composition of the gas-phase mixture, modelling the rate constants as Arrhenius processes. The particle model (Sweep) includes terms for inception, surface growth, coagulation and sintering; describing primary particles by their volume, connectivity and surface area. The properties of a particle Pi may be completely described by the variables of its type space. This model’s type space is represented by a finite dimensional vector, given by Pi ¼ P i ðp1 , . . . ,pn ,C,I,SÞ

ð2Þ

where i 2 f1, . . . ,Ng and N is the total number of particles in the system. Each particle Pi consists of n primary particles pj (j 2 f1, . . . ,ng), where N and n vary with time. C, I and S are matrices which describe the sintering between primaries (Sander et al., 2009). The primaries are described by their volume vj pj ¼ pj ðvj Þ

ð3Þ

and the sum of all primary volumes is equal to the total volume V of the particle V¼

n X

vj

ð4Þ

j¼1

Thus, particles are modelled by tracking the volume of every primary particle and the common surface between any two. A variety of sub-models are also included in the particle model which describe particle formation, growth and sintering. These are briefly outlined in the following paragraphs, however a detailed description is provided by Celnik (2007) and Sander et al. (2009). Inception: An inception event is modelled as a collision of two silicon hydrides and is dependent on their gas-phase concentration and the transition regime coagulation kernel (Lavvas et al., 2011; Pratsinis, 1988). Inception events increase the number of particles N in the type space by 1 moleculeþ molecule-P N þ 1 ðp1 ,C,I,SÞ

ð5Þ

The rate of inception is calculated using the transition regime kernel and the concentration of colliding species Rinception ¼ 12Zincep K tr N 2A C Sii Hj C Sik Hl

ð6Þ

where Zincep is the inception collision efficiency, K tr is the transition regime coagulation kernel, NA is Avogadro’s Number, and C Sii Hj and C Sik Hl are the concentrations of colliding silicon hydride species Sii Hj and Sik Hl respectively. The form of the transition kernel is given in Patterson et al. (2006a) and is dependent on the value of the Knudsen number. Surface growth: Sweep can consider surface growth as a condensation process (requiring the gas-phase species’ diameter) or a surface reaction process (requiring modified Arrhenius parameters). Surface growth alters the type space by increasing the volume of a selected primary pj of particle Pi by increment dv pj ðvÞ-pj ðv þ dvÞ Condensation is modelled as a collision process, with rate given by Celnik (2007) sffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi pkB T Rcond ¼ 2:2Zcond C Sii Hj þ dcol Þ ðd 2mSii Hj Sii Hj

ð7Þ

ð8Þ

where Zcond is the collision efficiency, kB is Boltzmann’s constant; mSii Hj & dSii Hj are the mass and diameter of the colliding species; and dcol is the collision diameter of the particle. In this work, the rate of a surface reaction process is proportional to the rate constant kSR and the surface area Sj of particle j RSR ¼ kSR Sj

I Y

C vi i

ð9Þ

i¼1

where vi is the forward stoichiometric coefficient of species i. The rate constant is given by typical modified-Arrhenius form   EA,SR ð10Þ kSR ¼ ASR T n,SR exp RT Coagulation: Coagulation is modelled as the collision and subsequent adhesion of two primary particles. Storing the common surface information in the C, I and S matrices allows for the connectivity of primaries to be tracked, and thus the collision diameter of the aggregate particle dcol to be calculated (Shekar et al., in press). A coagulation event may be represented by the following change in the type space: Pi þP j -Pk ðp1 , . . . ,pnðPi Þ ,pnðPi Þ þ 1 , . . . ,pnðPk Þ ,C,I,SÞ

ð11Þ

50

W.J. Menz et al. / Journal of Aerosol Science 44 (2012) 46–61

Coagulation processes are described by the general Smoluchowski coagulation equation (Celnik, 2007), again using the transition-regime coagulation kernel Rcoagulation ¼

x1 1 X 1X K tr ðxy,yÞnðxyÞnðyÞ K tr ðx,yÞnðxÞnðyÞ 2y¼1 y¼1

ð12Þ

where n(x) and n(y) are the numbers of particles of type x and y. Sintering: The exponential decay formula popularised by Koch & Friedlander (1990) is adapted in the present work to locally calculate the level of sintering between two neighbouring primary particles d½C i,j  1 ¼  ð½C i:j ½Si,j Þ tS dt

ð13Þ

where ½C i,j  and ½Si,j  are the i, jth elements of the sintering matrices C and S; and tS is the characteristic sintering time. The characteristic sintering time is usually a function of temperature and particle diameter and is discussed in Section 4.1. ¨ It is also important to highlight the difference in nomenclature between the present model and the model of Kormer et al. (2010b). Here, condensation is assumed to be the collision and coalescence of gas-phase species with a particle. ¨ In Kormer et al. (2010b), a distinction is made between collisions of gas-phase species and subcritical clusters with a particle: the former is labelled ‘surface growth’, while the latter ‘condensation’. This distinction cannot be made for the present model as the type space does not consider the formation of subcritical silicon clusters. The interaction of these processes with the ensemble of stochastic particles is implemented with a Linear Process Deferment Algorithm (Patterson et al., 2006). The particle population balance solver is coupled with a kinetic mechanism through the operator-splitting method described by Celnik (2007) and Celnik et al. (2007, 2009). 3.1. Numerical processes In order to accurately model a physical system, an AM typically requires a mathematical description of key processes and input parameters to generate results. The numerical processes and key parameters used in the present work are described in this section, where the latter are tabulated in Table 2. The operator-splitting code which couples the gas-phase chemistry and population balance is described in detail by Celnik (2007) and Celnik et al. (2007, 2009). Strang splitting is applied to reduce the splitting error, where the stochastic population balance equations are solved for one step followed by the deterministic gas-phase for one more step (Celnik et al., 2007). The timesteps Dt and number of splits Nsplits must be defined for MOPS, the latter describing the number of times which information is exchanged between gas- and particle-phase. MOPS also takes input of the basic process conditions such as temperature (T), component partial pressure (pSiH4 , pAr) and residence time (t). The conditions used here are those of ¨ Kormer et al. (2010a). A temperature gradient is also applied across the reactor, displayed in Fig. 2. At longer residence

Table 2 Summary of input parameters to MOPS. Parameter

Symbol

Value

Ref.

Nsplits

30 1.0  10  6 s 1  10  18 1  10  4 16,384 1 1  10  10

– – – – – – –

pSiH4 pAr

750 1C 0.08–0.42 s 1.0–3.2 mbar 22–24 mbar

¨ Kormer ¨ Kormer ¨ Kormer ¨ Kormer

Model parameters Kinetic Arrhenius parameters Surface reaction parameters Sintering parameters Inception and collision efficiency

AGP , nGP , EA,GP ASR , EA,SR AS ,ES ,dp,min Zincep , Zcond

Various This work This work 1.0

Swihart and Girshick (1999a) – – –

Diameter of species SiiHj

dSii Hj

Various

–

Mass of species SiiHj

mSii Hj

Various

–

Molecular weight of silicon Bulk density of silicon

Mw

28.1 g/mol 2.33 g/cm3

Perry et al. (1997) Perry et al. (1997)

Numerical parameters Number of splits Timestep Absolute error tolerance (Sweep) Relative error tolerance (Sweep) Number of stochastic particles Number of runs Maximum zeroth moment Process settings Initial temperature Residence time SiH4 partial pressure Ar partial pressure

Dt

ea er N SP L M0,max T

t

rSi

et et et et

al. al. al. al.

(2010a) (2010a) (2010a) (2010a)

W.J. Menz et al. / Journal of Aerosol Science 44 (2012) 46–61

1200

51

o 1100 C

Temperature,°C

1100

1000

900

800

700

600 10-6

10-5

10-4

10-3

10-2

10-1

Time, s Fig. 2. Temperature profile applied across the reactor for 0.08 s residence time.

times or lower temperatures, the profile is ‘stretched’ or ‘compressed’, maintaining the same initial and final temperatures ¨ and lineshape. This profile is identical to that used in Kormer et al. (2010b). Sprog requires input of the chemical mechanism and thermodynamic polynomials describing the silicon hydride species’ thermodynamic properties as a function of temperature. For this work, the thermodynamic polynomials and Arrhenius kinetic parameters were primarily sourced from the work of Swihart & Girshick (1999a). Absolute and relative error tolerances (ea and er respectively) were chosen to ensure the ODE system was converged. The particle population balance is solved by a stochastic particle algorithm first used in the modelling of formation of soot nanoparticles (Goodson & Kraft, 2002). It has since been further refined by Patterson et al. to incorporate a wider variety of processes and more efficient numerical methodology (Patterson & Kraft, 2007; Patterson et al., 2006a, 2006b). Two numerical parameters must be defined for proper operation of the stochastic algorithm: the number of runs L and stochastic particles NSP , which were chosen sufficiently high such that increasing them does not alter the ultimate PSD. Sweep additionally requires input of the molecular weight (M w ) and bulk density of silicon (rSi ) to calculate the volume added into the particle-phase with each inception or surface growth step. The mass and diameter of each silicon hydride species is required for evaluation of the coagulation kernel—these were obtained from the quantum mechanical calculations of Swihart & Girshick (1999b) and estimated from the structures reported in Swihart & Girshick (1999a). ¨ The AM was fitted to the IM of Kormer et al. (2010a) by adjusting the sintering or gas-phase/surface growth parameters using a detailed optimisation methodology, discussed in more detail in Section 4.3. 4. Model formulation The majority of the previous work surrounding the gas-phase and film growth mechanisms of silicon from silane have been conducted at atmospheric pressure and temperatures below 900 1C (Ho et al., 1994; Houf et al., 1993; Hogness et al., 1936; Kruis et al., 1994; Nijhawan et al., 2003; Purnell & Walsh, 1966; Swihart & Girshick, 1999). While expressions exist for extending the gas-phase kinetics to low-pressure regions (Koshi et al., 1991; Petersen & Crofton, 2003), it was shown in ¨ earlier work that typical kinetics still require scaling to match experimental predictions (Kormer et al., 2010). Thus, the ¨ conditions Kormer et al. (2010a) are poorly characterised in terms of the gas-phase. Furthermore, as discussed in Section 1, the previous formulations of the sintering kinetics (Kruis et al., 1993; Zachariah & Carrier, 1999) for silicon nanoparticles have never been verified in a particle population balance. For developing a new model for narrowly distributed silicon nanoparticle synthesis, this raises the question of which pieces of information may be translated to these process conditions. The present work proposes two models: one which assumes the gas-phase chemistry is correct and adjusts the sintering parameters (Model A); and another which adjusts the gas-phase, neglecting coagulation and sintering (Model B). 4.1. Model A In Model A, the popular kinetics of Swihart & Girshick (1999a) are chosen for the gas-phase chemistry. These kinetics were chosen as they include a detailed description of the decomposition of silane into three-dimensional polycyclic silicon clusters. The present work uses a reduced version of these kinetics, incorporating 131 silicon hydride species (up to Si13 ) and 1098 chemical reactions. The population balance modelling works using this mechanism typically assume a constant critical nucleus number N Si,min of 10 silicon atoms (Dang & Swihart, 2009; Nijhawan et al., 2003). However, this work will assume no such knowledge about the critical nucleus.

52

W.J. Menz et al. / Journal of Aerosol Science 44 (2012) 46–61

With so many gas-phase species, we must consider how to best model the particle inception and surface growth processes. Here, a general inception reaction of the following form was considered: Sii Hj þ Sik Hl -½Sii þ k  þ

jþl H2 2

ð14Þ

where the left-hand species Sii Hj and Sik Hl are in the gas-phase and the right hand species is a newly incepted particle (square brackets denoting a particle). Then, a critical nucleus number NSi,min was chosen; the minimum number of silicon atoms present in the smallest particle incepted. NSi,min falls in the range of 2–26 due to the number of silicon atoms present in the smallest and largest species of the gas-phase mechanism. All combinations of i and k whose sum was greater than or equal to NSi,min were subsequently assumed to form particles in a collision event. The model also requires treatment of surface growth. While experimental parameters exist for surface reaction (Ho et al., 1994), these only extend as far as Si3 species. Thus, surface growth was modelled as a condensation process (i.e. collisional), with all species permitted to condense. As Model A assumes the gas-phase chemistry is correct, the sintering functional must be chosen. The form of tS varies depending on the sintering mechanism assumed. For a grain-boundary diffusion (GBD) mechanism, there is a power of four dependence on particle diameter (Kobata et al., 1991)    d E ts ¼ As d4i,j T exp s 1 p,min ð15Þ RT di,j where As is the sintering pre-exponential, Es is the characteristic energy and di,j is the minimum diameter of two adjacent primary particles. dp,min is the diameter below which the particles sinter instantaneously, originally suggested by Tsantilis et al. (2001). The GBD kinetic is employed here at the recommendations of the work of Kruis et al. (1993). 4.2. Model B ¨ Model B follows a similar modelling philosophy to that employed by Kormer et al. (2010b). Here, a single-step mechanism is assumed for simplicity, given by the decomposition of silane into silylene and hydrogen kSiH4

SiH4 þM - SiH2 þ H2 þ M where ‘M’ is a collision partner. This reaction has a rate constant kSiH4 given by   E kSiH4 ¼ AGP T nGP exp A,GP RT

ð16Þ

ð17Þ

¨ An activation energy of 45 kcal/mol is used, consistent with the models of Kormer et al. (2010b) and Woiki et al. (1997). The rate of decomposition of silane attributed to only the gas-phase is therefore written as  dC SiH4  ¼ kSiH4 C SiH4 C M ð18Þ dt GP The SiH2 produced from this reaction is allowed to collide to incept particles into the population balance. The gas-phase decomposition rate AGP primarily effects the nucleation rate and is therefore one of Model B’s adjustable parameters. Surface growth must now also be addressed. Ho et al. (1994) and Houf et al. (1993) presented a mechanism for the reaction of silicon hydrides on a silicon surface. They proposed that silanes must proceed through an activation barrier before they could adsorb to the surface, while the interaction of silylenes and silenes with the surface was described by a collisional process. However, the parameters quoted in Ho et al. (1994) and Houf et al. (1993) are incompatible with the type space of the present work, and the gas-surface mechanism proposed in those works shows considerable divergence from experimental results at temperatures above 800 1C. A similar observation is made for the surface growth parameters of Peev et al. (1991), which diverge above 960 1C. Thus, Model B also adjusts the surface growth rate. Silane is assumed to react through an Arrhenius process, with rate given by Eq. (9) and rate constant by Eq. (10). The surface reaction rate for silane is adjusted with the parameters ASR,SiH4 , assuming an activation energy of 37.5 kcal/mol (Ho et al., 1994; Houf et al., 1993). The reaction of silylene (SiH2) with the surface is modelled as a barrierless process (EA,SR,SiH2 ¼ 0) with an adjustable pre-exponential factor ASR,SiH2 . The total rate of surface reaction on particle j for Model B is therefore given by RSR ¼ Sj ðkSR,SiH4 C SiH4 þ kSR,SiH2 C SiH2 Þ

ð19Þ

¨ Kormer et al. (2010b) observed that the best model correspondence was obtained when coagulation was neglected, due to the such low number concentration of particles. Model B follows in suit with this, and no coagulation process is included. The implications and limitations of this assumption are thoroughly discussed in Section 5.4. Sintering is therefore not accounted for in this model. Finally, Model B assumes that only 25% of the silane reacts to form particles. This ¨ is done to account for the growth of silicon on the walls of the reactor and is the same proportion as assumed in Kormer et al. (2010b).

W.J. Menz et al. / Journal of Aerosol Science 44 (2012) 46–61

53

Table 3 Comparison of the gas-phase and particle processes describing Models A and B. Process

Model A

Model B

Gas kinetics Inception

Swihart & Girshick’s mechanism (variable N Si,min )

SiH4 -SiH2 þ H2 (variable AGP , EA,GP ¼ 45:0 kcal=mol) SiH2 þ SiH2 -½Si2  þ 2H2 (constant Zincep ¼ 1:0)

Surface growth

j Sii Hj þ ½Six -½Six þ j  þ H2 (condensation with constant 2 Zcond ¼ 1:0)

l Sii Hj þ Sik Hl -½Sii þ k  þ j þ 2 H2 (constant Zincep ¼ 1:0)

Coagulation Sintering

Transition kernel GBD sintering (variable AS , ES , dp,min )

SiH4 þ ½Six -½Six þ 1  þ 2H2 (surface reaction with variable ASR,SiH4 , EA,SR,SiH4 ¼ 37:5 kcal=mol) SiH2 þ ½Six -½Six þ 1  þ H2 (surface reaction with variable ASR,SiH2 , EA,SR,SiH2 ¼ 0 kcal=mol) None Not applicable

Thus, Models A and B represent two different pathways in which the various sub-models describing nanoparticle synthesis can be assembled. The process of which they are composed are directly compared in Table 3. 4.3. Estimation of parameters Choice of appropriate model parameters is of key importance in obtaining theoretical results comparable to the experimental results which they seek to reproduce. For Model A, the sintering parameters are adjusted. This system’s parameter vector is therefore given by xA ¼ ðAs ,Es ,dp,min Þ

ð20Þ

Model B chooses instead to adjust the pre-exponential factors for the gas-phase chemistry and surface growth xB ¼ ðAGP ,ASR,SiH4 ,ASR,SiH2 Þ

ð21Þ

In order to find the optimal set of parameters for either model, a three-tiered strategy is applied. In the first step, we locate a region of parameters with approximate correspondence to experimental results by following the approach of Smallbone et al. (2011) and Shekar et al. (in press). Here, the parameter space is scanned with a low-discrepancy (Sobol) sequence. The four best parameter sets from this analysis are chosen, based on the objective function FðxÞ

FðxÞ ¼

N exp X

ð½/dpri Sexp /dpri Ssim ðxÞ2 þ ½Varðdpri Þexp Varðdpri Þsim ðxÞ2 Þ i i i i

ð22Þ

i¼1

¨ et al. (2010a), /dpri Sexp , Varðdpri Þexp , /dpri Ssim and where Nexp is the number of experimental cases from Kormer i i i sim Varðdpri Þi are the experimental and simulated mean primary diameter and variance of the primary diameter respectively. In the second step, further refinement towards the optimal set of sintering parameters is achieved by applying the simultaneous stochastic approximation algorithm (SPSA) (Hirokami et al., 2006), leading to the refined set of parameters xn that satisfy xn ¼ argminðFðxÞÞ

ð23Þ

x

In the final, third step, the optimal set of parameters is obtained through the Bayesian approach used by Mosbach et al. (accepted for publication). As before, the mean primary diameter and the variance of the primary diameter distribution are used as responses. An important ingredient for the Bayesian analysis is the choice of the prior distributions, and since the uncertainties on the responses are unknown, an Inverse-Wishart (Box & Tiao, 1974; Mosbach et al., 2011) prior is used for this. As in previous studies (Braumann et al., 2011; Mosbach et al., accepted for publication), a Metropolis–Hastings algorithm (Hastings, 1970; Metropolis et al., 1953) is employed for the estimation of the posterior distribution of the parameters. Since the algorithm requires many model evaluations, it is computationally prohibitive to use the detailed population balance model. Instead, we approximate the model responses locally around the point xn by surrogate models. Here, second order response surfaces are used, which are cheaper to evaluate. This analysis provides an estimate of the uncertainty of the model parameters. 5. Results 5.1. Model A To verify the integrity of Sprog in solving this expansive gas-phase mechanism, the kinetic model was first solved at the process conditions of Swihart & Girshick (1999a). The evolution of the concentration profiles was identical, so the system ¨ was subsequently solved at the baseline conditions of Kormer et al. (2010a) without the particle model coupled. The results of this calculation are presented in Fig. 3.

54

W.J. Menz et al. / Journal of Aerosol Science 44 (2012) 46–61

-8

SiHX

Si11+HX

Concentration, mol/cm

3

10

10

-10

10

-12

10

-14

10

-16

Si6HX

Si7HX

Si2HX

Si8HX Si3HX

Si4HX 10

-6

10

-5

Si9HX Si10HX

Si5HX 10

-4

10

-3

10

-2

10

-1

Time, s Fig. 3. Temporal evolution of the concentration of silicon hydride species according to the mechanism of Swihart and Girshick (1999a).

Table 4 Optimal parameters obtained from the estimation procedure for Model A. The upper/ lower bounds represent the 99% confidence interval. Parameter AS ES dp,min

Units 4

s/m K kJ/mol nm

Lower bound 1.5  10 440 4.8

12

Mode

Upper bound 12

9.4  10 450 6.5

1.1  1013 480 6.5

It is evident here that with only the gas-phase mechanism considered, one-silicon species (SiH4 and SiH2) are still present in large concentrations, even as the residence time is reached: this suggests that the reaction does not go to completion in this model. The concentration of ‘nuclei’ (arbitrarily labelled as species with eleven silicon atoms or more in Swihart & Girshick, 1999a) gradually increases as the higher species are irreversibly formed. In a fully coupled model, these nuclei would subsequently react to form particles. The particle model was then coupled to the gas-phase kinetics, and the parameter estimation procedure conducted. The optimal set of sintering parameters given in Table 4 was the result of fitting Model A to the six cases described in Table 1. The sintering rate was only weakly dependent on AS , as illustrated by the broad confidence bounds. The characteristic energy found through this analysis is a large departure from the unvalidated energies reported by Kruis et al. (1993) and Zachariah & Carrier (1999). When varying NSi,min, the primary particle size ranged from 26 to 40 nm at the base case, using the modal parameters of Table 4. The particles obtained for different values of NSi,min varied in characteristics; from highly spherical, narrowly distributed particles similar to those observed experimentally (NSi,min of 2–5) to broadly distributed and sintered agglomerates (NSi,min of 6–16). The latter phenomena are represented by the increase in primary diameter and variance, coupled with a decrease in collision diameter. Fig. 4 clearly demonstrates that the choice of critical nucleus number can strongly affect the PSD obtained. Despite this, ¨ there are multiple values of NSi,min which match the experimentally obtained (Kormer et al., 2010) primary diameter and variance. This raises the issue of model discrimination: given that several acceptable models have been found, how should an appropriate inception mechanism be chosen? To answer this question, we employ homogeneous nucleation theory to estimate a value of NSi,min which is appropriate, ¨ given the process conditions. Using the equations presented by Kormer et al. (2010b), the critical diameter may be calculated at the initial conditions of the reactor based-on a macroscopic properties of silicon particles. Here, measurements of the surface energy y are used in conjunction with the Kelvin Equation (Kodas & Hampden-Smith, 1999) to find the critical diameter. The smallest particle nucleus considered in Model A is a Si2 nucleus, for NSi,min ¼ 2. The corresponding gas-phase species to this is silene (Si2H4A). A quantum mechanical calculation at the B971/3-21G level in Gaussian 03 was carried-out on silene to estimate its molecular volume, from which the isotropic collision diameter (Halpern & Glendening, 1996) was determined. The results of these analyses are presented in Table 5. Under the assumptions made here, the critical diameter is less than the diameter of the smallest precursor nucleus. This corresponds to the collision-limited regime of homogeneous nucleation (Kodas & Hampden-Smith, 1999), where growth of clusters is not inhibited by a thermodynamic barrier. It is therefore reasonable to choose NSi,min ¼ 2 for the critical nucleus number, as the smallest nucleus number must be at least two for homogeneous nucleation (Kruis et al., 1994). This is consistent with the results of Wu et al. (1987), who noted that the critical nucleus size was of atomic dimension based on

W.J. Menz et al. / Journal of Aerosol Science 44 (2012) 46–61

55

Diameter measurement, nm

140 〈dcol〉

120 100 80 60 40

〈dpri〉

Exp.

σ(dpri)

Exp.

20 0 0

5

10

15 NSi,min, -

20

25

30

Fig. 4. Comparison of the mean and standard deviation of the primary diameter and the mean collision diameter at baseline process conditions. ¨ The experimental data of Kormer et al. (2010a) are given by the horizontal lines. Table 5 Comparison of critical diameters and surface energies. Arrows indicate the direction of calculation. Reference

NSi,min

Interpolation from Mezey and Giber (1982)

–

This work (Si2H4A)

2

g (N/m)

dcrit (nm) 0.20 QM

-

0.59

HNT

’

0.99 –

Table 6 Optimal parameters obtained from the estimation procedure for Model B. The upper/lower bounds represent the 99% confidence interval. Parameter

Units

Lower bound

Mode

Upper bound

AGP ASR,SiH4 ASR,SiH2

cm3/mol s 1/cm2 mol s 1/cm2 mol s

3.1  1010 5.7  1032 1.3  1029

5.1  1010 6.9  1032 3.0  1029

6.1  1010 9.1  1032 3.7  1029

macroscopic properties of silicon. Choosing a critical nucleus number of two also saves computational expense, allowing for cheaper model evaluations in the parameter estimation procedure. The response for Model A and the underlying mechanism by which it predicts particles to be formed are discussed in Section 5.3. 5.2. Model B The optimal parameters from fitting to the six experimental cases (Table 1) are given in Table 6. The gas-phase preexponential (AGP ) is four orders of magnitude lower (cf.7  1014) than that predicted by the mechanism of Ho et al. (1994) and Swihart & Girshick (1999a) when including the Tn term at 1100 1C. Indeed, it is also two orders of magnitude smaller ¨ (cf. 9  1012) than reported by Kormer et al. (2010b) when fitting to only the base case. The much slower gas-phase decomposition kinetic found through the optimisation method is necessary to suppress ¨ nucleation and promote the growth of silicon on the particles’ surface. It is lower than that found in Kormer et al. (2010b) because of fitting to the 420 ms residence time cases. In these cases, an even lower nucleation rate is needed in order to allow for more precursor to grow on the surface, thus creating larger particles. It also suggests that under these process conditions, silicon is removed from the gas-phase almost exclusively by the surface growth processes. The magnitude of the parameters ASR,SiH4 and ASR,SiH2 cannot be evaluated against previous silicon surface growth parameters, as there exists no literature data to which they may be directly compared. 5.3. Comparison of models With the modal parameters of Tables 4 and 6, both models showed some agreement with experimental PSDs, as illustrated in Fig. 5. Model A shows excellent correspondence for the base case, and good modal agreement for the 420 ms residence time cases. However, there PSDs predicted at increased-pressure conditions were far too narrow that those

56

W.J. Menz et al. / Journal of Aerosol Science 44 (2012) 46–61

o

o

T=1100 C, τ=420 ms, pSiH4=1.0 mbar

0.4

0.4

0.3

0.3 qo, 1/nm

qo, 1/nm

Base: T=1100 C, τ=84 ms, pSiH4=1.0 mbar

0.2 0.1 0 15

0.2 0.1

20

25 dpri, nm

30

0 20

35

0.4

0.3

0.3 qo, 1/nm

qo, 1/nm

0.4

0.2 0.1

40

0.1

20

30 dpri, nm

40

0 20

50

o

25

30 dpri, nm

35

40

o

T=900 C, τ=420 ms, pSiH4=1.0 mbar

0.4

0.4

0.3

0.3 qo, 1/nm

qo, 1/nm

35

0.2

T=1100 C, τ=84 ms, pSiH4=3.2 mbar

0.2 0.1 0 10

30 dpri, nm

T=1000 oC, τ=420 ms, pSiH4=1.0 mbar

T=1100 oC, τ=84 ms, pSiH4=2.0 mbar

0 10

25

0.2 0.1

20

30 40 dpri, nm

50

60

0 15

20

25

30

dpri, nm

¨ Fig. 5. PSDs for variation of process conditions. The experimental data of Kormer et al. (2010a) are given by crosses, and the optimal solutions for Models A and B given by the dashed lines and solid lines respectively.

¨ observed by Kormer et al. (2010a). This is also evident (to a lesser degree) for the width of the 420 ms cases. The cause of this phenomenon is explained in the following paragraphs, which examine the particle rates and PSD evolution in detail. The base case and PSD predicted by Model B is slightly too broad, with a calculated GSD of 1.10 as compared to 1.07 for the experimental data. As already discussed, this is partially a result of weighting the objective function equally to all cases: better fit could have been obtained, but would come at the cost of correspondence with other PSDs. The 420 ms cases show excellent agreement with the experimental results, with the exception of the 1100 1C which is slightly shifted to lower diameters. The cases at increased pressure are too narrow under Model B, but not to the extent as observed in Model A. The 3.2 mbar case shows very good agreement with the modal value of the experimental PSD. However, it is unlikely that the broadening observed in the experimental data can solely attributed to the inception or surface reaction process. That is, with increased pressure, coagulation could become important. Neglecting coagulation for these cases may therefore not

W.J. Menz et al. / Journal of Aerosol Science 44 (2012) 46–61

57

28

10

26

28

10

Inception Surface growth Coagulation

1024 10

22

10

20

10

18

10

16

Inception SiH4 Surface reaction SiH2 Surface reaction

1026 3

10

Particle rate, 1/cm s

Particle rate, 1/cm3s

be a reasonable assumption, and including coagulation with a sintering functional could allow for better model correspondence. This is further discussed in Section 5.4. The mechanism of particle growth for either model can be elucidated by comparison of the rates of the particle processes described in Table 3. This is presented in Fig. 6. For Model A (Fig. 6a), it is evident that nucleation and surface growth occur very quickly, allowing coagulation to control the the characteristics of the distribution. Typically, this would lead to a self-preserving PSD, however the adjusted sintering kinetics predict that aggregates with very narrowly distributed particles are formed. The rapid particle inception observed in this figure is in agreement with the results of Kruis et al. (1994), who report short time-lags before nucleation for the decomposition of silane at high temperatures. In contrast to the rapid nucleation and surface growth observed in Model A, Fig. 6b demonstrates that the particle processes occur at a much later time (after 10  2 s) in Model B. It may be observed here that inception and nucleation follow the increase in temperature imposed by the temperature profile (Fig. 2). After a brief burst of inception, surface ¨ reaction dominates, causing depletion of the precursor. This growth mechanism is identical to that proposed by Kormer et al. (2010b). The temporal evolution of particles in the reactor may also be analysed by investigation of the primary PSD throughout the simulation. Fig. 7a illustrates how the primary PSD changes with time for Model A. Here, it is evident that primaries are rapidly incepted and grown to several nanometres through condensation. Coagulation of primaries subsequently begins to broaden the distribution, after which the smaller primaries rapidly sinter to form large primaries, as shown by the shoulder in the 8 ms curve. Continuous coagulation and sintering narrow the primary PSD and causes aggregation of the primaries. This is typical of a collision-limited nucleation process, which exhibits a rapid sintering rate relative to the collision rate, causing a high degree of sphericity in primary particles (Kodas & Hampden-Smith, 1999). The rapid cooling at towards the end of the reactor massively increases the characteristic sintering time, preventing the spherical primaries within an aggregate from sintering into irregularly shaped particles.

24

10

1022 1020 18

10

16

10

-6

10

-5

10

-4

-3

-2

10

10

10

10

-1

10

-6

-5

10

10

-4

10

-3

-2

-1

10

10

Time, s

Time, s

Fig. 6. The temporal evolution of simulated particle rates for (a) Model A and (b) Model B at the base case.

0.25

0.35

84 ms 8 ms 0.8 ms 0.08 ms

0.3

84 ms 35 ms 25 ms 20 ms

0.2

q0, 1/nm

q0, 1/nm

0.25 0.2 0.15

0.15

0.1

0.1 0.05 0.05 0

0 0

5

10

15 20 dpri, nm

25

30

35

0

5

10

15 20 dpri, nm

Fig. 7. Temporal evolution of the primary PSD for (a) Model A and (b) Model B at the base case.

25

30

35

58

W.J. Menz et al. / Journal of Aerosol Science 44 (2012) 46–61

The mechanism of forming narrowly distributed particles by coagulation and sintering was also identified for the formation of titania particles by Heine & Pratsinis (2007), where primary PSDs narrower than self-preserving (GSD¼1.4) were obtained by accounting for localised sintering within aggregate particles. This highlights that there are at least two potential mechanisms for growth of narrowly distributed particles: this coagulation-sintering mechanism, and the ¨ surface-growth mechanism recognised by Kormer et al. (2010b) and Tsantilis & Pratsinis (2004). The rapid nucleation shown in Figs. 6a and 7a is observed for all cases, with the entire precursor typically being consumed by 0.05 ms. Since inception occurs so rapidly, only the number concentration is affected by pressure, to which the adjusted sintering kinetic is insensitive to. This causes the overly narrow PSDs observed in Fig. 5. As suggested by the comparison of particle rates, a completely different growth mechanism is observed in Model B. In Fig. 7b, we see the ‘seed’ particles being incepted at a much later time (20 ms as compared to o 0:08 ms for Model A). The surface growth process subsequently dominates after this: the PSD uniformly grows—as coagulation is neglected—shifting to larger particle sizes with no change in distribution width. ¨ The models may also be evaluated by considering the number concentration (M0). In the work of Kormer et al. (2010b), 0 a number concentration ‘in the range of 109 #=cm3 was reported for the base case. In this work, Model A predicts a number concentration of 3.0  109 #/cm3, and Model B a concentration of 5.0  109 #/cm3. These are both consistent with experimental and theoretical predictions of previous work. It is interesting to note that despite assuming 100% of the silane reacts to form particles, Model A still produces a number concentration with reasonable correspondence to experimental results. The number concentration quoted for Model A represents the concentration of aggregates of highly spherical primaries, rather than of the individual ‘subparticles’ themselves. There are therefore two factors affecting M0 for Model A, acting in opposite directions. 5.4. Critical analysis of models The results of these two models indicate that from a purely theoretical point of view, it is possible to synthesise narrowly distributed nanoparticles through two completely different mechanisms. However, it is important to understand whether these parameters exist as reasonable physical quantities by comparing them with available experimental data. Model A shows generally good agreement for the base and 420 ms cases, but it is unable to account for variations in silane partial pressure. It is unlikely that the sintering parameters should include a pressure dependency, which suggests that there is some fundamental process in this model which is not being adequately represented. Furthermore, if the characteristic sintering time predicted by Table 4’s parameters at 1100 1C (the base-case temperature) are compared with previous literature kinetics, we observe a difference of many orders of magnitude and a much stronger dependency on particle diameter (Fig. 8). These sintering times are clearly unphysical and are in contradiction to the MD calculations of Zachariah & Carrier (1999). It is therefore proposed that the kinetics of Swihart & Girshick (1999) and Ho et al. (1994) are invalid at the process conditions of this study. As such, the ‘error’ in these kinetics is reflected in unphysical sintering rates. While the kinetics surrounding the gas-phase decomposition of silane are fairly well established (Petersen & Crofton, 2003; Swihart & Girshick, 1999a), the results of Model B could suggest that under some specific process conditions, a fundamental change in the decomposition process occurs. This would be associated with heterogeneous growth becoming favoured over the gas-phase reaction. The accuracy of this statement could be assessed by conducting quantum mechanical calculations on the decomposition and heterogeneous group process. Model B could therefore provide insight into the appropriate order of magnitude of the kinetic parameters for these process conditions. Its better quantitative and 10

2

100 10-2

τS , s

10-4 10-6 -8

10

-10

10

-12

10

This work Kruis et al., 1993 Zachariah & Carrier, 1999

-14

10

-16

10

0

5

10

15

20

25 30 dpri, nm

35

40

45

50

Fig. 8. Comparison of the characteristic sintering time at 1100 1C for the optimal parameters in Table 4, the GBD model of Kruis et al. (1993) and the solid-state diffusion model of Zachariah & Carrier (1999).

W.J. Menz et al. / Journal of Aerosol Science 44 (2012) 46–61

o

o T=1100 C, τ=84 ms, pSiH4=1.0 mbar

0.25

0.2

0.2

0.15

0.15

qo, 1/nm

qo, 1/nm

0.25

0.1 0.05 0 15

59

T=1100 C, τ=84 ms, pSiH4=3.2 mbar coag. on coag. off exp.

0.1 0.05 0

20

25

30

dpri, nm

35

40

10

20

30

40

50

60

dpri, nm

Fig. 9. PSDs generated by Model B with coagulation turned-off (solid line) and -on (dashed line) for the base and 3.2 mbar cases. The experimental data are shown by crosses.

qualitative agreement with experimental results indicates that the mechanism predicted by this model is representative of the real physical processes occurring in the centreline of a tubular reactor. ¨ The results of Model B confirm the assertions of Kormer et al. (2010b), Koshi et al. (1991) and Wu et al. (1987) that the decomposition of silane can in some cases be approximated as a single-step process. It is possible, however, that better fit could be obtained by including higher silanes in the gas-phase mechanism and limiting the nucleation rate based-on a dynamic calculation of the critical nucleus size. Again, this would require greater investigation into the gas-phase kinetics via quantum-mechanical calculations and is outside the scope of the present study. As raised in Section 5.3, the validity of the no-coagulation assumption may be questionable for some cases. Fig. 9 shows how the PSDs predicted by Model B are affected by coagulation, using the transition kernel and the GBD sintering kinetic of Kruis et al. (1993). It is shown here that the base case remains relatively unchanged, with GSD increasing from 1.10 to 1.16. This emphasises that the assumption of no-coagulation is reliable for the base case process conditions. Further, when using the sintering kinetic of Kruis et al. (1993), Fig. 8 indicates that sintering should occur very rapidly 1373 K. This causes any particle which coagulate to coalesce into spheres almost instantaneously, implying that the use of a spherical particle ¨ model (in the manner of Kormer et al., 2010b, for example) is also sound at these process conditions. At 3.2 mbar, the GSD increases from 1.10 with no coagulation to 1.20 with coagulation and sintering turned-on. The PSD with coagulation included shows much better correspondence with the experimental data, which reported a GSD of 1.30. It is evident here that coagulation is necessary to achieve the broadening observed in Figs. 5 and 9 for the increasedpressure cases. Thus, this system clearly sits on the threshold at which coagulation of particles becomes an important process. Finally, the issue of the loss of silicon to reactor walls must be addressed. No set of parameters for Model A with appropriate fit to experimental PSDs could be found when assuming a reduced initial concentration of silane in a similar way to Model B. The fact that this model is unable to account for the significant proportion of precursor lost to the walls further questions the validity of Model A. 6. Conclusions This paper has presented the development of a robust fully-coupled gas-phase and particle model for the synthesis ¨ of narrowly distributed silicon nanoparticles along the centreline of the hot-wall reactor of Kormer et al. (2010a). Two models describing this system were proposed; one using the chemical mechanism of Swihart & Girshick (1999a) while adjusting the sintering parameters (Model A), and the other adjusting the gas-phase chemistry and surface reaction rates while neglecting coagulation and sintering (Model B). A sophisticated estimation methodology was employed to find parameters and their error bounds for Models A and B. Both models were able to obtain approximate fit to experimental data for six of the cases at different process conditions. In Model A, the critical nucleus size was varied to identify an appropriate nucleus size for calculations. It was found that a nucleus number of two gave best-fit to experimental results and was consistent with arguments made through homogeneous nucleation theory. This model gave excellent fit to the base-case, but overly-narrow PSDs for cases at lower temperatures and higher pressures, attributed to the coagulation-sintering growth mechanism identified. However, the parameters found through the estimation process predicted unphysical sintering rates and departed significantly from previous literature results. This caused Model A to be rejected. Model B adjusted the gas-phase and surface-growth rates while neglecting coagulation. Superior agreement with ¨ experimental results as compared to Model A and the model of Kormer et al. (2010b) was obtained. Formation of narrow

60

W.J. Menz et al. / Journal of Aerosol Science 44 (2012) 46–61

PSDs was attributed to short nucleation periods and high surface growth rates. It was found that enabling coagulation at higher pressures allowed for better correspondence with experimental PSDs. It was also suggested that this model could be improved by developing a new gas-phase mechanism from first principles via quantum-mechanical calculations, tailored specifically to these process conditions. At this stage, a constant critical nucleus size of two silicon atoms is assumed in both models. This is an accurate assumption at the initial process conditions, however, could become unreliable as the precursor is depleted. Finally, the model does not incorporate any radial coordinate and currently only seeks to understand what occurs through the centreline. Coupling it to a computational fluid dynamics code could substantially improve its robustness and potentially allow for better modelling of wall-losses.

Acknowledgements W.M. acknowledges financial support provided by the Cambridge Australia Trust. The authors thank members of the Computational Modelling Group for their guidance and support. The authors also gratefully acknowledge the funding of the Deutsche Forschungsgemeinschaft through the Cluster of Excellence ‘‘Engineering of Advanced Materials’’ at the ¨ Friedrich-Alexander University Erlangen-Nurnberg. References Adamczyk, A., Reyniers, M., Marin, G., & Broadbelt, L. (2009). Exploring 1, 2-hydrogen shift in silicon nanoparticles: Reaction kinetics from quantum chemical calculations and derivation of transition state group additivity database. The Journal of Physical Chemistry A, 113, 10933–10946. Adamczyk, A., Reyniers, M., Marin, G., & Broadbelt, L. (2010). Kinetics of substituted silylene addition and elimination in silicon nanocluster growth captured by group additivity. ChemPhysChem, 11, 1978–1994. Bhandarkar, U., Swihart, M., Girshick, S., & Kortshagen, U. (2000). Modelling of silicon hydride clustering in a low-pressure silane plasma. Journal of Physics D, 33, 2731–2746. Box, G., & Tiao, G. (1974). Bayesian Inference in Statistical Analysis. Addison-Wesley: . Braumann, A., Man, P., & Kraft, M. (2011). The inverse problem in granulation modeling—Two different statistical approaches. AIChE Journal, 57, 3105–3121. Celnik, M., 2007. On the numerical modelling of soot and carbon nanotube formation. Ph.D. thesis. University of Cambridge. Celnik, M., Patterson, R., Kraft, M., & Wagner, W. (2007). Coupling a stochastic soot population balance to gas-phase chemistry using operator splitting. Combustion and Flame, 148, 158–176. Celnik, M., Patterson, R., Kraft, M., & Wagner, W. (2009). A predictor–corrector algorithm for the coupling of stiff ODEs to a particle population balance. Journal of Computational Physics, 228, 2758–2769. Dang, H., & Swihart, M. (2009). Computational modeling of silicon nanoparticle synthesis—II. A two-dimensional bivariate model for silicon nanoparticle synthesis in a laser-driven reactor including finite-rate coalescence. Aerosol Science and Technology, 43, 554–569. Frenklach, M., & Harris, S. (1987). Aerosol dynamics modeling using the method of moments. Journal of Colloid and Interface Science, 118, 252–261. Friedlander, S. (2000). Smoke, Dust and Haze. Oxford University Press: . Giunta, C., McCurdy, R., Chapple-Sokol, J., & Gordon, R. (1990). Gas-phase kinetics in the atmospheric pressure chemical vapor deposition of silicon from silane and disilane. Journal of Applied Physics, 67, 1062–1075. Goodson, M., & Kraft, M. (2002). An efficient stochastic algorithm for simulating nano-particle dynamics. Journal of Computational Physics, 183, 210–232. Grosschmidt, D., Bockhorn, H., Goodson, M., & Kraft, M. (2002). Two approaches to the simulation of silica particle synthesis. Proceedings of the Combustion Institute, 29, 1039–1046. Halpern, A., & Glendening, E. (1996). Estimating molecular collision diameters using computational methods. Journal of Molecular Structure, 365, 9–12. Hastings, W. (1970). Monte carlo sampling methods using Markov chains and their applications. Biometrika, 57, 97–109. Heine, M., & Pratsinis, S. (2007). Polydispersity of primary particles in agglomerates made by coagulation and sintering. Journal of Aerosol Science, 38, 17–38. Hirokami, T., Maeda, Y., & Tsukada, H. (2006). Parameter estimation using simultaneous perturbation stochastic approximation. Electrical Engineering in Japan, 154, 30–39. Ho, P., Coltrin, M., & Breiland, W. (1994). Laser-induced fluorescence measurements and kinetic analysis of Si atom formation in a rotating disk chemical vapor deposition reactor. The Journal of Physical Chemistry, 98, 10138–10147. Hogness, T., Wilson, T., & Johnson, W. (1936). The thermal decomposition of silane. Journal of the American Chemical Society, 58, 108–112. Houf, W., Grcar, J., & Breiland, W. (1993). A model for low pressure chemical vapor deposition in a hot-wall tubular reactor. Materials Science and Engineering B, 17, 163–171. Kobata, A., Kusakabe, K., & Morooka, S. (1991). Growth and transformation of TiO2 crystallites in aerosol reactor. AIChE Journal, 37, 347–359. Koch, W., & Friedlander, S. (1990). The effect of particle coalescence on the surface area of a coagulating aerosol. Journal of Colloid and Interface Science, 140, 419–427. Kodas, T., & Hampden-Smith, T. (1999). Aerosol Processing of Materials. Wiley: New York. ¨ Kormer, R., Jank, M., Ryssel, H., Schmid, H., & Peukert, W. (2010a). Aerosol synthesis of silicon nanoparticles with narrow size distribution—part 1: Experimental investigations. Journal of Aerosol Science, 41, 998–1007. ¨ Kormer, R., Schmid, H., & Peukert, H.W. (2010b). Aerosol synthesis of silicon nanoparticles with narrow size distribution—part 2: Theoretical analysis of the formation mechanism. Journal of Aerosol Science, 41, 1008–1019. Koshi, M., Kato, S., & Matsui, H. (1991). Unimolecular decomposition of SiH4, SiH3F, and SiH2F2 at high temperatures. Journal of Physical Chemistry, 95, 1223–1227. Kraft, M., & Mosbach, S. (2010). The future of computational modelling in reaction engineering. Philosophical Transactions of the Royal Society A, 368, 3633–3644. Kruis, F., Kusters, K., Pratsinis, S., & Scarlett, B. (1993). A simple model for the evolution of the characteristics of aggregate particles undergoing coagulation and sintering. Aerosol Science and Technology, 19, 514–526. Kruis, F., Schoonman, J., & Scarlett, B. (1994). Homogeneous nucleation of silicon. Journal of Aerosol Science, 25, 1291–1304. Lavvas, P., Sander, M., Kraft, M., & Imanaka, H. (2011). Surface chemistry and particle shape: Processes for the evolution of aerosols in titan’s atmosphere. The Astrophysical Journal, 728, 80. Li, X., He, Y., Talukdar, S., & Swihart, M. (2003). Process for preparing macroscopic quantities of brightly photoluminescent silicon nanoparticles with emission spanning the visible spectrum. Langmuir, 19, 8490–8496.

W.J. Menz et al. / Journal of Aerosol Science 44 (2012) 46–61

61

Metropolis, N., Rosenbluth, A., Rosenbluth, M., Teller, A., & Teller, E. (1953). Equation of state calculations by fast computing machines. The Journal of Chemical Physics, 21, 1087–1092. Mezey, L., & Giber, J. (1982). Calculation of surface free energies and characteristics of surface segregation of solids. Surface Science, 117, 220–231. Mosbach, S., Braumann, A., Man, P., Kastner, C., Brownbridge, G., Kraft, M. Iterative improvement of Bayesian parameter estimates for an engine model by means of experimental design, accepted for publication. Nijhawan, S., McMurry, P., Swihart, M., Suh, S.-M., Girshick, S., Campbell, S., & Brockmann, J. (2003). An experimental and numerical study of particle nucleation and growth during low-pressure thermal decomposition of silane. Journal of Aerosol Science, 34, 691–711. Patterson, R., & Kraft, M. (2007). Models for the aggregate structure of soot particles. Combustion and Flame, 151, 160–172. Patterson, R., Singh, J., Balthasar, M., Kraft, M., & Wagner, W. (2006a). Extending stochastic soot simulation to higher pressures. Combustion and Flame, 145, 638–642. Patterson, R., Singh, J., Balthasar, M., Kraft, M., & Norris, J. (2006b). The linear process deferment algorithm: A new technique for solving population balance equations. SIAM Journal on Scientific Computing, 28, 303–320. Peev, G., Zambov, L., Shanov, V., & Tserovski, L. (1991). Kinetics of polysilicon film growth by thermal decomposition of silane. Semiconductor Science and Technology, 6, 281. Perry, R., Green, D., & Maloney, J. (1997). Perry’s Chemical Engineers’ Handbook (7th ed.). McGraw-Hill: New York. Petersen, E., & Crofton, M. (2003). Measurements of high-temperature silane pyrolysis using SiH4 IR emission and SiH2 laser absorption. Journal of Physical Chemistry A, 107, 10988–10995. Pratsinis, S. (1988). Simultaneous nucleation, condensation, and coagulation in aerosol reactors. Journal of Colloid and Interface Science, 124, 416–427. Purnell, J., & Walsh, R. (1966). The pyrolysis of monosilane. Proceedings of the Royal Society of London. Series A: Mathematical and Physical Sciences, 293, 543–561. Sander, M., West, R., Celnik, M., & Kraft, M. (2009). A detailed model for the sintering of polydispersed nanoparticle agglomerates. Aerosol Science and Technology, 43, 978–989. Shekar, S., Sander, M., Riehl, R., Smith, A., Braumann, A., Kraft, M., Modelling the flame synthesis of silica nanoparticles from tetraethoxysilane. Chemical Engineering Science, doi:10.1016/j.ces.2011.06.010. In press. Smallbone, A., Brownbridge, G., Phadungsukanan, W., Kraft, M., Johansson, B. Automated IC engine model development with uncertainty propagation, SAE Paper SAE 2011-01-0237. Swihart, M., & Girshick, S. (1999a). Thermochemistry and kinetics of silicon hydride cluster formation during thermal decomposition of silane. Journal of Physical Chemistry B, 103, 64–76. Swihart, M., & Girshick, S. (1999b). Ab initio structures and energetics of selected hydrogenated silicon clusters containing six to ten silicon atoms. Chemical Physics Letters, 307, 527–532. Talukdar, S., & Swihart, M. (2004). Aerosol dynamics modeling of silicon nanoparticle formation during silane pyrolysis: A comparison of three solution methods. Journal of Aerosol Science, 35, 889–908. Tsantilis, S., & Pratsinis, S. (2004). Narrowing the size distribution of aerosol-made titania by surface growth and coagulation. Journal of Aerosol Science, 35, 405–420. Tsantilis, S., Briesen, H., & Pratsinis, S. (2001). Sintering time for silica particle growth. Aerosol Science and Technology, 34, 237–246. Vepˇrek, S., Schopper, K., Ambacher, O., Rieger, W., & Vepˇrek-Heijman, M. (1993). Mechanism of cluster formation in a clean silane discharge. Journal of The Electrochemical Society, 140, 1935–1942. West, R., Shirley, R., Kraft, M., Goldsmith, C., & Green, W. (2009). A detailed kinetic model for combustion synthesis of titania from TiCl4. Combustion and Flame, 156, 1764–1770. Woiki, D., Catoire, L., & Roth, P. (1997). High-temperature kinetics of Si-containing precursors for ceramic processing. AIChE Journal, 43, 2670–2678. Wong, H., Li, X., Swihart, M., & Broadbelt, L. (2004). Detailed kinetic modeling of silicon nanoparticle formation chemistry via automated mechanism generation. Journal of Physical Chemistry A, 108, 10122–10132. Wu, J., Nguyen, H., & Flagan, R. (1987). A method for the synthesis of submicron particles. Langmuir, 3, 266–271. Zachariah, M., & Carrier, M. (1999). Molecular dynamics computation of gas-phase nanoparticle sintering: A comparison with phenomenological models. Journal of Aerosol Science, 30, 1139–1152.