Modeling of chemical vapor deposition of silane for silicon production

Article Cite This: Ind. Eng. Chem. Res. XXXX, XXX, XXX−XXX

pubs.acs.org/IECR

Modeling of Chemical Vapor Deposition of Silane for Silicon Production in a Spouted Bed via Discrete Element Method Coupled with Eulerian Model Pongpawee Chanlaor,†,‡,§ Sunun Limtrakul,*,†,‡,§ Terdthai Vatanatham,†,‡,§ and Palghat A. Ramachandran∥ †

Department of Chemical Engineering, Faculty of Engineering, Kasetsart University, Bangkok 10900, Thailand Center of Excellence on Petrochemical and Materials Technology, Department of Chemical Engineering, Faculty of Engineering, Kasetsart University, Bangkok 10900, Thailand § Center for Advanced Studies in Industrial Technology, Faculty of Engineering, Kasetsart University, Bangkok 10900, Thailand ∥ Department of Energy, Environmental & Chemical Engineering, Washington University in St. Louis, St. Louis, Missouri 63130, United States

Ind. Eng. Chem. Res. Downloaded from pubs.acs.org by DURHAM UNIV on 09/03/18. For personal use only.

‡

ABSTRACT: A CFD−DEM coupling model was applied to the silicon production by silane chemical vapor deposition in a spouted bed to investigate local phase movements, mass transfer, heat transfer, and the deposition rate of silicon. Silicon particle growth by heterogeneous deposition on the seed particle and fines formation in the gas phase were included. In addition, the scavenging effect was also considered in the particle growth mechanism. The spouted bed has advantages to promote heterogeneous deposition time due to a larger solid holdup and higher gas residence time in the annulus region. Furthermore, the effects of operating conditions, i.e., inlet gas temperature, silane composition, and gas velocity, on the reactor performance were investigated in detail, and some guidelines are provided for the choice of operating conditions. Verification of solid movement results with the experimental data and correlations show good agreement. In addition, the axial profile of bed voidage agrees well with correlation.

1. INTRODUCTION Solar energy is an alternative energy source in which the energy is generated by using solar cells, a concept known for over a century and now getting increased attention due to focus on being green. The dominant material used in photovoltaic cells is silicon, which is required in six9s quality (99.999999%). Polycrystalline silicon production is the starting point of the production of wafers over which the cells are fabricated. This part of the process takes up approximately 20% of capital cost of solar cells,1 and there is a great incentive to reduce the cost using fundamentally based engineering tools. The chemical process for producing poly silicon is chemical vapor deposition (CVD) from precursors; trichlorosilane, SiHCl3 (TCS), and monosilane, SiH4 (silane), are most commonly used. The precursor is decomposed by thermal decomposition in the case of silane, while it is reacted with hydrogen in the case of TCS. Silicon product will be then deposited on a starting substrate surface to achieve high-purity polysilicon. Two types of reactors are commonly used for the process: (i) a Siemens reactor and (ii) a fluidized bed reactor. The workhorse of today’s industry is the Siemens reactor, where either TCS or silane is deposited on silicon rods. Several improvements have been made to the design, resulting in a substantial reduction in energy consumption.2 However, the need for maintaining © XXXX American Chemical Society

large temperature differences between the deposition surfaces and the reactor wall, as well as the limited deposition surface area, cause a limit to theoretical reactant yield and lowering of the energy consumption. TCS has proven successful for these reactors among several reasons because some reverse reactions will remove loosely bound solid silicon and thereby continuously ensure high quality on the deposited material. The drawback is that the products of some of these reactions may form gaseous species that are not decomposed at the reactor operating temperature, and this in turn causes a lower reactant gas yield. The energy requirements are also high in the Siemens reactor because a large amount of energy is wasted in the process of heating silicon rods to reach very high temperatures. Thus, further reduction in the energy consumption in order to reduce production costs is difficult.3 An alternative to the Siemens reactor is the fluidized bed reactor (FBR), which has excellent heat and mass transfer, suitable for chemical vapor deposition. In a fluidized bed reactor silane is commonly used as the reactant. This process requires a Received: Revised: Accepted: Published: A

June 22, 2018 August 10, 2018 August 22, 2018 August 22, 2018 DOI: 10.1021/acs.iecr.8b02794 Ind. Eng. Chem. Res. XXXX, XXX, XXX−XXX

Article

Industrial & Engineering Chemistry Research lower reaction temperature, i.e., 600−800 °C, and gives a high production rate. Silicon particles are used as substrate and grow along the process. The particle product has acceptable purity and is more convenient to process in the next steps of solar cell production. However, the complex interaction of the hydrodynamic behavior of gases and particles, heat transfer, mass transfer, and chemical reaction in the reactor makes the design and scale-up of these reactors a difficult issue. In silane pyrolysis, the chemical vapor deposition of silicon occurs in two major pathways. The first pathway is homogeneous decomposition that leads to the formation of new solid fines of silicon which are not useful for further processing. The second one is heterogeneous decomposition, in which the reactant gas gets decomposed and the product is deposited on a seed particle surface. The seeds grow with time and produce larger sized particles which are useful for further processing. It is therefore important to promote heterogeneous growth on the seed particles and reduce the fines formed by the gas-phase reaction. The flow regime of the column can also be difficult to control if there is a wide distribution of particle size. The modeling of the reactor, in order to understand the various interconnected phenomena, is important for the design, scaling up of this reactor, and reducing the energy costs. Hence both phenomenological4,5 as well as computational fluid dynamics (CFD) based models have been proposed. Most CFD based models have used the two-fluid flow model (Euler−Euler model) to investigate the phenomenon of silicon production fluidized bed reactors.6−10 The discrete element method (DEM) or discrete particles model is an alternative way of modeling a fluid bed reactor. No study using this tool to model fluid bed pyrolysis of silane has been reported, and this is the main focus of this present work. A discrete element simulation, first proposed by Cundall and Strack,11 tracks the motion of individual particles while the gas phase is treated as a continuum. Hence this model provides more realistic results than the two-fluid flow model, particularly in the system with large particle and dense bed, i.e., spouted bed. This model is also known as an Euler−Lagrangian model and has been applied to monitor multiphase flow in fluidized beds.12−16 The DEM model can be combined with heat transfer,17,18 mass transfer, and chemical reaction.19,20 This combined model can predict the hydrodynamics, concentration in gas, concentration in solid phases, and local temperature throughout the beds. The advantage of the DEM model is that it provides the tracking of particles in the system and, correspondingly, the calculation of growth and size distribution of the seed particles. The growth of a particle takes place by heterogeneous reaction and also by the capture of the fine particles by scavenging, and these effects can be more accurately tracked in the DEM setup than in a two-fluid model. Eulerian continuum model cannot directly provide particle tracking, particle growth, and size distribution. Unless it is coupled to detailed population balance equation, for example, see Passalacqua and co-workers’ work.21 However, the CFD−DEM coupling model has advantages that it can predict discrete particle growth and size distribution directly and also included the scavenging effect. Thus, DEM-based model is a proper method for applying in a chemical vapor deposition system. The local temperatures in each zone and on the particles lead to different reaction rates, and the particle holdup in the system is also not uniform. These effects can be captured more effectively with DEM, and strategies for wall temperature profiling and control of the reactor can be proposed based on the computational model.

A spouted bed arrangement is used in the simulation in this work. This is a modification of the fluid bed wherein the inlet gas is introduced by a nozzle somewhat above the bed. The particles are fluidized in the region above this nozzle, while in the bottom part, the solids circulate as a spout. Some advantages of this design are as follows: (i) larger solid holdup in the annulus region which permits heterogeneous reaction and fine capture, (ii) smaller gas residence time and bubble size is reduced decreasing fines formation, and (iii) more intimate contact with larger particles and fines promoting scavenging of fines. A quantitative prediction of some of these effects can be accomplished by the DEM-based model, which is suitable for investigating behavior of large particles in dense area of a spouted bed. The paper provides some information on the efficiency of the spouted bed arrangement. This work developed a combined model of discrete element method (DEM) coupled with Eulerian model to investigate the reactor performance of silicon production by silicon chemical vapor deposition in a spouted bed. The combination of discrete element method and Eulerian model is usually called CFD− DEM coupling method. The combined model included both heterogeneous deposition on seed particles and homogeneous fines formation reaction, as well as heat transport effects providing unmeasurable local phase movements, mass transfer, heat transfer, and the deposition rate of silicon in detail. This local information is unavailable. However, it is important for analysis, design, and scale up of reactors especially chemical vapor deposition reactors with highly nonisothermal.

