Computer Simulation of Flow Processes in Fluidized Bed Reactors

The paper outlines and exemplifies a multi-dimensional multi-phase Computational Fluid Dynamics (CFD) model for the various processes that occur in fluidized bed reactors. The model is based on the Eulerian description of the two phases: gas and particles. This means that separate conservation equations are set up for these phases. Calculations are shown for some examples that include flow patterns in circulating beds.


The problem
Many different processes influence the total performance characteristics of fluidized bed reactors. Among these are the gas/particle flow and distribution, chemical reactions, and heat transfer. All these processes occur simultane· ously and are intimately coupled to each other. The performance of the reactor is strongly dependent on the interactions between these processes and the geometrical design of the reactor. Because of the complexities of such problems, there are no analytical methods available that can describe in detail what is occurring inside a fluidized bed reactor. In each new reactor it is very often the geometrical details that will determine the performance of the final design of the reactor. Until recently the design of fluidized bed reactors has therefore been based on experimental work done on a laboratory and pilot-plant scale, and on extrap· alation of these data to full scale by the application of various empirical scaling techniques.
These design techniques may be very expensive, and give insufficient and inaccurate information about the performance. New computer techniques for fluid flow processes have now made it possible to develop new design methods based on mathematical modelling of the basic processes that occur in the actual geometrical model of the reactor. To apply these methods, TMIH

96
the basic conservation equations of flow, chemical reactions, and heat transfer as well as the coupling between them and the coupling between the two phases, must be solved for a multi-dimensional calculation domain.

Previous work
Gidaspowll, has reviewed the CFD models that had evolved up until about 1985. The models at that time had demonstrated the capability to simulate some of the crucial characteristics such as bubble formation and development in bubbling fluidized beds. The models assumed a two-phase gas-particle flow with given apparent viscosities in the gas and in the particle phases. Shih et aUl, introduced a model that takes account of the multi-sized particle problem and showed application to a sedimentation problem. The work of Tsuo and Gidaspow 3 l, showed that the model was able to calculate the flow patterns in circulating fluidized beds. This model had to assume an apparent viscosity for the particle phase to be able to predict the experimentally observed phenomena. The works of Ding and Gidaspow 4 l, and Ma and Ahmadi 5 l, propose to couple the apparent viscosities in the two phases by introducing turbulence models for each of the phases. A comprehensive model for two-phase flow, heat transfer and chemical reactions in a fluidized bed gasifier has been presented by Gidaspow et al.6l_

Objectives of paper
The present paper will present a two-phase multi-dimensional model for the processes that occur in fluidized beds for catalytic chemical reactions. The model will be verified against experimental data for isothermal flow fields in various systems.

Fluid dynamics model
The governing equations are deduced based on the Eulerian concept. This means that the phases share the space and interact with each other. All equations are such that all volume fractions may take values between zero and one.
Based on the works of Ding and Gidaspow 4 l, and Ma and Ahmadi 5 l, the following flow model for gas particle flow may be formulated using the following dependent variables: ag, ap, Uj,g, uj,p, kg, kp,cg, Cp and p. Here the variables have the index g for the gas phase and index p for the particle phase. The variables are volume fractions, velocity components, turbulent kinetic energy, dissipation of turbulent kinetic energy and the pressure. Conservation equations written in Cartesian tensor notation may be formulated as:

Mass balances:
Gas: ( 1 ) Particles: ( 2) These equations are valid for catalytic reactions with no mass transfer between the phases.
Where <l>s is the form factor of the particles.
Shear stresses are related to the gradients of velocity components and a turbulent viscosity as: Particles: . Where ap,m is the maximum solid volume fraction for a randomly packed bed. ap,m is taken to be 0.65. Solid stress is given by: Turbulent kinetic energy: Gas: ( 11) Particles: (12) Here, the effective transport coefficients are related to the turbulent viscosities as: Particles: The dissipation in the particle phase is given by an algebraic equation as: (15) where r is the restitution.
Here, the effective transport coefficient is related to the turbulent viscosity as: The heat transfer is taken account of by solving the two energy equations for the two phases where heat is generated by the exothermic (homogeneous) gas phase chemical reactions and then transported between the phases.

