PARTICLE TECHNOLOGY AND FLUIDIZATION Modeling of an Aerosol Reactor for Optimizing Product Properties Eric Bain Wasmund Inco Special Products, Mississauga, ON, Canada L5K 2L3 Shadi Saberi Inco Technical Services Limited, Mississauga, ON, Canada L5K 1Z9 Kenneth S. Coley Dept. of Materials Science and Engineering, McMaster University, Hamilton, ON, Canada L8S 4L7 DOI 10.1002/aic.11181 Published online April 16, 2007 in Wiley InterScience (www.interscience.wiley.com). Important and unique metal powders are made industrially by a variety of vapor condensation processes in tube reactors. Often, however, the fundamental mechanisms for particle formation and growth are still not well understood. In this article, a computational fluid dynamics (CFD) model was developed to examine a tube reactor’s internal flow characteristics. The model identified a massive zone of fluid recirculation in the top half of the reactor. In-situ sampling from an experimental reactor under the same conditions revealed a large increase of aerosol particle size corresponding to the region of recirculation. A first principles mass balance model based on chemical kinetics and aerosol physics was developed for this system which showed that the average particle size grew monotonically with time in the reactor. On the basis of this firmly established link between residence time and particle size, a new reactor geometry was proposed to produce a ‘‘plug-flow’’ velocity profile with a narrower particle size distribution. A CFD model was used to prototype the new configuration, and then this new reactor design was tested experimentally to confirm that the design objective was achieved. This work shows the potential synergies between first principles models for process understanding and CFD models for process prototyping and optimization. Ó 2007 American Institute of Chemical Engineers AIChE J, 53: 1429–1440, 2007 Keywords: aerosols, particle technology, reactor analysis, computational fluid dynamics Introduction Aerosol tube reactors for making nickel and iron powders by the thermal decomposition of their metal carbonyls have been used commercially for more than 80 years. One of the pioneers in this field, Alwyn Mittasch originally disclosed Correspondence concerning this article should be addressed to E. B. Wasmund at ewasmund@inco.com. Ó 2007 American Institute of Chemical Engineers AIChE Journal June 2007 the process to the world in German patent 500,692.1 Today, many thousands of tons of specialty powders are made by companies such as Inco Limited, International Specialty Products and BASF. The technique is very simple and flexible. On the basis of empirical study, practitioners have learned to control the morphology, purity, and particle size so that these products can find their way into diverse applications such as magneto-rheological fluids, powder metallurgy, EMI shielding, catalysts, sintered battery plaques, and binders for hard metal tools. Until recently, the system has been Vol. 53, No. 6 1429 Figure 1. Cross-section of the reactor tube configuration for the Case A scenario. [Color figure can be viewed in the online issue, which is available at www.interscience.wiley.com.] treated as a ‘‘black box’’. Inputs such as flowrates and composition, dopants and wall temperatures have been empirically correlated with the attributes of the resulting powder, without understanding or proving the fundamental mechanisms that are in effect. This approach is suboptimal, since it does not allow the application of engineering first principles such as reactor design and fluid dynamics. In this article, we will show how computational fluid dynamics (CFD) and chemical kinetics can be applied to derive insights that allow for the optimization of an aerosol tube reactor, in particular, for making powders with a narrow size distribution. Experimental The model system that we are investigating is a 5-cm diameter, 30-cm long hot-wall tube reactor for decomposing nickel carbonyl vapor carried in CO. The decomposition reaction for nickel carbonyl is shown as Eq. 1. Nickel carbonyl becomes thermally unstable at temperatures greater than about 1208C at standard conditions. It can decompose directly to form fine nuclei particles, or it can react heterogeneously onto surfaces. Once particles have been formed, they can continue to grow by coagulation and reorganize by sintering. Each of these phenomena can in principle occur simultaneously, which challenges our understanding of how to optimize the reactor design. NiðCOÞ4 $
NiðsÞ þ 4 CO DG ðkJ=molÞ ¼ 155:20  0:3964 TðKÞ DH
ðkJ=molÞ ¼ 140 (1) This reactor, its operating conditions, and the particle sampling technique, have been described in detail elsewhere2,3 and will be briefly reviewed here, please refer to Figure 1 which illustrates the reactor in cross-section. The outside wall of the reactor is heated with five independently controlled 500 W resistance heating zones, which are evenly spaced along the reactor axis as shown in Figure 1. Each 1430 DOI 10.1002/aic heating zone is 50 mm along the axis. The power in each zone is the manipulated variable for controlling the outside wall temperature, which is measured with a contact thermocouple at the center of each zone. The top of the first heating zone is ~22 mm below the reactor feed-gas inlet. Feed-gas is fed into the top of the reactor through a tube nozzle. Two reactor configurations will be discussed here, which are described in Table 1. For Case A, the reactor outside wall temperature is controlled between 575 and 6508C. The total inlet flow-rate is 3 standard liters per minute (slpm), consisting of 0.15 slpm of Ni(CO)4 and 330 ppm of NH3, the balance being CO. The reactor configuration for Case A is as shown in Figure 1. For Case B, the reactor outside wall temperature is controlled between 575 and 7258C, and the inlet flow-rate is 18 slpm, consisting of 0.38 slpm of Ni(CO)4 and 330 ppm of NH3, the balance being CO. For reasons that will be evident later, the reactor configuration is modified by replacing the existing straight inlet nozzle with a slowly expanding cone nozzle that opens up to one half of the reactor diameter. This is the so-called ‘‘well fared entry’’ that will be described later. The reactor cross-section for the Case Table 1. Flow, Temperature, and Energy Values for the Case A and Case B Scenarios Parameter Case A Case B Measured outside wall temps for each zone (8C) Measured center-line temps for each zone (8C) Measured flowrate, (slpm) Nozzle diameter (mm) Measured Ni(CO)4 (v/v%) Calculated energy requirements for heating and reaction (Watts) Calculated contribution of reaction to total (%) Measured actual heat consumed from furnaces (Watts) 625/575/575/ 575/650 350/410/430/ 460/450 3 3.1 4.8 725/575/575/ 575/650 330/360/360/ 370/340 18 22.0 2.1 Published on behalf of the AIChE 50 210 25 16 1110 1350 June 2007 Vol. 53, No. 6 AIChE Journal Figure 2. Cross-section of the reactor tube configuration for the Case B scenario. In this case, the inlet nozzle is a slowly expanding tube, designed to create the ‘‘well-fared’’ entry conditions for minimizing recirculation. [Color figure can be viewed in the online issue, which is available at www.interscience.wiley.com.] B experiments is shown as Figure 2. Finally, one experiment is run where the conditions of Case B are repeated, except that the reactor is turned upside down and fed from the bottom, all other conditions being the same. This experiment will be designated as Case B up-flow. After the reactor, the gas consists of CO and CO2, some unreacted Ni(CO)4 and nickel aerosol particles. This gas is quenched with 90 slpm of nitrogen, passes through a fabric filter and the composition is analyzed by an on-line Fourier Transform Infra-red spectrometer. The gas composition of inlet and outlet for each experiment is shown in Table 2. Forty-four grams of powder are collected on the gas filter from the Case A experiment, 61 g from the Case B experiment, and 61 g from the Case B up-flow experiment. These samples of the final product are so-called ‘‘integral results’’ since they represent the summation of all powders that have gone through the reactor irrespective of path. In particular, the volume size distribution is measured by laser light scattering using a Malvern Mastersizer 2000. The details are described in Appendix A. The tube reactor shown in Figures 1 and 2 is equipped with five water-cooled thermophoretic sampling ports (shown on the right-hand side of Figures 1 and 2). Each sampling port is located at the center of each of the five heating zones. To gain an experimental approximation of the internal temperatures at the center-line of the reactor for each of the experimental conditions, the following method is used. The experimental flow-rate and outside wall temperatures as shown in Table 1 are replicated for a flow of nitrogen, which has nearly identical thermal and transport properties as CO. At steady state, a 1.5-mm unshielded Type K thermocouple AIChE Journal June 2007 Vol. 53, No. 6 junction is passed through the first sampling port and the center-line gas temperature is recorded. This method is repeated sequentially to obtain center-line temperatures for the axial coordinates corresponding to the next four samplers. Because the thermocouple is unshielded, the junction may absorb radiation from the tube walls and could represent a higher temperature than the surrounding gas. This makes the thermocouple measurement more representative of the temperature of a radiation-absorbing solid, such as the nickel aerosol particles. This measurement also does not take into consideration the cooling effect of the nickel carbonyl reaction that would be present during the actual experiment. Later in this article, this system is modeled with and without consideration of the effect of the reaction on the energy balance to show the sensitivity of the reaction on the temperature profiles in the reactor. The primary purpose of the thermophoretic sampling ports is to obtain spatially resolved ‘‘snapshots’’ of the aerosol particle growth inside the reactor. The use of thermophoresis to Table 2. Stream Compositions (Measured and Calculated) for Case A and Case B Measured Gas Compositions (v/v%) Case A Case B Inlet Ni(CO)4 (measured) Inlet CO (from mass balance) Inlet N2 (from mass balance) Outlet Ni(CO)4 (measured) Outlet CO (from mass balance) Outlet N2 (from mass balance) Outlet CO2 (measured) 4.8 92 3.1 0.02 97 2.8 0.15 2.1 97 0.6 0.02 >99 0.6 0.07 Published on behalf of the AIChE DOI 10.1002/aic 1431 Figure 3. Samples obtained from the micro-grids attached to the second sampler under Case A conditions. 3A and 3B are both from the middle of the tube, but measured at 5000 and 60,000 respectively. 3C is extracted from the half-way position (middle) of the tube and 3D is extracted from near the wall, both at 5000. Figures 3A, C, D are measured using SEM and the background is the copper grid. Figure 3B is measured using STEM and the background is carbon film. obtain unbiased aerosol particle samples has been pioneered by Dobbins and Megaridis4 and the method used in this study has been explained in greater detail elsewhere.2,3 The main advantage of thermophoresis is that the aerosol collection mechanism is relatively insensitive to particle size, when compared with Brownian diffusion, which is strongly temperature dependent. To obtain thermophoretic samples from inside the reactor, the following technique is used. A sampler with six TEM micro-grids is located inside of each sampling arm prior to the experiment. The TEM micro-grids are 3-mm diameter copper grids with deposited holey carbon (manufactured by SPI Supplies Division of Structure Probe, West Chester, PA). Throughout the experiment, a small 20 ml/min bleed of nitrogen flows through each sampling arm into the reactor to maintain isolation of the sampling arm, except during deployment. During the experiment, each sampler is sequentially injected into the gas flow for a controlled period of time, always less than about 1 s. During deployment, the TEM micro-grids are parallel to the axial flow of gas, and are positioned such that they obtain aerosol samples that correspond to the inside wall, half-way (middle) position and center-line of the reactor (2 TEM micro-grids for each radial sampling position, which allows for comparison of each sam1432 DOI 10.1002/aic ple with a replicate). Calculations have been performed to show that the thermophoretic velocity induced by these conditions is at least 10 times the velocity induced by Brownian diffusion.2 In this study, we are focusing on the variation in particle growth as a function of the radial position inside the tube. The TEM microgrids from each of the 15 spatial locations inside the reactor have been examined for the Case A experiment and are reported elsewhere.2 As an example of the strong variation of particle size with radial position at some positions inside the reactor, Figure 3 shows example micrographs from the second sampler at the center-line (3A, B), at the half-way (middle) position of the tube (3C), and at the tube wall (3D). Figures 3A, C, D are obtained by a scanning electron microscope (JEOL 6400) at 5000 magnification. In these images, the particles are on a background of the copper grid. Figure 3B is obtained from the same field used to generate Figure 3A, using a scanning transmission electron microscope at 60,000 (JEOL 2010F). In this case, the particles are on a background of the holey carbon, to allow for electron beam transmission. The micro-grids represented in Figures 3C, D have also been examined using the STEM, but small particles as shown in Figure 3B were not observed. These micrographs suggest that at the center-line, the nickel Published on behalf of the AIChE June 2007 Vol. 53, No. 6 AIChE Journal Figure 4. Cumulative size distributions from the second sampler under the Case A conditions. At the center-line the size distribution is bimodal and the d50 is 0.5 mm (based on 8910 objects). Samples at the halfway (middle) position have a d50 of 0.6 mm (based on 1030 objects) while samples near the wall have a d50 of 0.9 mm (based on 1190 objects). This suggests that the product size distribution is a strong function of radial position inside the reactor. particles have a bimodal distribution, and that the particle size distribution, at least qualitatively, is larger near the wall. To convert these spatially resolved microgrid samples into particle size distributions, SEM images are obtained from 10 representative and repeatable sectors of the grid. These images are taken at 2000 or 5000. In the case where the presence of finer particles has been identified, the imaging is performed at 60,000 using the STEM. Each of the 10 SEMs is run through an image analysis software program (Image Pro Plus Version 4.5 by Media Cybernetics, Silver Spring Md, USA) to calculate and record the mean diameter for each object. A roundness metric for each object is also calculated, which quantifies the object’s nearness to the shape of a circle. Through an analysis of many micrographs, it was concluded that all particles less than about 0.7 mm are roughly spherical (having a roundness metric less than 2.5) and their volume can be adequately described by the standard formula for a sphere (p/6 d3) where d is the mean diameter. Larger particles more closely resemble filaments or agglomerates. For these nonround particles, many particles were studied and the relationship between mean diameter and volume was estimated as shown in Eq. 2. The details for this estimation are shown in Appendix B. Under the Case A conditions, more than half of all particle studied are nonspherical, following the criterion above, after the first sampler position.2 v½m3  ffi 1:4  104 m0:5 d½m2:5 (2) The list of mean diameters for all objects obtained from each grouping of SEMs can be used to calculate number-based size distributions corresponding to each position in the reactor that was sampled. Using the scaling law of V ¼ p/6 d3 for particles with a mean diameter less than 0.7 mm and Eq. 2 for particles with a mean diameter greater than 0.7 mm, the number distributions are converted to volume-based cumulative size distributions. It should also be noted that at the upper end of these size distributions, there may be a reduced representation of the larger particles in the overall sample because of decreased thermophoretic mobility. To address this, Talbot’s equation for estimating the thermophoretic veAIChE Journal June 2007 Vol. 53, No. 6 locity is used to calculate a size dependent correction factor, which is applied to each segment of the distribution.22 The resulting size distributions for the three radial positions from the second sampler are shown in Figure 4. The size distribution corresponding to the center-line is estimated from 700 objects measured at 5000 magnification and 5000 objects measured at 60,000 magnification, the half-way (middle) position was estimated from 1030 objects at 5000 magnification and the wall position was estimated from 1200 objects at 5000 magnification. The solid-lines in Figure 4 are from fitting the experimental size data to a log-normal size representation. This approach is used to obtain particle size metrics for the samples obtained at all points in the reactor. A test is made to verify that the size distributions obtained by the combination of thermophoretic sampling and image analysis are comparable to other size measurement techniques such as laser light scattering. At the position of the last sampler (near the exit of the reactor), the size distributions at each of the three radial coordinates are measured by SEM with 10 images each at 5000 magnification. The number size distributions at the center-line, half-way position (middle), and wall position are estimated based on counting 1510, 1380, and 1200 objects, respectively. It is assumed that this sampler is close enough to the exit of the reactor that the particles obtained there will not grow appreciably before being quenched. If this is the case, then these size distributions can be compared with the size distribution obtained by the Malvern laser light scattering from the product sample. This representation of the size distribution is obtained using the product sample and the method described in Appendix A. In making this comparison, it is important to recognize that the Malvern volume size distribution is derived by assuming that all particles are spherical, so when converting the experimental sample number distributions to volume distributions, the volume weighting also uses this assumption, rather than using Eq. 2 for particles greater than 0.7 mm. This comparison is shown in Figure 5. Figure 5. A comparison of the cumulative size distributions by sampling and by laser light scattering. The sampled size distributions are taken from each of the three radial sampling positions at the last sampler (before exit) under the Case A conditions. The laser light scattering distribution is taken from an analysis of the reactor product (integral product) using a Malvern Mastersizer 2000 (see Appendix A for details). In both cases, the size distributions are scaled by assuming that the particles are spherical. The sampled size distributions from the center-line to the wall are calculated by counting 1380, 1510, and 1200 objects, respectively. Published on behalf of the AIChE DOI 10.1002/aic 1433 Modeling—introduction The radially-differentiated size distributions shown in Figures 3 and 4 suggest that the model reactor is not behaving according to the expectations of a plug-flow reactor. As a result, a CFD model of this reactor is built, to derive an understanding of the possible contribution of variable fluid effects on the formation of particles of different sizes within the same axial position of the reactor. For fluid systems, CFD is a well developed approach for predicting the internal states in a complex flow by numerically solving the conservation equations at many positions within the domain of interest.5 The modeling of aerosol synthesis processes, however, is more complicated by the fact that there are two distinct phases: the particle phase, which is solid or liquid, and the host phase, which is a gas. Each of these has a separate set of constitutive equations. To make modeling even more challenging, changes in each of these phases should influence the other phase. Some simplifications are required to be able to model reacting multiphase flows. One approach is to retain the complex description of the reactor flow, and simplify the description of the aerosol dynamics. Johannessen et al.6,7 and Schild et al.8 have studied aerosol material synthesis in complex flows. Johannnessen, for example, 6 studied the formation of Al2O3 by the oxidation of Al-tri-sec-butoxide in an oxidizing diffusion flame. The fluid dynamics of the flame were modeled using Fluent 4.4, a commercial fluid dynamics code. Simplified population and mass balances were written to describe the evolution of the aerosol number and area concentration. These were embedded into each of the cells within the CFD model and solved along with each of the conservation equations for the fluid phase. The aerosol calculations were solved using fluid properties from the converged solution for each cell. In this way, the conditions of the fluid were allowed to exert an influence on the calculation of the aerosol states. The evolving condition of the aerosol however was not fed back into the fluid calculations; the coupling between the aerosol equations and the fluid equations was one-way. The justification for this simplification was that the volume fraction of the aerosol phase is very small, typically 105. Equations that include the effect of a dispersed phase on fluid properties, such as Einstein’s equation for fluid viscosity, confirm that particles do not exert a significant influence on the fluid phase in this range of volume fraction.22 In the case where the volume fraction of aerosol is high enough to affect the fluid characteristics, the situation is further complicated because the fundamental relationships that govern the interaction between the solid and gas phases are not known. Another approach is to simplify the description of the reactor, as is done in the classical discipline of reactor design, where real-world continuous reactors are idealized as one or a combination of plug flow (PF) and continuously stirred tank (CST) reactors. In this case, more complexity can be introduced into the aerosol dynamic equations, such as the inclusion of surface reactions. One example of this approach was utilized by Pratsinis and Spicer9 who studied the simultaneous effects of nucleation, surface reaction and coagulation of particles made from the oxidation of TiCl4 in flames. The effect of the surface reaction was included, but the size distribution was assumed to be monodisperse at all times. 1434 DOI 10.1002/aic Kruis et al.10 also used these simplifications for modeling the effect of coagulation and sintering for particles where the shape was treated as ‘‘fractal-like’’. Instead of assuming a monodisperse population, the shape of the size distribution can be modeled more realistically using a log-normal approach, or a sectional approach. These latter two approaches have been compared by Xiong and Pratsinis11 for the oxidation of TiCl4 in a tube reactor. They reported that the size distributions for both models reduced to a log-normal form when the residence time in the reactor was sufficiently long for coagulation to dominate the process. Kruis et al.10 showed that the results of a monodisperse and log-normal model are about the same for tracking particle dynamics in a process dominated by coagulation. In the present work, we choose a sequential approach, where we begin by modeling the internal flow fields of our reactor using CFD. That information is used to generate a residence time for the trajectory at the center-line of the reactor. Considered in isolation, an element of fluid near the center-line conforms reasonably well to the plug-flow type. From there, we describe the aerosol dynamics by a system of simultaneous ordinary differential equations (ODEs), and compare this to data obtained from in-situ sampling along the center-line. From an understanding of the dominant mechanisms of aerosol growth, changes in the way the reactor is configured can be proposed to optimize the properties of the product. And finally, these decomposer changes are evaluated experimentally. Modeling—computational fluid dynamics The nickel carbonyl tube reactor CFD simulation for this system is developed using CFX-4.5 from ANSYS.12 A full three-dimensional model of the reactor is created in the geometry-modeling module of the ICEM CFD preprocessor software (ICEM CFD is a product of ANSYS). The reactor geometry was modeled as a straight tube, fed from the top. Any contribution or effect from the sampling port purges is neglected. The geometry is imported to HEXA, the hexahedral mesher module, and a hexahedral mesh of about 600,000 cells is created. The mesh is read into the multiblock solver CFX-4.5 from ANSYS. The model solution is calculated with the following assumptions. The reactor is at steady-state in time, in other words, the state of the fluid at each position inside the reactor is time invariant. The flow behavior is transitional between turbulent and laminar. Although the flow-rate is low, and would normally be laminar when fully developed under isothermal conditions, turbulence is introduced from the high momentum of the inlet jet and the instability created by natural convection at the reactor wall. The turbulence is described mathematically with a Menter Low-Reynolds number k-omega turbulence model.12 The fluid is treated as nonisothermal (diffusivity, thermal conductivity, viscosity, heat capacity, and density are all functions of temperature) and buoyancy is included. The fluid is at all times a binary mixture of the carrier carbon monoxide (CO) and nickel carbonyl. The mixture properties are calculated using the ideal mixing law for density, Wilke’s model for mixture viscosity and thermal conductivity,13 and the Chapman-Enskog model for diffusivity of nickel carbonyl in CO.13 A simplified chemical reaction model, which fits Published on behalf of the AIChE June 2007 Vol. 53, No. 6 AIChE Journal Figure 6. The radial temperature profile at the position of the first sampler under the Case A conditions. The solid line is the prediction for a CFD model of a nonreacting flow. The hatched line is a CFD model that includes the enthalpy of reaction. The diamond points are experimental data for an inserted thermocouple junction in a nonreacting flow-rate. The thermocouple has been inserted from the lefthand side of this figure, and the asymmetry in the measured profile is cause by the purge flows used for each of the samplers. the data of Jones14 is included. Each steady-state run took about 4–10 CPU h on 2.8 GHz Intel Pentium 4 PC’s, with 1.0–1.8 GB RAM. Figure 6 shows a comparison of the radial temperature profile for the Case A experiment at the position of the first sampler. The solid and hatched lines are the model predictions with and without the effect of considering the enthalpy of reaction. In the case where the reaction is included, the temperature is a maximum of 60 K less than the prediction that does not include the reaction. In other words, neglecting the chemical reaction accounts for a maximum local temperature deviation of 60 K under the Case A Conditions, where the volume percentage of decomposing nickel carbonyl is about 5%. The diamonds represent thermocouple measurements in a nonreacting flow, discussed in the Experimental section. For these measurements, the thermocouple is fed through the sampling port from the left-hand side of Figure 6. The experimentally derived radial temperature profile is clearly asymmetric, especially near the center of the reactor, being as much as 100 K hotter on the side where the thermocouple is inserted. It is explained in the experimental section that a small purge of inert gas is maintained through each sampling port to keep it clean. It is hypothesized that the asymmetry in the experimental temperature profile is caused by forced convective heat transfer from these purge flows. The resulting velocity profiles at the position of each sampler for the Case A conditions are shown in Figure 7. The resulting velocity profiles at the same position for the Case B conditions are shown in Figure 8 and Case B up-flow in Figure 9. Although the reactor geometry was modeled as a full three-dimensional object, the predicted velocity profiles under Case A, Case B, and Case B up-flow conditions are symmetric through the central axis of the reactor. The velocity profiles in Figure 7 indicate that there is a large flow inversion in the top half of the reactor under Case A conditions. This phenomenon is caused by two factors; the expansion of the inlet jet in the proximity of the tube wall (studied under isothermal conditions by Becker15), and free convection from the heated walls. A flow inversion is also present, although to a much lesser extent in the Case B simuAIChE Journal June 2007 Vol. 53, No. 6 lations (Figures 8 and 9). The significant reduction in internal recirculation is largely because of the increased diameter of the inlet nozzle. Strong spatial velocity gradients inside the reactor are not ideal from the perspective of reactor design. This feature, combined with the parabolic nature of the velocity field, leads to a broad residence time distribution. This is because particles of the size shown in Figure 3 effectively track the convective flow of the host fluid, so that a strong velocity gradient across the reactor diameter will create a vastly different family of residence times for the product. We hypothesize that the strong velocity gradients inside the reactor are contributing to the differences in particle size shown in Figures 3 and 4. Particles near the reactor center-line have experienced a direct path through at high velocity, while particles in the middle have spent more time and particles near the wall are moving more slowly, with the possibility of being returned upstream by the flow inversion, a kind of internal recycle. This undesired feature of the flow should lead to a broadening of the distribution of all properties of the final Figure 7. Axial velocity profiles calculated using the CFD model for the Case A conditions. The initial velocity profile is parabolic with recirculation near the wall. As a consequence, particles made near the center-line will travel through the reactor much more quickly than product near the wall. Published on behalf of the AIChE DOI 10.1002/aic 1435 A simple mass and population balance can be used to describe the formation and growth of aerosol particles by nickel carbonyl decomposition. These can be written as a system of coupled ODEs, shown as Eqs. 3–6 that describe the interaction between nickel carbonyl concentration (CNC, moles/g gas), particle number concentration (N, particles/g gas) and solid volume concentration (f, m3/g gas). In this model, the particle size is assumed to be monodisperse at all points inside the reactor. The decomposition reaction (Eq. 1) is known to proceed by two pathways; a homogeneous reaction, and a heterogeneous reaction. The reaction constants for each of these reactions as a function of temperature have been taken from the reaction mechanisms of Chan16 and Carlton and Oxley,17 respectively. The mass balance for depletion of nickel carbonyl from the gas-phase is shown as Eq. 3. The first term on the RHS of Eq. 3 is the homogeneous contribution of the overall reaction (from Chan) and the second term is the heterogeneous contribution (from Carlton and Oxley). In this mass balance, it is assumed that the surface area of each particle can be approximated by p/4 dp2, which is not strictly true since the particles are not perfect spheres, as has been noted in the experimental section. For all of the kinetic expressions, we have used the same nomenclature as the original authors, which are reproduced in the notation section at the end of this article. Equation 5 is a mass balance for the solid phase. It is the same as the mass balance for the gas-phase, using the molar volume of the Figure 8. Axial velocity profiles calculated using the CFD model for the Case B conditions in down-flow mode. product that depend directly on residence time, such as particle size. Modeling—a simplified plug-flow model To test the hypothesis that a highly variable residence time inside the reactor can broaden product properties such as the size distribution, it is important to show that the fundamental mechanisms of particle growth predict that particles will continue to grow as they spend more time in the reactor. To do this, a simplified reactor model is proposed, following the methodology of Pratsinis and Spicer9 and compared against the experimental data obtained by thermophoresis for the Case A conditions. For comparing against the experimental results of Case A, it should be noted that the predicted internal velocity profiles shown in Figure 7 do not closely conform to either the CST or PF reactor type for the purposes of developing a simple reactor model. Nevertheless, there is one part of the reactor that does behave according to the plug-flow type, and that is the flow at the center-line of the reactor. If we assume that there is no radial mixing of gases at the center-line, then the central velocities shown in Figure 7 can be integrated to give estimates of the residence time at the center-line of the reactor. This allows the sampler positions in the reactor to be paired up with their corresponding residence times at the center-line. 1436 DOI 10.1002/aic Figure 9. Axial velocity profiles calculated using the CFD model for the Case B conditions in upflow mode. Published on behalf of the AIChE June 2007 Vol. 53, No. 6 AIChE Journal solid phase as the proportionality constant to convert moles to volume. The population balance, Eq. 4, is used to describe the evolution of N, the number of particles per gram of the gas. Following Johannesen,6 an algebraic transformation from particles per volume to particles per gas weight has been made to account for the nonisothermal conditions of the reactor. The first term on the RHS accounts for the creation of new particles by the homogeneous reaction, and the second term accounts for the consumption of particles by collisions. The coagulation constant b is calculated from the Fuchs interpolation formula,18 which includes the effect of temperature, particle size, and fluid properties of the aerosol throughout the entire range of particle size. The comparison results from Figure 5 suggest that particle fragmentation can be neglected from this number balance. Equation 6 transforms the volume and number concentrations to particle For spherical particles, it should be ffiffi qffiffiffisize. 6f simply dp ¼ p N . However, as discussed in the experimental section, the majority of the particles after the first sampler on a volume basis are more accurately described by the power-law relation shown as Eq. 2. For this reason, Eq. 6 is derived from the experimentally derived power-law relation shown in Eq. 2. 3 dCNC k2 ðTÞCNC ¼   pdp2 Nr2g rw0 C2NC dt 1 þ b2 ðTÞCCO RT (3) dN k2 ðTÞNava CNC =4 1 ¼  bðT; m; dp ÞN 2 rg dt 1 þ ðb2 ðTÞ=4ÞCCO RT 2 (4) df dCNC v1 ¼  dt dt rffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi f 2:5 dp ¼ 1:4  104 N June 2007 Vol. 53, No. 6 This was calculated using Eqs. 3–6 for the Case A conditions. Experimental data (one test plus a replicate) are shown from the Case A conditions. Reasonable agreement between the model and the experimental data was obtained for reactor conditions where the particle size distribution was roughly log-normal (after the second sampler). Results (5) (6) This system of ODEs has been solved using the ODE15S solver of MATLAB version 6.1.19 The center-line temperature estimates required for this model are taken from the internal center-line temperatures shown in Table 1 and the starting conditions at the reactor inlet are taken from the inlet flowrate and composition details described in Tables 1 and 2. The trajectory of particle number concentration N and particle diameter dp are shown in Figure 10, along with two sets of experimental data obtained from the center-line under the Case A conditions. Following the rationale of Pratsinis et al.20 the monodisperse plug-flow model is compared with experimental data by considering the experimental diameter of number average volume. The model and experimental data together confirm that particles grow in size monotonically with time inside the reactor after initial nucleation, at first quickly by a surface reaction, and then more slowly by Brownian coagulation. It should be mentioned that this plugflow analysis assumes implicitly that there is negligible radial conduction of heat, particles or gas species from some thin core of fluid at the center-line. For this reason, repeating a plug-flow model of this type for other stream-lines within the reactor would require additional assumptions about the degree of mixing and dispersion within the reactor. AIChE Journal Figure 10. The model trajectory for number concentration and particle size as a function of residence time for the first principles model assuming plug-flow fluid characteristics. In-situ sampling data obtained from experiments run under the Case A conditions reveal that this reactor does not behave ideally from the point of view of making products with narrow ranges of properties. In particular, under the Case A conditions, we have shown that the size distribution of particles being produced is strongly influenced by the position inside the reactor (Figures 3 and 4). Thermophoretic sampling at the end of the reactor and comparison of the particle size distribution of the product measured by light scattering has shown that the in situ sampling method is capable of obtaining good representations of the particles within the reactor. This suggests that the radial variation in size distributions is a real feature of this type of reactor. Comparison of the center-line sampled data with a simplified plug-flow model derived from first principles shows reasonable agreement between the model and the experimental data, and this is verified by running a replicate experiment under the Case A conditions -this is shown in Figure 10. The trajectory and evolution of the particles through the reactor has been described in greater detail elsewhere.3 To summarize here, there is an initial nucleation event that creates a high number concentration of small particles. After nucleation, the particles continue to grow by a surface condensation reaction, and finally by coagulation and sintering. Because the first principles model has a number of limiting assumptions; it cannot describe the particle distribution and it is only valid for a plug-flow trajectory, it agrees reasonably well with the experimental data, suggesting the fundamental particle growth mechanisms that have contributed to the observed experimental results. After the nickel carbonyl gas has been depleted from the reactor, the particles will continue to grow by Brownian coagulation. If we ignore Published on behalf of the AIChE DOI 10.1002/aic 1437 Table 3. Values of the Plug-Flow Comparison Index (Eq. 7, Normalized) at Each of the Axial Positions Shown in Figures 8 and 9 for Case B Down-Flow and Up-Flow Axial Distance from Inlet, cm 5 10 15 20 25 Case B Down-Flow 1.83 7.79 8.31 9.33 9.93      101 102 102 102 102 Case B Up-Flow 1.50 2.75 5.55 7.37 8.46      101 102 102 102 102 dynamic particle fragmentation in the reactor, this means that particles that spend more time in the reactor have the opportunity to grow to larger sizes. The complementary computational fluid dynamic model of the reactor under Case A conditions reveals that the current design of the reactor creates the conditions for a broad residence time distribution. With an inlet nozzle/ inside diameter <0.1, and buoyancy induced from the hot walls, a strong upward flow of gas is created, which opposes the forced downward convection of the reactor flow. By comparing the Case A simulation with the Case B simulation where the inlet nozzle/ inside diameter ¼ 0.5, it is clear that the intensity of the recirculation has been reduced (Figures 7 and 8). According to Becker, under isothermal conditions, recirculation should subside completely as the ratio of inlet nozzle/ inside diameter approaches 0.5.15 Unfortunately, redesign of the nozzle is not enough to completely rehabilitate the internal flow to make it more like a plug-flow. This is because of the natural convection that is induced by the hot walls. To remedy this situation, it is proposed that the reactor can be turned upside down to make the velocity profiles conform more closely to the plug-flow type. In this case, the inlet profile will still be parabolic (as shown in Figure 9), however the buoyancy near the hot wall will actually increase the velocity where it is low, which will tend to flatten the overall velocity profile. To check this hypothesis, the CFD simulation for Case B can be compared with simulation for Case B up-flow, where identical conditions are maintained except for orientation of the reactor and feed inlet position. In Figure 8, which shows the axial velocity profiles for Case B (down-flow) the inlet velocity profile is parabolic, it is low near the walls, and reaches a maximum near the center. The negative velocity near the wall that is the result of buoyancy from the hot walls exacerbates this effect; the velocity near the walls becomes less, resulting in a greater difference between the velocity at the wall and velocity near the center. In Figure 9, which shows the equivalent axial velocity profiles for Case B up-flow, the initial velocity profile is also parabolic. In this case, the buoyancy-induced convection from the hot walls actually increases the velocity near the walls, where it is low. As a result, the superposition of the buoyancy driven flow and the natural parabolic flow of the inlet create an artificial PF, where the velocity profile is approximately flat. To compare the closeness of the reactor flow to the ideal ‘‘plug-flow’’, a comparison index can be invoked to quantify the difference between two flow profiles based on the mini1438 DOI 10.1002/aic mization of variation of the axial velocity distribution. The quantity to be minimized is the summation over the radius of the deviations between the local velocity and the mean velocity -the minimum of this quantity corresponds to the condition where the velocity profile is flat, and all of the fluid elements in the flow field have the same residence time. Each of the contributions to this summation should be weighted by the corresponding mass flux. From the principle of continuity, the mass flux is proportional to the axial velocity multiplied by the square of the radius. The comparison index, which should be minimized, is therefore i¼i max X 2 jni jðri2  ri1 Þjni  navg j (7) i¼1 where ni and ri are the axial velocity and tube radius for the ith element of the summation. If the velocity profile is symmetric about the center of the tube, then the summation can be over one half of the tube diameter. The native version of this comparison index will have units of distance to the fourth power divided by time squared. To make this index dimensionless, the summation shown in Eq. 7 can be divided by the square of the average velocity and the square of the full radius, which is the convention, used in the accompanying comparison table. For two velocity profiles with all other conditions being equal, the plug-flow characteristics will be best for the profile with the smaller value of this comparison index. Table 3 shows this comparison index for the Case B ‘‘upflow’’ and ‘‘down-flow’’ velocity profiles shown as Figures 8 and 9. As expected, the ‘‘up-flow’’ configuration has a profile that is mathematically closer to the PF ideal at each of the five sampled positions within the reactor. To verify the impact of this configuration change, three experiments were run at the Case B conditions; the first with the standard down-flow configuration and the remaining two with the upflow configuration. The size distribution of the accumulated product is measured by laser light scattering using a Malvern Mastersizer 2000 as described in Appendix A. The volume size distributions are shown as Figure 11. These size distribu- Figure 11. Size distributions of powder, measured by laser light scattering (Malvern Mastersizer 2000) made under the Case B condition (one experiment in down-flow and one experiment in up-flow with one replicate). Published on behalf of the AIChE The main effect of the up-flow configuration is to eliminate the coarse shoulder of the size distribution. June 2007 Vol. 53, No. 6 AIChE Journal tions provide some evidence that the more uniform velocity profile of the reactor in the up-flow configuration has created a product with a narrower size distribution. In particular, the coarse tail of the size distribution greater than 4 mm has been completely eliminated using the up-flow configuration. A similar comparison of size distributions by thermophoretic sampling is not available. Professor Sotiris Pratsinis of the ETH provided many fruitful discussions around aerosol modelling and we are indebted to him for many aspects of the first principles model. Randy Shaubel of Inco Technical Services Limited ran the experimental lab reactor. Bill Nowosiadly of Inco Technical Services Limited conducted many of the SEM and TEM sessions that resulted in Figure 4. We would also like to express our gratitude to the staff at Inco Technical Service Limited. Notation Conclusions In-situ sampling by thermophoresis is an excellent way to obtain spatially resolved samples of aerosol particles, and when combined with high magnification microscopy and image analysis, this technique can be used to chart the genesis and growth of particles inside an aerosol reactor. A higher level of process knowledge can be obtained than the usual procedure of treating the reactor as a black-box and trying to understand internal dynamics through an analysis of the integral product properties. At the same time, thermophoretic sampling is laborious, costly and requires special equipment to be successful. In this study, thermophoretic sampling is used to obtain data that can be used to check the accuracy of a simplified first principles model, which includes mass and number balances developed from existing equations for nickel carbonyl chemical kinetics and aerosol physics. The experimental and model results are in reasonable agreement, which validates the first principles-based equations used in the simplified model. In particular, it indicates that particles grow monotonically with time in the reactor, and that particle fragmentation does not seem to have a major impact under these conditions. A more sophisticated model would be required to accurately describe the expression of these equations under the complex flow profiles of a real reactor. A CFD analysis that only considers the system as a single fluid phase confirms the extreme velocity gradients that exist inside a conventional nickel carbonyl tube reactor. In particular, there is a massive zone of recirculation in the reactor. Radial samples obtained from an axial position approximately coincident to this zone confirm that the particles have grown in size qualitatively in proportion to their residence time. This provides a link between simplified fluid models and the particle size distribution. It is also shown that the velocity gradients can be greatly improved by increasing the nozzle diameter of the reactor. Finally, it is shown that the velocity gradients in the radial direction inside a hot-wall tube reactor can be reduced by operating the reactor in an ‘‘up-flow’’ configuration. In this case, the inlet velocity profile is approximately parabolic and the additional upward velocity near the walls contributes to make the overall profile closer to the plug-flow type. A plugflow comparison index is developed to quantitatively show that this is the case. Several experiments are run under the Case B conditions in both the conventional down-flow mode and the proposed up-flow mode. An analysis of the product size distribution shows that a coarse shoulder is removed from the size distribution, suggesting that the proposed improvement in the reactor configuration has been successful in narrowing the product size distribution. This result shows how CFD modelling can be used to prototype new aerosol reactor configurations, to optimize product properties. AIChE Journal Acknowledgments June 2007 Vol. 53, No. 6 d ¼ mean diameter (m) dp ¼ diameter for idealized monodisperse population (m) CCO ¼ carbon monoxide concentration (mol/m3) CNC ¼ nickel carbonyl concentration (mol/g gas) k2, b2 ¼ kinetic constants from Chan,16 also defined in Ref. 2 N ¼ particle number concentration (particles/g gas) Nava ¼ Avagadro’s number 6.02  1023 ri ¼ reactor radius at position i (m) 1 rw ¼ surface rate constant from Carlton and Oxley,17 also defined in Ref. 2 R ¼ ideal gas law constant b ¼ Coagulation constant (m2/sec) f ¼ Particle volume concentration (m3/g gas) ni ¼ radial velocity at position i (m/sec) navg ¼ average radial velocity across radius (m/sec) rg ¼ gas density (kg/m3) Literature Cited 1. Mittasch A. German Patent 500692, 1925. 2. Wasmund EB. Powder making in an aerosol tube reactor. PhD Thesis. Hamilton: McMaster University, 2005. 3. Wasmund EB, Coley KS. In-situ sampling uncovers the dynamics of particle genesis and growth in an aerosol tube reactor. J Mater Sci. 2006;41:7103–7110. 4. Dobbins RA, Megardis CM. Morphology of flame generated soot as determined by thermophoretic sampling. Langmuir. 1987;3:254–259. 5. Ferziger JH, Peric M. Computational Methods for Fluid Dynamics (2nd edition). Berlin: Springer, 1999. 6. Johannessen T, Pratsinis SE, Livbjerg H. Computational fluid-particle dynamics for the flame synthesis of alumina particles. Chem Eng Sci. 2000;55:177–191. 7. Johannessen T, Pratsinis SE, Livbjerg H. Computational analysis of coagulation and coalescence in the flame synthesis of titania particles. Powder Technol. 2001;118:242–250. 8. Schild A, Gutsch A, Muhlenweg H, Pratsinis SE. Simulation of nanoparticle production in premixed aerosol flow reactors by interfacing fluid mechanics and particle dynamics. J Nanoparticle Res. 1999;1:305–315. 9. Pratsinis SE, Spicer PT. Competition between gas phase and surface oxidation of TiCl4 during synthesis of TiO2 particles. Chem Eng Sci. 1998;53:1861–1868. 10. Kruis FE, Kusters KA, Pratsinis SE. A simple model for the evolution of the characteristics of aggregrate particles undergoing coagulation and sintering. Aerosol Sci Technol. 1991;19:514–526. 11. Xiong Y, Pratsinis SE. Gas phase production of particles in reactive turbulent flows. J Aerosol Sci. 1991;22:637–655. 12. CFX 4.5 User’s Manual, 2002. 13. Reid R, Prausnitz, J, Sherwood T. The Properties of Gases and Liquids (3rd edition). New York: McGraw-Hill, 1977. 14. Jones L. Formation of nickel powder by decomposition: flow measurements within a laboratory-scale nickel carbonyl powder decomposer. PhD Thesis. Swansea: University College of Wales, 1994. 15. Becker HA. Concentration fluctuations in ducted jet mixing. PhD Thesis. Boston: Massachusetts Institute of Technology, 1961. 16. Chan R. The heterogeneous and homogeneous decomposition of nickel carbonyl. PhD Thesis Toronto: University of Toronto, 1961. 17. Carlton HE, Oxley JH. Kinetics of the heterogeneous decomposition of nickel tetracarbonyl. Am Inst Chem Eng J. 1967;13:86–91. 18. Seinfeld JH. Atmospheric Chemistry and Physic of Air Pollution. New York: Wiley, 1998. Published on behalf of the AIChE DOI 10.1002/aic 1439 19. Matlab Users Manual, Version 6.1. Natick Massachusetts: The Mathworks, 2000. 20. Pratsinis SE, Arabi-Katbi OF, Megaridis CM, Morrison PW, Tsantilis S, Kammler HK. Flame synthesis of spherical nanoparticles. Mater Sci Forum. 2000;8:511–518. 21. Yazicioglu AG, Megaridis CM. Measurement of fractal properties of soot agglomerates in laminar coflow diffusion flames using thermophoretic sampling in conjunction with transmission electron microscopy and image processing. Combust Sci Technol. 2001;171:71–87. 22. Friedlander SK. Smoke, Dust and Haze, 2nd edition. New York: Wiley, 2000. Appendix A Method for measuring particle size using the Malvern Mastersizer 2000 Samples are taken from the filter bag after the reactor described in the experimental section. Approximately 0.5 g of nickel powder is dispersed in water using a small amount of surfactant and external ultrasonic agitation. The sample is added into a flow-through water bath. The particle size distribution is calculated by the instrument. Appendix B Estimating the scaling relationship between mean particle diameter and volume using image analysis From 13 of the micro-grids described in the experimental section, 10 SEM images are taken from each grid to charac- 1440 DOI 10.1002/aic terize the particle population. These populations each contained about 500–1000 objects. The image analysis software described in the experimental section is used to sort the objects in each population based on a roundness metric. The objects with roundness metric greater than 2.5 are deemed to be nonround, in other words, their volume could not accurately be described by the formula for a sphere. Nonround particle represent from 5 to 15% of the population on a number basis, in general, their frequency increases with distance from the reactor inlet. To estimate, the volume of each of these nonround particles, the image analysis software is used to calculate a skeleton length using the thinning process filter. Then, the volume of the particle is calculated by assuming that the particle consisted of a flat, branched cylinder. The details can be found in Section 5.2 of Ref. 2. Although this calculation neglects particle volume projecting out of the viewing plane, it is a much closer approximation of the particle volume than the assumption of sphericity. Others have estimated the ‘‘2D–3D correction factor’’ is about 1.24.21 For each population, the calculated volume is regressed against the mean diameter. The average pre-exponential and exponent for each of these regressions is as shown in Eq. 2. This is the basis for weighting nonspherical particles when converting from number to volume size distributions throughout this work. Manuscript received Aug. 23, 2006, and revision received Mar. 5, 2007. Published on behalf of the AIChE June 2007 Vol. 53, No. 6 AIChE Journal