2. MATHEMATICAL MODELS The mathematical model used to predict the reactor performance in terms of silicon production consists of three sets of governing equations: momentum equations for particle and gas flows, energy equations for gas and solid phases, and mass equations for gas and solid phases. The particle movement was investigated by the discrete element method which was proposed by Cundall and Strack,11 known as the Lagrangian model. This model has been applied variously to studies of flows, heat transfer, and mass transfer. Meanwhile, the flow, heat transfer, and mass transfer in the gas phase were studied by an Eulerian model. The heat and mass transfer of the particle phase were developed by focusing the group of particles in a fluid cell. The details are explained below. 2.1. Equations of Particle Motion. The details of particle movement in terms of solid velocity in a fluidized bed can be predicted by Newton’s second law of motion. DEM model monitors individual particle movement. It provides the location of particles in the bed leading to the distribution of solid holdup, which is related to the force acting on the particle, consisting of particle contact force and force from surrounding fluid. Each particle moves with acceleration from translational and rotational motions as given in eqs 1 and 2, respectively. The acceleration of a particle (a⃗) is caused by the total force acting on the particle while the angular acceleration of the particle (α⃗ ) depends on the torque (T⃗ ) and the moment of inertia of the particle (I). a⃗ =

α⃗ = B

⃗ Ftotal + g⃗ mp

(1)

T⃗ I

(2) DOI: 10.1021/acs.iecr.8b02794 Ind. Eng. Chem. Res. XXXX, XXX, XXX−XXX

Article

Industrial & Engineering Chemistry Research The total force is composed of a particle contact force (F⃗ C) and a force from surrounding fluid (F⃗ D) which can be written as ⃗ = FC⃗ + FD⃗ Ftotal

summation of all these contact forces. This is written in eqs 9 and 10: N

(3)

⃗ = Fc,nor

The velocities, written as the function of acceleration, after the time step Δt are shown as vs⃗ = vs,old ⃗ + a ⃗Δt

(4)

ωs⃗ = ωs,old + α⃗Δt ⃗

(5)

N

⃗ = Fc,tan

⃗ ij ∑ Fc,tan, (10)

j=1

where i is the particle under consideration and j refers to other adjacent particles and/or walls and N the total number of adjacent particles/walls in the neighborhood. The force from the surrounding fluid (F⃗ D) can be calculated by Gidaspow’s drag model. ÅÄ πd p3 ÅÅÅ β ÅÅ r direction: FDr = (ug,⃗ r − vs,⃗ r ) ÅÅ 6 ÅÅ (1 − εg) Ç ÑÉ ∂p ÑÑÑÑ + (1 − εg) ÑÑ ∂r ÑÑÑ (11) Ö Ä Å πd p3 ÅÅÅ β ÅÅ z direction: FDz = (ug,⃗ z − vs,⃗ z) Å 6 ÅÅÅ (1 − εg) Ç ÉÑ ∂p ÑÑÑÑ + (1 − εg) ÑÑ ∂z ÑÑÑ (12) Ö

Figure 1. Models of contact forces between two spherical particles in the normal direction (left) and in the tangential direction (right).

where FDr and FDz are the force exerted on the group of particles in the r and z directions, respectively. 2.2. Equations of Gas Motion. The gas movement in terms of gas velocity (u⃗g) and gas volume fraction (εg) in a fluidized bed reactor was predicted by an Eulerian model. The model consists of the terms of convection and momentum transfer between gas and solid, as shown in eq 13:

The discrete element method monitors the particle movement in three dimensions. The particle contact force is calculated based on the overlapping of particles according to a soft particle model. Once particles are in contact with each other, the model allows the particles to overlap. The distance of this overlap is particle displacement (δ), which is the key to calculating particle contact force by the concept of spring (stiffness), dash-pot (dissipation), and slider (friction). The contact force is a function of the stiffness coefficient, damping coefficient, and particle displacement. A particle itself has both normal and tangential contact forces, as shown in eqs 6 and 7, respectively: ÷◊÷ ̇ ⃗ − η δnor ⃗ = −kstif δnor Fc,nor damp (6)

∂(εgρg ug⃗ )

+ (ug⃗ ·∇)εgρg ug⃗ = −∇εgp + Fg⃗ → s (13) ∂t where F⃗ g→s is the drag force of fluid on particles and is given by Fg⃗ → s = β(vs⃗ − ug⃗ )

(14)

The momentum coefficient β depends on the void fraction of the bed.22 The eq 15 is Ergun’s equation valid for ε ≤ 0.8, while eq 16 is Wen and Yu’s model valid for ε > 0.8. The combined model valid for the entire range of holdup is called Gidaspow’s drag model. ÄÅ ÉÑ (1 − εg)μg (1 − εg) ÅÅÅ ÑÑÑ ÅÅ150 ÑÑ β= + 1.75 ρ ε | v − u | ⃗ ⃗ ÑÑ g s g g 2 Å ÑÑ dp d pεg ÅÅÅÅ ÑÖ Ç (15) for ε ≤ 0.8

(7)

⃗ and δtan ⃗ are the particle displacements in the normal where δnor ÷◊÷ and tangential directions, respectively, δ ̇ the relative particle velocity of the contact particles, kstif the stiffness of the spring, and ηdamp the damping coefficient. If the value of contact force in the tangential direction in eq 7 is greater than the friction force in the normal direction, the contact force in the tangential direction will be replaced by the force from friction: ⃗ = −μ |Fc,nor ⃗ | Fc,tan f

(9)

j=1

where Δt is the time step increment. To calculate the contact force between two spherical particles, the soft particle model is applied by using the concept of spring, dash-pot, and friction slider as first proposed by Cundall and Strack and as shown in Figure 1.

÷◊÷ ̇ ⃗ − η δtan ⃗ = −kstif δtan Fc,tan damp

⃗ ij ∑ Fc,nor,

β=

3 |vs⃗ − ug⃗ |ρg (1 − εg) −2.7 CD εg 4 dp

for ε > 0.8 (16)

The drag coefficient, CD, can be calculated by

(8)

C D = 24(1 + 0.15Re 0.687)/Re

where μf is the friction coefficient. In the reactor, each particle is in contact with many neighboring particles and/or walls. Thus, the total contact force involved for considering particle i is obtained by the

C D = 0.43

for Re > 1000

for Re ≤ 1000

(17) (18)

where Re is the Reynolds number C

DOI: 10.1021/acs.iecr.8b02794 Ind. Eng. Chem. Res. XXXX, XXX, XXX−XXX

Article

Industrial & Engineering Chemistry Research Re =

|vs⃗ − ug⃗ |ρg εgd p μg

where Cps is the heat capacity of the solid particle. As mentioned previously, the solid temperature (Ts) in eq 21 is the average temperature of those of the group of particles in the fluid field (as described in Figure 2). 2.4. Species Transport Model. In the reactor, the gas phase is composed of silane, hydrogen, and silicon fines. The conservation equation of mass in the gas phase uses the concept of Eulerian model, and the continuity equation of the gas phase is written as

(19)

The spouted bed is divided into small fluid cells as shown Figure 2. Axial symmetry coordinates were applied to the fluid

∂(εgρg ) ∂t

+ ∇·(εgρg ug⃗ ) = ṁ gs

(22)

where ṁ gs is the interphase mass transfer from gas to solid. The species transport model is related to convection, conduction, and reaction terms: ∂(εgCi) ∂t