98
These equations are: Gas: Particles: (19) Here, the effective transport coefficients are related to the turbulent viscosities as: ( 20) where Xtam,g is the laminar conductivity of the gas and Cp,g is the specific heat of the gas.
The volumetric heat transfer coefficient hvpg is calculated as: The basis for calculating the volumetric heat transfer coefficient is based on correlations of the Nusselt number NP in the various flow regimes characterized by the Reynolds number and volume fraction of gas according to:

Chemical reaction model
Since chemical reactions are only taking place in the gas phase, only species conservation for this phase is needed. These are: Gas: Here, ri is the net bulk gas-phase reaction rate for component j, and ri.c is the gas phase reaction occurring on the surface of the catalyst particle for component j. The effective gasphase transport coefficient is calculated as: DYj,g is the molecular diffusivity of the gasphase component with mass fraction Yj,g·

Solution Procedure
Solution of the set of partial differential equations given above is performed by finitedomain methods. The calculation method is taken from the work of Spalding 7 l. The calculation domain is divided into a finite number of main grid points where pressure, densities, void fractions, turbulent quantities, enthalpies and mass fractions of the chemical species are stored. The velocity components are, on the other hand, stored at grid points located midway between the main grid nodes. The relevant conservation equations are integrated over control volumes surrounding the relevant grid points in space and over a time interval. This integration is performed using upwind differencing in space and implicit differencing in time. The result of this is a set of non-linear algebraic equations which are solved by the application of the well known Tri-Diagonal Matrix Algorithm (TDMA) used along the various coordinate directions. For the void fraction and velocity equations, a point iteration method is used. Special care has been taken to solve the coupling of the momentum equations with the continuity equations. For this we use the so-KONA No. 10 (1992) called IPSA procedure proposed by Spalding 7 l. The Partial Elimination Algorithm (PEA)?l, is used to de-couple the drag force between the phases in the momentum equations.

Simulations
Parts of the model given above are compared against two isothermal non-reacting gas-particle flow situations. In the present paper, focus has been directed towards calculating the flow fields and particle distributions and comparing the predictions with data from literature and data generated in our own laboratory using Laser-Doppler Anemometry (LDA). The first case is related to a core annulus type flow in an 11-m-high, 0. 305-m diameter tube of Bader et al 8 l, and the second is related to our own laboratory setup of a 1-m-high, 0. 03-m diameter circulating fluidized bed.

CASE 1: Large-scale tube of Bader et al 8 l
The particles used in this study had a diameter of 76 ,urn, a density of 1714 kg/m 3 and a sphericity of 1.0. The tube was initially empty and the flow was started at the inlet with flow velocities of 0. 228 m/ s for the particles and 3. 7 m/s for the gas. The particle volume fraction at the inlet was fixed at 0. 25. The calculation was performed using a Cartesian description with 76 grid points along the height and 30 grid points across the diameter. It took about 14 seconds to fill the the empty tube. Fig. 1 gives an overview of the results of stream lines for gas, particles and distribution of volume fraction of particles at 16 seconds after start-up. Please note that the contours are expanded in the transverse direction to improve visualisation of the details. It was generally found that the flow did not reach a true steady state situation, but had a cyclic behaviour. The general trend was that the gas and particle flow had an upward flow in the centre of the tube and a downward flow along the wall. Bader et al 8 l, had measured profiles across the tube diameter of both volume fraction of gas (porosity) and particle velocity. Fig. 2 and Fig. 3 show comparisons between measured and predicted profiles at a height of 9. 1 m from the inlet. The predicted profiles were determined by averaging over 5 seconds. Fig. 2 shows that the porosity is well predicted, whereas the particle velocity profile in Fig. 3 shows some discrepancies between Measured and predicted profiles of particle velocity as a function of distance across the tube diameter at 9 .l m from the inlet. experiments and predictions. The general trend seems, however to be well predicted, namely upward flow in the centre and downward flow along the wall.