field. Each fluid cell consists of gas phase contacting a group of solid particles. The particle movement is monitored in a threedimensional system. To couple discrete element method for particle movement with an Eulerian model for gas movement, the solid information monitored in three-dimensional coordinates was averaged into information on the group of particles in a two-dimensional axisymmetric system. Thus, the momentum, heat, and mass transports between a group of particles and gas phase in a fluid cell are considered. 2.3. Conservation of Energy. The energy balance in the gas phase is in terms of convection, conduction, reaction, and energy transfer between gas and a group of solids and heat of reaction in the gas phase:

∂t

(23)

where Ci is the concentration of species i in the gas phase, Di,m diffusivity of species i in the gas medium, and Ri the reaction rate of species i. Silicon fines are considered to constitute a pseudohomogeneous phase. This phase is dealt as a gaseous species based on aerosol in this work. 2.5. Chemical Vapor Deposition of Silicon. The precursor gas is decomposed via two major pathways, one involving a homogeneous reaction and the other a heterogeneous reaction. The first reaction is the heterogeneous decomposition on the particle surface providing the desired solid silicon product, while the second reaction produces fines particles which are not useful for further processing. However, the silicon fines can be captured by silicon particles according to a scavenging effect as shown in Figure 3.

Figure 2. Spouted bed reactor with fluid field as a two-dimensional model.

̂ Tg) ∂(εgρg Cpg

+ ∇·(εgug⃗ Ci) = ∇·(∇(εgDi ,mCi)) + R i

̂ ug⃗ Tg) + ∇·(εgρg Cpg

= ∇·(εgλg ∇Tg) +

6(1 − εg)h dp

(Ts − Tg) + Q Rxn,hom (20)

where Cpg is the heat capacity of the gas phase, λg the thermal conductivity of the gas, h the heat transfer coefficient, dp the particle diameter, and Tg the gas temperature. Silane decomposition in the gas phase is a homogeneous reaction to produce silicon fines; it is a weakly exothermic reaction, and its heat of homogeneous reaction is represent by QRxn,hom. The solid temperature is represented by the average temperature of a group of particles in the fluid cell. Thus, in the energy equation for the solid phase, accumulation of energy arises from the energy transferred from the gas phase to the solid phase and the heat of heterogeneous reaction in the solid phase (QRxn,het): ∂((1 − εg)ρs CpsTs) ∂t

=

6(1 − εg)h dp

(Tg − Ts) + Q Rxn,het

Figure 3. Schematic diagram of silicon chemical vapor deposition reaction.

The rates of both heterogeneous and homogeneous reactions used in this work have been proposed by Lai and co-workers.4 Meanwhile, the homogeneous kinetic rate constants as a function of temperature was proposed by Hogness et al. (1936)23 as follows: k hom = 2 × 1013 exp( −26000/Tg) s−1

The heterogeneous reaction rate constant is as follows:

(24) 24

k het = 2.79 × 106 exp(− 19530/Ts) m 3 gas/(m 2 particle· s)

(21)

(25) D


Article

Industrial & Engineering Chemistry Research Furthermore, Lai and co-workers4 proposed the scavenging coefficient (ksc), which is related to a single large-particle collection coefficient (E) and particle diameter: ksc = EU

3(1 − εg) 2d pεg

s −1 (26)

where εg and U are the bed voidage and the superficial velocity of gas, respectively. The single large-particle collection coefficient represents a sum of particle capture efficiency from impaction, interception, and diffusion. In the fluidized bed reactor, the mechanism that dominates the particle capture is diffusion. Hence, the single large-particle collection is related to the Peclet number (Pe) in terms of mass transfer diffusion: E = 2Pe−2/3

(27)

where the Peclet number, Pe, is defined as Pe =

d pU D

. In the

preceding, the diffusivity coefficient, D, is defined by D =

kBzTg

Figure 4. Spouted bed reactor geometry.

,

Table 1. Silicon Seed Particle Properties

3πμf d pf −23

where kBz is the Boltzmann constant, 1.38066 × 10 J/(molecule·K). With information about silicon produced from the gas phase on the particles, the CFD−DEM model can provide silicon mass deposited from both heterogeneous reaction and scavenging effect in a local fluid cell (eq 28). This is used to evaluate a new particle diameter in that cell as shown in eq 29:

properties initial particle diameter (mm) particle density (kg/m3) heat capacity at 293 K (J/kg·K) thermal conductivity at 293 K (W/m·K)

ji (1 − εg) zy MSi,deposited = jjjj6 εgk hetCSiH4 + εgkscCSi zzzzVcellM wSi dt j z dp k {

MSi,0 + (MSi,deposited /n p,cell) zyz ji 3 zz d p,new = 2jjjj × z ρs k 4π {

2.0 2329 712 148

differential equations were solved by finite difference method. Semi-implicit method for pressure linked equation (SIMPLE), developed by Patankar (1980),25 was applied. The program was coded in FORTRAN, in which the calculation domain was divided into cells. The size of a fluid cell is necessarily smaller than the macroscopic motion of particles but larger than individual particle size. In this work, a size of 9.5 × 16.4 mm (Δr × Δz) was used. The time step cannot be too large to have stability of the calculation. Tsuji et al. (1993)15 suggested that the time step (Δt) should follow the below equation for stable calculation: π Δt = m /kstif (30) 5

(28)

1/3

(29)

3. NUMERICAL METHOD Silane decomposition in a spouted bed is the investigated system. The reactor is cylindrical and has the dimensions of 0.152 m I.D. and 1.0 m height, with a conical gas plenum with an inclined angle of 60°, as shown in Figure 4. The properties of silicon seed particles used in the simulation are shown in Table 1. The two-dimensional differential equations of momentum, species balance, and energy balance for gas phase were considered in an Eulerian framework while the three-dimensional particle motion is considered in a Lagrangian framework. The particle motion equation provides the movement tracking of individual particle. In Eulerian−Lagrangian flow mapping, each fluid cell was considered to contain a group of particles interacting with the fluid as shown Figure 2. When solving algebraically, interaction between gas and particles was included in the particle motion calculation. Solid and gas motion equations of the Eulerian−Lagrangian framework were simultaneously solved to provide particle positions, particle velocities, and gas velocities. The voidage in each fluid cell was obtained from the information on particle positions (see Figure 2). The energy and species equations of solid phase written for the particle group in a fluid cell, combined with the differential equations of energy and species equations of gas phase, were simultaneously solved to give locally averaged concentrations and temperatures in both gas and particle phases. The

where m is the mass of a particle. This simulation was carried out by using the time step of 2 × 10−5 s.

4. PARAMETERS ESTIMATION AND SIMULATION CONDITIONS Spouted beds have been widely investigated for use in various physical and chemical processes. Many studies have focused on developing mathematical models to predict the hydrodynamics parameters, e.g., minimum spouted velocity, maximum spoutable bed height, and spouted diameter. 4.1. Minimum Spouted Velocity (Ums). The minimum spouted velocity is a critical parameter for the selection of operating gas velocity. The correlation that was used in this work were provided by Mathur and Gishler (1955):26 Ä ÑÉ1/2 1/3Å ij d p yzij Di yz ÅÅÅÅ 2gH0(ρs − ρg ) ÑÑÑÑ j z ÑÑ Ums = jj zzjjj zzz ÅÅÅ ÑÑ j Dc zj Dc z ÅÅ ρg ÑÑÖ k {k { ÅÇ (31) where Di and Dc are the inlet (orifice) diameter and column diameter, respectively, and H0 is the static solid bed height. 4.2. Maximum Spoutable Bed Height (Hm). The maximum spoutable bed height reflects a limitation on stable E


Article

Industrial & Engineering Chemistry Research

same in both works, as shown in Table 3. Figure 5 shows the axial solid velocities in the spouted, annulus, and fountain zones obtained from this work and the experimental data. Figure 5a shows that in the spouted zone, solid particles move upward with very high velocity in the central zone but very low near annulus zone. In the annulus zone, the particles flow downward with relatively low velocity, as seen in Figure 5b. The radial profiles in the fountain zone seen in Figure 5c illustrate that the particles move upward near the center but move downward near the wall. Comparison shows that the simulation results from this work agree well with the experimental results for all three zones. The average errors of simulation results are 8.4%, 6.0%, and 8.7% in the spouted, annulus, and fountain zones, respectively. 5.2. Comparison Empirical Model for Spouted Bed with the Simulation Results. This section shows the comparison of the spouted diameter, the superficial gas velocity in the spout, and the annulus regions from the simulation results with known correlations. The spouted diameter (ds) is another important parameter which is essential in order to understand the flow pattern of fluid in the reactor bed. The correlation for spouted diameter in this work was proposed by McNab (1972):30 ÄÅ É ÅÅ (ρ U )0.5 D0.68 ÑÑÑ ÅÅ g ÑÑ c ÑÑ ds = 2.0ÅÅÅ ÑÑ 0.41 ÅÅ ρ ÑÑ ÅÅÇ b (33) ÑÖ

spouting behavior. It is related to the amount of material that can be processed in the spouted bed. The correlation was proposed by Malek and Lu (1965)27 in SI units as follows: ji D zy Hm = 418jjjj c zzzz j dp z Dc k {

0.75

ij Dc yz −1.2 −2 jj zz ρ ϕ jj D zz s k i{ 0.4

(32)

where ϕ is the shape factor of particles. This correlation is suitable for 0.01 < dp/Dc < 0.02, a range in which the ratio used in this work falls (dp/Dc = 0.013). The minimum spouted bed velocity and maximum spoutable bed height were estimated via eqs 31 and 32, respectively. Simulation conditions and parameters are shown in Table 2. Table 2. Simulation Conditions and Parameters gas phase inlet mole fraction silane:nitrogen minimum spouted velocity (Ums) (m/s) superficial gas velocity (=1.30Ums) (m/s) inlet gas temperature (K) initial bed temperature (K) wall temperature (K) solid phase particle loading, particles (initial loading 4.87 kg) stiffness (N/m) friction coefficient coefficient of restitution maximum spoutable bed height (Hm) (cm) static bed height (H0) (cm)

0.2:0.8 1.103 1.434 473 873 973 500000 800 0.3 0.9 34 28

where ρb is the bulk density. A two-region model for calculating the velocity in the spouted bed was proposed by Mathur and Lim (1974).31 Smith et al. (1982)28 studied the reactor model so as to investigate the superficial gas velocity in the annulus region32 and spouted region.33 From these studies, the correlations are provided as below: ÄÅ É 3Ñ ÅÅ ij yz ÑÑÑÑ ÅÅ z a zz ÑÑ Uz (z) = Umf ÅÅÅ1 − jjj1 − zz ÑÑ j ÅÅ H m k { ÑÑÖ (34) ÅÇ

5. MODEL VALIDATION Hydrodynamics in a spouted bed obtained from our model is validated with two sets of data, experimental data and data from the empirical model of Smith et al. (1982).28 The dimension of the investigated spouted bed are provided in Figure 4. Comparison between our simulation work and the experimental work of Ren and co-workers (2011)29 is carried out based on the conditions listed in Table 3.

Uzs(z) =

Table 3. Conditions and Parameters Used in Verification Case

(35)

where Uaz,Usz are the superficial gas velocity in the annulus and spouted regions, respectively, and Ac, Aa, and As are cross sectional area of column bed, annulus region, and spouted region, respectively. The calculation of spouted diameter via eq 33 gives 38.0 mm, while the simulation results of spouted bed diameter give 36.0 mm. This shows that the simulation result agrees quite well (5.6% relative error) with the correlation used for spouted bed diameter. Figure 6 shows the comparison of gas velocity in the spouted region and that the annulus region obtained from correlations with those obtained by simulation. For the spouted region, the correlation for superficial gas velocity used was eq 34, while the results for the gas velocity in the annulus region were obtained by eq 35. It is found that, relative to the results of this work, the mean relative error of the axial gas velocity in spouted and annulus regions are 5.8% and 4.7%, respectively. In addition, Smith and co-workers28 showed the correlation for bed voidage in a spouted region as a function of bed height (εs) as shown in eq 36.

properties particle diameter (mm) particle density (kg/m3) gas density (kg/m3) gas viscosity (Pa/s) initial bed height (mm) minimum spouted velocity (Ums) (m/s) superficial gas velocity/minimum spouted velocity number of particles

Ac A U0 − a Uza(z) As As

2.9 1450 1.205 1.81 × 10−5 325 1.15 1.30 217552

5.1. Comparison of Radial Distribution of Particle Velocity with Experimental Results. Verification of the simulation results was carried out by comparing the radial profiles of axial solid velocities at various bed heights in three zones of the spouted bed with the experimental results obtained by Ren and co-workers (2011).29 Both experiments and simulations were carried out in a spouted bed with the same geometry and size. Particles of 2.9 mm in diameter and a superficial gas velocity of 1.30Ums were used in the comparison. All conditions are the F


Article


Figure 5. Comparison of solid velocities obtained from experiments and simulation in (a) spouted zone, (b) annulus zone, and (c) fountain zone.

can only with difficulty describe the solid information in a threedimensional system. Therefore, the axially sliced quantities in the central region of the solid particle phase are shown in a twodimensional domain. The radial profiles that are provided in the results were averaged axially within the range 0.1 < z/Z < 0.3 to represent information in the dense particle area. 6.1. Key Hydrodynamic Features. Figure 8a shows the particle distribution in the spouted bed as a function of time. It was found that there are very low solid holdup in the central region and high solid holdup in the annulus zone due to feeding reactant at the central bottom of the bed, as seen in Figures 8a and 9. Above these two zones, a fountain was found where the particles are lifted in the spouted region and create the fountain because of high gas velocity in the central zone. Figures 10 and 11 show the axial solid velocity and the vector plots of for both gas and solid, which confirm that the solid particles move upward with very high velocity in the central zone and turn down in the fountain zone. The particles fall downward in the annulus with downward velocity, especially near the bottom of the bed (z/Z = 0.1), as seen in Figure 11, leading to circulation flow. Particles in annulus zone behave like a moving bed, which increases the residence time of solid to contact with the reactant gas. This in turn induces the heterogeneous reaction giving the desired solid silicon solid product. In the central zone, the gas velocity is high, leading to short gas residence time. This then reduces silicon fines formation in the gas phase. A spouted bed is a good design to minimize fine formation and maximize silicon solid product. A spouted bed with short gas residence time in the core zone and large residence time of solid particles in the annulus shows the above advantages. 6.2. Temperature Distributions. The temperature distributions of solid and gas are shown in parts a and b of Figure 8, respectively. Initially, temperatures of both gas and solid particles are 873 K. At later times, in the core zone, the gas temperature is lower than the solid temperature due to the cooling effect from the feed. The gas and solid particle temperature in the annular region of the bed are not different from one another and still high because of less cooling effect from the feed and also because the flow in the annular region is much smaller than in the spouted zone, as seen in Figure 12. These effects together led to a higher heterogeneous reaction rate in the annular zone. 6.3. Silane and Fines Profiles. Silane conversion presented in the terms of unconverted silane concentration in the gas phase at different times is shown in Figure 13a. In the center zone of the bed, the silane concentration is the highest, and it decreases along the radial position. In the central region, the solid holdup is very low and the temperature is low, leading to lower heterogeneous and homogeneous reaction rates. Thus, low

Figure 6. Comparison of averaged axial gas velocity with empirical model in (a) spouted region and (b) annulus region.

ε s (z ) = 1 − β s

z H

(36)

where β is the spouted porosity constant, which is equal to 0.2.28 The axial profile of bed voidage in the spouted region obtained by this work was compared with that calculated by the correlation as shown in Figure 7. The results of this work were in s

Figure 7. Comparison of axial profile of bed voidage in spouted region with empirical model.

the range of 0.87−0.92 and agreed well with that from the correlation, fell in the range of 0.86−0.93. The mean realtive errors of simulation results is 0.99%.

6. RESULTS AND DISCUSSION The solid particle movement was monitored in a threedimensional domain, while the fluid flow, temperature, and gas concentration distributions were calculated in a two-dimensional axisymmetric domain. The average temperature of groups of particles in a fluid cell was obtained and shown in a twodimensional axisymmetric domain. The particle visualization G


Article


Figure 8. Distributions of solid temperature (a) and gas temperature (b) at various times.

Figure 12. Averaged radial profiles of gas and solid temperature.

Figure 9. Averaged radial profile in the dense zone of solid holdup.

residence time and high temperature were found. Therefore, a high heterogeneous reaction rate takes place here, leading to high conversion. Figure 13b shows distribution of silicon fines concentration which are produced by homogeneous decomposition of reactant gas. In this work, the silicon fines are considered to form a pseudohomogeneous phase, which is dealt with as a gaseous species in the fashion of an aerosol. The silicon fines concentration is low in the particle-dense phase, where the heterogeneous reaction is dominant. The fines concentration in the freeboard, where the solid silicon holdup is low, is high with high homogeneous reaction rate, as shown in Figures 13b and 14a. The axial profile confirms that the silicon fines were produced more in the freeboard area due to there being less solid content and high gas temperature in the freeboard zone, allowing for more homogeneous decomposition to form fines. To reduce the homogeneous reaction in the freeboard region, a cooler wall temperature should be maintained to suppress the homogeneous reaction rate there. Figure 14b shows the radial profile of silicon fines concentration at different bed heights. The concentration is low in the spouted region (central area) and higher in the annular zone. Although there is low solid content in the spouted region, fines are formed less due to there being a very high gas velocity, leading to low gas residence time. In addition, the temperature in the core zone is low here. Thus, silicon fines concentration is lower in this zone. On the other hand, in the annular zone, even though the solid holdup is high, the homogeneous reaction rate is also high because of the high wall temperature. Hence, more silicon fines are formed in the annulus area than in the spouted zone. The above results indicate that a spouted bed with short gas residence time in the core zone and large residence time of solid particles in the annulus shows advantages for reducing

Figure 10. Vector plot for gas (left) and particles (right).

Figure 11. Radial profile of averaged solid axial velocity.

conversion was found in this zone. On the other hand, in the annular zone, where there is high solid loading, long solid H


Article


Figure 13. Distributions of unconverted silane concentration (a) and silicon fine concentration (b) at various times.

Figure 14. Dimensionless silicon fines concentration in axial position (a) and radial position (b).

Figure 15. Silicon production rates by radial position (a) and axial position (b).

fine formation and inducing the formation of silicon solid product. Figure 15 shows the profiles of silicon production rates in the radial and axial directions. The radial profile of silicon solid production (solid line) is higher in the annulus region than that in the spouted zone. This confirms that more heterogeneous reaction takes place in this area. On the other hand, in the spouted zone, where there is very low solid content, the heterogeneous reaction is very sparse. The axial profile of the solid silicon production rate demonstrates that it is higher in particle dense area and approaches zero in the freeboard. This indicates that the heterogeneous reaction takes place in the particle-dense area. Figure 15 also shows comparison of the fine silicon and solid silicon productions. The solid silicon production is higher than fine silicon production for all reactor radial positions. In the spouted region, the fine silicon production rate (dash line) is less than that of the solid silicon, although there is less solid content.

This is due to the gas temperature being cooler than the solid temperature, caused by the cooling effect of feed gas (see Figure 12). In the annular region, the solid silicon production is obviously higher than the fines production. This is due to high solid content inducing more heterogeneous reaction in this area. Figure 15b shows the comparison of the axial profiles of silicon production rates. The solid silicon production is higher than fines production in the particle-dense area, while the fine product is obviously produced more in the freeboard area. This indicates that the heterogeneous reaction is dominant in the particle-dense area, and the homogeneous reaction takes place in the freeboard region, where there is less solid content. To minimize the homogeneous reaction rate in the freeboard, the operation suggestion is that the wall temperature should be cooler. Figure 16 shows the overall mass balance for silicon decomposition. The heterogeneous reaction provides 76.12% of silicon I


Article


Table 4. Silicon Solid and Fines Production Rates, Particle Growth Rate, and Standard Deviation of Particle Size Distribution at 15 s for Different Operating Conditions operating conditions

solid Si production rate (mg/s)

Feed Gas Temperature 300 K 134.2 473 K 147.1 673 K 161.3 Inlet Silane Concentration 20% 147.1 30% 218.4 40% 290.6 Superficial Gas Velocity 1.30Ums 147.1 1.52Ums 175.9 1.63Ums 205.1

Figure 16. Overall mass balance indicating the silicon production rates of various pathways for 20% feed silane composition.

product deposited on the particles’ surface. On the other hand, 23.87% of the silane is decomposed by the homogeneous reaction to form fines. As the freeboard is excluded from the results, the fines formation is reduced to 12.44%, which it is in the range of values reported by Lai and co-workers (1986).4 This indicates that the temperature in the freeboard region should be cooled down to suppress the fines formation in that region. Only 0.1% of fines are captured by scavenging effect to form a solid silicon product. This is due to large initial size of silicon particles. The scavenging coefficient that was calculated by the equation from Lai’s work depends on the seed particle. At the same operating condition, when the particle size decreased by 50%, the scavenging effect is increased 2.24 times. Comparing to other works, the initial particle size in this work is quite large, leading to low scavenging coefficient. This shows that the heterogeneous reaction is the dominant mechanism of particle growth, while the scavenging of silicon product has only a slight effect on the particle growth. This is probably because of the initial size of particles in this work is large, leading to low particle surface area. Thus, a low degree of scavenging effect on the particle surface was found. 6.4. Particle Size Distribution. Figure 17a shows the distribution of the silicon particle size change as (dp,t − dp,0). The

fines production rate (mg/s)

particle growth rate (nm/s)

standard deviation of PSD (nm)

37.00 55.91 66.21

21.18 22.14 23.19

32.3 13.7 7.7

55.91 97.79 147.5

22.14 33.07 44.01

13.7 20.6 28.1

55.91 38.15 30.59

22.14 22.12 22.06

13.7 22.6 20.9

longer is the operating time, the broader is the particle distribution, indicating nonuniform deposition on the solid particles. 6.5. Effect of Inlet Gas Temperature. In this section, the effect of inlet gas temperature on the silicon production is discussed. The inlet gas composition consists of 20% silane, and the inlet gas velocity is 1.434 m/s, while inlet gas temperatures of 300, 473, and 673 K are used in the studies. The wall temperature was relatively high and kept constant at 973, whereas the solid bed temperature was initialized at 873 K. Table 4 provides the effect of inlet gas temperature on silicon production. It can be seen that the silicon production rate is increased with feed gas temperature. Figure 18 shows the average radial profiles of gas temperature, solid temperature, fines production rate, and solid silicon production rate in the dense area at 15 s of operation. The gas temperature profiles show that the core zone is mostly affected by the inlet gas temperature, as seen in Figure 18a. When the feed temperature is 300 K, the gas temperature in the core zone, where the reaction rate can be very low, is low. The lowest feed gas temperature gives the largest cooling effect of feed gas. However, the solid temperature in the core zone is less affected by the cool gas feeding, as it does not change much with inlet temperature despite the very low gas feed temperature of 300 K, as seen in Figure 18b, because the initial solid particle temperature is relatively high (850 K). The cool gas requires a long time to cool down the solid temperature. Both gas and solid temperatures in the annular zone changes less with the feed temperature (see Figures 18a,b) because the dense zone with 60−65% of total volume is close to the wall and is therefore effective in maintaining high solid temperature. This is the advantage of using a spouted bed, rather than another kind of fluidized bed, as it here promotes the heterogeneous reaction. Figure 18c shows that the silicon fines production rate increases with the feed gas temperature. At very low temperature, 300 K, the core zone gas temperature is too low, that is, below 500 K, and the reaction sharply drops and approaches zero. The silicon solid production rate also increases with the inlet gas temperature, as seen in Figure 18d. However, in the core zone, the silicon solid production rate slightly increases with the feed temperature due to the presence of dilute solid concentration and becomes more different in the annular zone due to high solid content and temperature in that region. The axial profiles of averaged gas temperature covering both solid bed and freeboard are shown in Figure 19a. Meanwhile, Figure 19b shows the solid temperature as a function of reactor

Figure 17. Particle size distribution (a) and histogram of particle size distributions (b) at various times.

coupling model can predict the size of individual particles in the spouted bed reactor. The particle diameter increases over the time via chemical vapor deposition of silane. Figure 17a shows that the mixture contains different sizes of particles with good mixing. The averaged particle growth over time of operation shows a growth rate about 22.14 nm/s (see Table 4). Figure 17b shows the histogram of particle size change distributions at 10, 15, and 60 s. The particle size distribution is not uniform because of the nonuniformity of gas, solid movement, concentration, and temperature in the system. Moreover, it can be noted that the standard deviation of particle size distribution at 10, 15, and 60 s are 7.6, 13.7, and 36.0 nm, respectively. This means that the J


Article


Figure 18. Radial profiles of averaged gas temperature (a), solid temperature (b), fines production rate (c), and solid silicon production rate (d) for various inlet gas temperatures.

Figure 19. Axial profiles of averaged gas temperature (a), solid temperature (b), fines production rate (c), and solid silicon production rate (d) for various inlet gas temperatures. K


Article

Industrial & Engineering Chemistry Research bed height in the dense zone. In the dense bed area, small temperature differences for various feed temperatures were found because the high wall temperature effect is dominant. In addition, the heat conduction through the particle solid bed from the wall is large, making the effect of feed temperature on the gas and solid temperatures in the dense bed zone to be less significant. On the other hand, the gas temperature in freeboard area drops significantly when the feed temperature is too low (at 300 K). With a low feed gas temperature of 300 K, less fines formation occurs. The cooling effect from the feed reduces the freeboard temperature significantly compared to those at higher feed temperatures (e.g., 473 and 673 K). The axial position of silicon fines production shows that the inlet gas temperature has an effect on the fines production rates, as seen in Figure 19c. The results show that, in the particle dense bed zone near the feed zone, the fines formation is low due to there being a low gas volume fraction. In the freeboard area, in contrast, results show that the fines are generated more. It can be seen that when the gas temperature is too low, at 300 K, the fines formation is the lowest. This condition provided a cooling of feed gas through the freeboard area, which helps to reduce the temperature in this zone and prevents the fines generation. Low temperature in the freeboard can suppress the fines formation, which is an undesired reaction. Figure 19d shows the axial profile of silicon solid product rate in the dense bed zone. It can be seen that the silicon solid production rates increase with the feed temperature. Thus, in actual operation, the feed temperature should be properly chosen in order to suppress fines formation and enhances silicon solid production. Table 4 provides the effect of inlet gas temperature on particle growth rate as well as standard deviation of the particle size distribution. The results show that the growth rate increases with inlet gas temperature, leading to the particle size, too, increasing with the feed temperature. With higher feed temperature the particle size distribution is narrower, as shown in Figure 20a. In contrast, at low feed temperature (300 K), the cooling effect is dominant, and nonuniformity of solid temperature in the bed was found. The temperature at the central zone is much different from that at the wall and from the initial solid bed temperature, leading to yet more nonuniformity. As a result of all this, the local production rate is not uniform, and this ends in a broader distribution of particle size. In the case of 300 K operation, the standard deviation is 32.3 nm (see Table 4). On the other hand, at higher feed temperature (i.e., 673 K), the temperature of the solid phase is more uniform, and consequently the particle size distribution is narrower and has a lower standard deviation (i.e., 7.7 nm, see Table 4). 6.6. Effect of Inlet Gas Composition. In this section, the effect of inlet gas composition on silicon production is discussed. Throughout this discussion, the inlet gas temperature is 473 K and inlet gas velocity is 1.434 m/s. The inlet silane composition, however, is varied from 20 to 40%, with nitrogen as a carrier gas. Figure 21 shows the gas and solid temperature profiles as a function of axial positions for various inlet gas compositions. It can be observed that the gas and solid temperatures in the dense area were not affected by different gas compositions. The gas temperature changes significantly with the silane composition only in the freeboard zone (z/Z > 0.3), as seen in Figure 21a. Because the reaction of silane decomposition is exothermic, the gas temperature increases with bed height, leading to very high gas temperature in the freeboard. Indeed, the temperature in the freeboard is above 850 K in all cases. The reaction rate rises up suddenly when the bed temperature is higher than 850 K. This

Figure 20. Particle size distribution at 15 s of operation for various inlet gas temperatures (a), inlet gas compositions (b), and superficial gas velocities (c).

Figure 21. Axial profiles of averaged gas temperature (a) and solid temperature (b) for various inlet gas compositions.

causes a very high reaction rate in the freeboard area, where the gas temperature is very high. The reaction rate depends on the silane concentration causing difference in heat generation. Thus, the effect of silane composition in the freeboard, where the temperature is above 850 K, on the gas temperature is significant. The operation suggestion is that the wall temperature in the freeboard zone should be cooled down to prevent too high a gas temperature in this zone by suppressing fines formation. However, the temperature in the dense zone does not affect by the inlet gas composition because the reaction rate at below 850 K slightly increases with temperature. L


Article


Figure 22. Radial profiles of averaged silicon solid (a) and silicon fine (b) production rates for various inlet gas compositions.

Figure 23. Axial profiles of averaged silicon solid (a) and silicon fines (b) production rates for various inlet gas compositions.

The silicon production rates are strongly affected by inlet gas composition. Increasing the fraction of inlet reactant gas enhanced both silicon solid and fines production rates, as seen in Table 4. Figure 22 shows the radial profiles of both silicon production rates for various inlet gas compositions. The results show that both solid silicon and fines productions in the central region are low, and the effect of gas composition is minimal. This is due to low gas residence time, which causes low conversion in this zone. On the other hand, the solid silicon and fines products increase with reactant composition in the annular zone. Figure 23 shows the silicon production rates as a function of axial position for various inlet silane compositions. The solid silicon production rate increases with inlet silane composition as seen in Figure 23a. Figure 23b shows that the silicon fines production rate significantly increases with silane composition in the freeboard region. The selectivity, defined as the mole of silicon solid product per the mole of total silicon product, for inlet reactant composition of 20%, 30%, and 40% in the dense area are 72.5%, 69.1%, and 66.3%, respectively. They show that the selectivity of silicon solid production is decreased with inlet silane composition. The rate of homogeneous reaction is more sensitive to the reactant composition than that of heterogeneous reaction. Thus, increasing reactant composition increases the homogeneous rate, leading to lower selectivity of solid silicon product. The effect of inlet silane composition on the particle growth rate, average particle size, and standard deviation of particle size distribution are shown in Table 4. As silane composition is increased, the particle growth rate increases. Figure 20b shows that the particle size distribution is broader as the inlet gas composition is increased, which is confirmed by the standard deviation increasing from 13.7 to 28.1 nm as the inlet silane

composition is increased from 20% to 40%. This is due to there occurring more nonuniformity of reaction rate for a higher inlet composition (see Figure 22). 6.7. Effect of Gas Superficial Velocity. In this section, the effect of gas superficial velocity on silicon production is discussed. Throughout this section, the inlet gas composition is 20% and inlet gas temperature is 473 K. To study the impact of hydrodynamics on the reactor performance, the weight of solid particles per feed flow rate (weight−time) is fixed as the gas feed velocity is varied. To keep the weight−time constant at 46.86 s, the particle loading is increased from 4.88 kg (500000 particles) to 6.83 kg (700000 particles) while the superficial gas velocity is increased from 1.30Ums to 1.63Ums. Figure 24 shows the distributions of solid temperature, gas temperature, silane concentration, silicon fines concentration, and particle size distribution for various superficial gas velocities. Results were monitored at 15 s of operation time. At very low gas velocity, the solid loading is low to keep constant weight time. Thus, the solid particles are easily spouted. This in turn causes lower solid holdup in the annular zone, as shown in Figure 25. On the other hand, using a higher gas velocity leads to higher solid holdup (see the case of 1.52Ums), as higher solid loading is required to keep constant weight time. However, when the feed velocity is very high (1.63Ums), the effect of solid loading is more dominant, leading to a less expanded bed than that in the case of 1.52Ums. It thus happens that the solid holdup for 1.63Ums is in between those of 1.52 and 1.30Ums. At low feed velocity (1.30Ums), the solid holdup is low, leading to low solid silicon production rate (see Figure 26e). In this case, the cooling effect is also low, and the solid and gas temperatures see less of a cooling down effect (see Figure 26c). Hence, the fines formation is higher due to there being a high gas temperature, as seen in Figures 26c,f, and 27b. M


Article


effect of medium solid temperature and solid holdup, the solid silicon production is the highest of all studied cases, as seen in Figures 26 and 27. Thus, in practice, the feed velocity and the solid loading should be chosen to get optimum solid silicon production. The effect of inlet superficial gas velocity on the particle growth rate and standard deviation of particle size distribution are shown in Table 4. The particle growth rate is less significantly changed with gas velocity. Figure 20c shows the particle size distribution for various superficial gas velocities. The results show that low superficial gas velocity (1.30Ums) gives a narrow distribution. Increasing gas velocity to 1.52Ums, however, gives a broader distribution. At a very high feed velocity, the particle size distribution is narrower. The distribution should correspond to the solid holdup. Higher solid holdup tends to have more nonuniformity. At medium feed velocity, the solid holdup is the highest, thus the distribution is broader, with a standard deviation of 22.6 nm. In contrast, at very low velocity, 1.30Ums, the solid holdup is very low and thus the standard deviation is low also (13.7 nm).

7. CONCLUSIONS In this paper, we have examined the application of an Euler− Lagrangian model to assess the hydrodynamics and performance of silicon pyrolysis in a spouted bed. The discrete element method was used for tracking the motion of individual particles, while the gas phase was treated as a continuum. This coupling model was used in the simulation of the reactor, which included both heterogeneous and homogeneous reactions of silane decomposition, as well as heat transport effects. Hydrodynamics feature such as solid holdup profile, particle circulation, and gas distribution were shown. The results indicate that the coupled model was able to reasonably describe particle growth in silane decomposition. The spouted bed reactor induced the heterogeneous reaction in the annular zone by allowing larger solid holdup there, and the smaller gas residence time in the central region reduced the fines formation in the system. Verification of the solid movement results with the results from experiments and correlations show good agreement. In addition, the axial profile of bed voidage agrees very well with correlation. In this work, effects of three parameters were studied in detail: inlet gas temperature, inlet gas composition, and superficial gas velocity. The main conclusions are given as follows. (1) Inlet gas temperature does not significantly affect the growth rate in a spouted bed as keeping high wall temperature, but lowering of feed gas temperature is beneficial in terms of suppressing fines formation. However, if the inlet gas temperature is too low, the gas temperature is much different from that at the wall, and the solid bed temperature in turn causes a broad particle size distribution. (2) Increasing inlet gas composition enhanced both solid silicon and fines production rate. In addition, it significantly increased fines formation in freeboard due to there being a very high temperature. Thus, in practice, the wall temperature in the freeboard region should be cooled down to prevent a too high gas temperature and in order to suppress the fines generation. This suggests that the wall temperature should be profiled in an optimum manner and could be the focus of further study. (3) The feed velocity should be optimum to spout the solid particles and to provide a cooling effect to the feed gas to

Figure 24. Distributions of solid temperature (a), gas temperature (b), silane concentration (c), silicon fines concentration (d), and particle size (e) for various superficial gas velocities.

Figure 25. Radial profiles of averaged solid holdup for various superficial gas velocities.

At a higher feed velocity (1.52Ums), the solid loading is required to be higher so as to keep a constant weight−time, leading in turn to higher solid holdup. However, at higher feed velocity, the solid temperature is cooled down to a greater extent, leading lower gas and solid temperature, as seen in Figures 26a−d. In the production of silicon solid, the solid holdup is dominant, while in the fines production, the effect of gas temperature is important. At this condition, higher solid holdup leads to higher production of solid silicon and lower gas temperature leads to lower fines production as seen in Figure 27a,b. At a very high feed velocity (1.63Ums), the solid loading is also very high. The bed expansion is therefore not very high because of higher loading. Thus, it happens that solid holdup for this case is slightly lower than that for the case of 1.52Ums but higher than 1.30Ums. Here, the gas temperature is similar to that in 1.52Ums, as seen in Figure 26c, while the solid temperature is not cooled down much because of the high loading of solid. Thus, the solid temperature is still hot and same as that in 1.52Ums. With the N


Article


Figure 26. Averaged profiles for various superficial gas velocities: solid temperature as a function of radial (a) and axial position (b), gas temperature as a function of radial (c) and axial position (d), radial profiles of solid silicon production (e) and fines production (f).

Figure 27. Averaged axial profiles of silicon solid (a) and silicon fines (b) production rates for various superficial gas velocities.

■

suppress the fines formation. Increasing superficial gas velocity at a constant weight time leads to decreased fines formation and increased solid silicon product. However, increasing the gas velocity the broader solid distribution is obtained.

AUTHOR INFORMATION

Corresponding Author

*Tel.: +66-2797-0999 ext 1210. Fax: +66-2561-4621. E-mail: [email protected]. O


Article

Industrial & Engineering Chemistry Research Uaz Usz v Vreactor z

ORCID

Sunun Limtrakul: 0000-0002-3340-7240 Notes

The authors declare no competing financial interest.

■

superficial gas velocity in annulus region (m/s) superficial gas velocity in spouted region (m/s) solid velocity (m/s) volume of reactor (m3) axial position (m)

Greek Letters

ACKNOWLEDGMENTS This work was financially supported by the Kasetsart University Research and Development Institute (KURDI), the Faculty of Engineering of Kasetsart University, the Centre of Excellence on Petrochemical and Materials Technology (PETROMAT), and the Centre for Advanced Studies in Industrial Technology. We are pleased to acknowledge Prof. Dr. Yutaka Tsuji and Assistant Prof. Dr. Toshihiro Kawaguchi for their collaboration.

α⃗ β δ⃗ ε ηdamp λg μf ρ ϕ ω⃗

■

NOMENCLATURE a⃗ acceleration vector (m/s2) A cross sectional area (m2) Ci Concentration of gas species i (kmol/m3) CD drag coefficient (−) Ĉ pg heat capacity of gas mixture (kJ/kg·K) Cps heat capacity of solid (kJ/kg·K) dp particle diameter (m) d̅pf averaged fines particle diameter (m) ds spouted diameter (m) Dc column diameter (m) Di inlet (orifice) diameter (m) Di,m diffusivity coefficient of species i in the gas mixture (m2/s) E single large particle collection coefficient (−) F⃗ g→s drag force of fluid on particles (N) F⃗ total total force acting on the particle (N) F⃗ C contact force (N) F⃗ D fluid force (N) g⃗ gravity acceleration (m/s2) h heat transfer coefficient (W/m2·K) Hm maximum spoutable bed height (m) I moment of inertia of the particle (kg·m2) kB Boltzmann constant, 1.38066 × 10−23 J/(molecule·K) khet kinetics rate constant of heterogeneous reaction (m3 reactor/m2 particle·s) khom kinetics rate constant of homogeneous reaction (1/s) ksc scavenging coefficient (1/s) kstif stiffness of spring (N/m) ṁ gs interphase mass transfer from gas to solid (kg/m3s) mp mass of particle (kg) Msi mass of silicon product (kg) Mwi molecular weight of species I (g/mol) p pressure (Pa) Pe Peclet number QRxn,het heat of reaction of heterogeneous reaction (kJ/mol) QRxn,hom heat of reaction of homogeneous reaction (kJ/mol) r radial position (m) Ri reaction rate of species i (kmol/m3·s) Re Reynolds number (−) t time (s) Tg gas temperature (K) solid temperature (K) Ts T⃗ torque (N.m) u fluid velocity (m/s) Umf minimum fluidization velocity (m/s) Ums minimum spouted velocity (m/s) U superficial velocity (m/s)

angular acceleration vector (rad/s2) momentum transfer coefficient (kg/s·m3) particle displacement (m) void fraction (−) damping coefficient (N·s/m) gas thermal conductivity (kW/m·K) fluid viscosity (kg/m·s) density (kg/m3) shape factor of particle Angular momentum (rad/s)

Subscripts

b f g mf nor p s Si tan

■

Bulk phase friction Gas phase Minimum fluidized bed Normal direction particle solid phase silicon tangential direction

REFERENCES

(1) Bamett, A. M.; Wenger, H. Solar Electric Power in the Mainstream. 2003 IEEE Power Engineering Society General Meeting (IEEE Cat. No. 03CH37491), Toronto, July 13-17, 2003; IEEE, 2003; Vol. 3, pp 1330−1333, DOI: 10.1109/PES.2003.1267342. (2) del Coso, G.; Tobias, I.; Canizo, C.; Luque, A. Temperature Homogeneity of Polysilicon Rods in a Siemens Reactor. J. Cryst. Growth 2007, 299 (1), 165. (3) Ramos, A.; Filtvedt, W. O.; Lindholm, D.; Ramachandran, P. A.; Rodríguez, A.; Del Cañizo, C. Deposition Reactors for Solar Grade Silicon: A Comparative Thermal Analysis of a Siemens Reactor and a Fluidized Bed Reactor. J. Cryst. Growth 2015, 431, 1. (4) Lai, S.; Dudukovic’, M. P.; Ramachandran, P. A. Chemical Vapor Deposition and Homogeneous Nucleation in Fluidized Bed Reactors: Silicon from Silane. Chem. Eng. Sci. 1986, 41 (4), 633. (5) Pina, J.; Bucala, V.; Schbib, N. S.; Ege, P.; de Lasa, H. I. Modeling a Silicon CVD Spouted Bed Pilot Plant Reactor. Int. J. Chem. React. Eng. 2006, 4, 1306. (6) Balaji, S.; Du, J.; White, C. M.; Ydstie, B. E. Multi-Scale Modeling and Control of Fluidized Beds for the Production of Solar Grade Silicon. Powder Technol. 2010, 199, 23−31. (7) Cadoret, L.; Reuge, N.; Pannala, S.; Syamlal, M.; Coufort, C.; Caussat, B. Silicon CVD on Powders in Fluidized Bed: Experimental and Multifluid Eulerian Modelling Study. Surf. Coat. Technol. 2007, 201, 8919−8923. (8) Liu, S. S.; Xiao, W. De. Numerical Simulations of Particle Growth in a Silicon-CVD Fluidized Bed Reactor via a CFD-PBM Coupled Model. Chem. Eng. Sci. 2014, 111, 112. (9) Reuge, N.; Cadoret, L.; Caussat, B. Multifluid Eulerian Modelling of a Silicon Fluidized Bed Chemical Vapor Deposition Process: Analysis of Various Kinetic Models. Chem. Eng. J. 2009, 148 (2−3), 506. (10) White, C. M.; Ege, P.; Erik Ydstie, B. Size Distribution Modeling for Fluidized Bed Solar-Grade Silicon Production. Powder Technol. 2006, 163 (1−2), 51. (11) Cundall, P. A.; Strack, O. D. L. A Discrete Numerical Model for Granular Assemblies. Geotechnique 1979, 29, 47. (12) Limtrakul, S.; Chalermwattanatai, A.; Unggurawirote, K.; Tsuji, Y.; Kawaguchi, T.; Tanthapanichakoon, W. Discrete Particle P


Article

Industrial & Engineering Chemistry Research Simulation of Solids Motion in a Gas-Solid Fluidized Bed. Chem. Eng. Sci. 2003, 58, 915−921. (13) Limtrakul, S.; Rotjanavijit, W.; Vatanatham, T. Lagrangian Modeling and Simulation of Effect of Vibration on Cohesive Particle Movement in a Fluidized Bed. Chem. Eng. Sci. 2007, 62, 232−245. (14) Limtrakul, S.; Thanomboon, N.; Vatanatham, T.; Khongprom, P. DEM Modeling and Simulation of a down-Flow Circulating Fluidzed Bed. Chem. Eng. Commun. 2008, 195, 1328−1344. (15) Tsuji, Y.; Kawaguchi, T.; Tanaka, T. Discrete Particle Simulation of Two-Dimensional Fluidized Bed. Powder Technol. 1993, 77, 79. (16) Tsuji, Y.; Tanaka, T.; Ishida, T. Langrangian Numerical Simulation of Plug Flow of Cohesionless Particles in Horizontal Pipe. Powder Technol. 1992, 71, 239. (17) Kaneko, Y.; Shiojima, T.; Horio, M. DEM Simulation of Fluidized Beds for Gas-Phase Olefin Polymerization. Chem. Eng. Sci. 1999, 54, 5809. (18) Li, Z.; Janssen, T. C. E.; Buist, K. A.; Deen, N. G.; van Sint Annaland, M.; Kuipers, J. A. M. Experimental and Simulation Study of Heat Transfer in Fluidized Beds with Heat Production. Chem. Eng. J. 2017, 317, 242. (19) Limtrakul, S.; Boonsrirat, A.; Vatanatham, T. DEM Modeling and Simulation of a Catalytic Gas-Solid Fluidized Bed Reactor: A Spouted Bed as a Case Study. Chem. Eng. Sci. 2004, 59, 5225−5231. (20) Zhuang, Y. Q.; Chen, X. M.; Luo, Z. H.; Xiao, J. CFD-DEM Modeling of Gas-Solid Flow and Catalytic MTO Reaction in a Fluidized Bed Reactor. Comput. Chem. Eng. 2014, 60, 1. (21) Passalacqua, A.; Fox, R. O.; Garg, R.; Subramaniam, S. A Fully Coupled Quadrature-Based Moment Method for Dilute to Moderately Dilute Fluid-Particle Flows. Chem. Eng. Sci. 2010, 65 (7), 2267. (22) Ergun, S. Fluid Flow through Packed Columns. Chem. Eng. Prog. 1952, 48, 89−94. (23) Hogness, T. R.; Wilson, T. L.; Johnson, W. C. Thermal Decomposition of Silane. J. Am. Chem. Soc. 1936, 58, 108. (24) Iya, S. K.; Flagella, R. N.; DiPaolo, F. S. Heterogeneous Decomposition of Silane in a Fixed Bed Reactor. J. Electrochem. Soc. 1982, 129, 1531−1535. (25) Patankar, S. V. Numerical Heat Transfer and Fluid Flow, 1st ed.; Hemisphere Series on Computational Methods in Mechanics and Thermal Science; Taylor & Francis, 1980. (26) Mathur, K. B.; Gishler, P. E. A Technique for Contacting Gases with Coarse Solid Particle. AIChE J. 1955, 1, 157. (27) Malek, M. A.; Lu, B. C. Y. Pressure Drop and Spouted Bed Height in Spouted Beds. Ind. Eng. Chem. Process Des. Dev. 1965, 4 (1), 123. (28) Smith, K. J.; Arkun, Y.; Littman, H. Studies on Modeling and Control of Spouted Bed reactorsI: Reactor Modeling. Chem. Eng. Sci. 1982, 37 (4), 567. (29) Ren, B.; Zhong, W.; Jin, B.; Yuan, Z.; Lu, Y. Computational Fluid Dynamics (CFD)-Discrete Element Method (DEM) Simulation of Gas-Solid Turbulent Flow in a Cylindrical Spouted Bed with a Conical Base. Energy Fuels 2011, 25 (9), 4095. (30) McNab, G. S. The Pattern of Change of Spout Diameter in a Spouting Bed. Br. Chem. Eng. Process Technol. 1972, 17, 532. (31) Mathur, K. B.; Lim, C. J. Vapour Phase Chemical Reaction in Spouted Beds: A Theoretical Model. Chem. Eng. Sci. 1974, 29, 789. (32) Mamuro, T.; Hattori, H. Flow Pattern of Fluid in Spouted Beds. J. Chem. Eng. Jpn. 1968, 1, 1. (33) Grbavcic, Z. B.; Vukovic, D. C.; Zdanski, F. K.; Littman, H. Fluid Flow Pattern, Minimum Spouting Velocity and Pressure Drop in Spouted Beds. Can. J. Chem. Eng. 1976, 54, 33.

Q


Modeling of chemical vapor deposition of silane for silicon production

Recommend Documents