CASE 2:
Laboratory-scale circulating fluidized bed The vertical particle flow velocity in the lab.
-scale fluidized bed reactor was determined by a laser Doppler anemometer delivered by Dantec. This system made it possible to measure both positive and negative particle velocities. Velocity profiles across the tube diameter were measured at three different heights above the ·distributor plate. The average particle diameter was about 55 J.lm and the density was 1600 kg/m 3 • The particles that were transported from the reactor were separated in a cyclone and directed back into the fluidized bed just above the distributor plate. The gas superficial velocity was 1. 7 m/s. The calculation was performed using a Cartesian description. The calculation domain included both the vertical tube and the separator and return tube. 102 grid points were used in the axial direction and 26 points were used in the radial direction. 75 50 25 II 75 50 Fig. 4 shows the predicted contours of volume fraction of particles, particle and gas streamlines. The figure shows that the gas and particles are flowing upwards in the centre of the tube and have regions with negative velocities. The figure also shows a significant slip between the gas and particle flow. This is especially seen in the separation section of the flow system. Fig. 5 shows a comparison between the measured and predicted profiles of the vertical particle velocity at three different heights above the distributor plate. It is seen that the overall agreement is good, i.e. the flow is upward in the centre of the tube and is directed downwards along the wall. However, the negative velocity is overpredicted.

Discussion
The comparisons between measurements and predictions in the two cases given above are reasonable. However, more work is needed to improve the model. This is particularly so for the turbulence model that determines the effective viscosity for the two phases. It must also be mentioned that the simulations were performed with a Cartesian geometry description, whereas the geometry is axisymmetric. Calculations should therefore also be performed using an axisymmetric geometry description.

Concluding Remarks
A comprehensive multi-dimensional CFD model for gas/particle flow, heat transfer and chemical reactions in fluidized bed reactors has been proposed. The model is incorporated in a 2D computer code and some initial verification simulations of isothermal flows have been carried out. The predictions give reasonable agreement with data from literature and our own experimental data. However, more model development is still needed especially with regard to turbulence modeling. In addition to this, the model needs to be tested against experimental data where chemical reaction and heat transfer is included.

Acknowledgements
The present work is financially supported by N orsk Hydro a.s. and the Royal Norwegian Council for Scientific and Industrial Research (NTNF). The authors are grateful to the following MSc students that have contributed to the work: E. Lode, H. Aadland and ]. Stokke.

c
: constant Cp,g : specific heat for gas Cpg : gas -particle friction coefficient Cct : drag coeficient C,,C,C,": constants in turbulence model dp : Greek : radial distribution function : solid stress modulus : gas and particle enthalpy : volumetric heat transfer coefficient : heat of reaction for reaction no. k : kinetic energy of turbulence in gas and particle phase : Nusselt number : pressure : gas phase reaction rate in the bulk phase and on the surface of the partides : restitution factor : Reynolds number : specific surface area of particles :time : gas and particle temperature : Lagrangian time scale of turbulence : j-component of velocity for gas and particles : velocity vector for gas and particles : coordinate direction in i-direction : mass fraction of component j in the gas phase ag,ap gas and particle volume fraction r ¢g,r <!>,p effective turbulent transport coefficient for scalar variable <I> for the gas and particle phase au Kroenecker delta cg,cp dissipation of kinetic energy of turbulence x molecular thermal conductivity J1 1 ,g,J1t,p turbulent viscosity of gas and particles pg,{Jp density of gas and particles Tu,g, -ru.P stress tensor in gas and particle phases <I>s form factor of particles Subscripts g :gas phase lam : laminar p : particle phase t : turbulent