The Discrete Phase Modelling Governing the Dynamics of Biomass Particles Inside a Fast Pyrolysis Reactor

: The influences of several important biomass pyrolysis process parameters such as the biomass feedstock type, flow pressure, biomass feeding rate and biomass particle size play an important role to ensure an efficient pyrolysis process. Unfortunately, the determination of these parameters can be cumbersome and often requires the method of trial and error. As a result, our work discusses the idea of the application of discrete phase modelling (DPM) in the fast pyrolysis process so that the optimum value of these essential parameters can be determined numerically. The numerical test demonstrated in this paper involving the fast pyrolysis of wood indicates that the application of DPM in the simulation is feasible for obtaining the initial prediction of the optimum process parameters.


INTRODUCTION
At the fundamental level, pyrolysis is the thermal degradation process of biomasses in the absence of oxygen. [1][2][3][4][5][6] Pyrolysis process produces two distinct classes of products, which are the primary products and secondary products. Conventionally, the primary pyrolysis products are classified into condensable and non-condensable products. For instance, tar is a condensable product whereas volatiles and char are non-condensable products. Also, the condensable products are often classified as liquids, whereas the non-condensable mainly consists of CO, CO 2 , H 2 and carbons. 1,2,5 Pyrolysis of biomass is one of the utmost effective technologies used to supply biofuel. Furthermore, it is a simple, cheap and promising technology which include bio-oils and bio-chars as the two substantial pyrolysis products. These products have their influential values which can be utilised in a viable way. An extensive experimental study concluded that the transformation of waste biomasses into beneficial and valuable energy products by pyrolysis is an appealing option. [20][21] Unfortunately, the process conditions, as well as the raw materials used, have varying effects on the yield of each product. Often, these process conditions are determined through experience via experimental observation or trial and error. Thus, a computational model of the pyrolysis process is required to enable the prediction of such optimum process conditions which will maximise the process yield.
Above all, this work focuses on the numerical modelling of biomass particles dynamics in a typical fast pyrolysis reactor. In general, fast pyrolysis involves the process of rapid heating of organic materials in the absence of oxygen. Specifically, such a process occurs inside an entrained flow reactor (EFR). 1 In this work, a computational fluid dynamics (CFD) computer program was used to account the tedious mathematical manoeuvring in the numerical analysis of this model. The EFR model was first designed using CAD software, and then the mesh of the model was generated using ANSYS meshing program (see Figure 1). The generated mesh of the EFR model was used by ANSYS Fluent for the numerical calculation of the solution to the problem. In this work, we proved that the application of discrete phase modelling for the fast pyrolysis process is capable of predicting the optimal process parameters such as the particle size and shape and the biomass particle residence time. Here, we define the particle residence time as the amount of time where a biomass particle flows from the EFR inlet to the outlet.
Most importantly, the residence time plays a vital role in ensuring efficient pyrolysis process. Correspondingly, the optimum residence time is always preferred such that the entire biomass particle volume is fully pyrolysed without leaving the unreacted fraction of the biomass particle flowing through the EFR outlet. Inevitably, a few process parameters may also affect the composition and temperature of gases products, such as the oxygen/steam-to-fuel ratio, particle size, biomass feeding rate and injection velocity. The particle residence time is sensitive to the changes in any of these process parameters. In particular, particle size influences the residence time required to complete the pyrolysis process. 7,8 Here, a smaller biomass particle size enables the pyrolysis process to produce a higher gaseous product yield. Such a process behaviour is due to longer residence time which will allow more conversion of biomass. 9 Papadikis et al. demonstrated that smaller particle sizes could produce a higher amount of syngas. 4 Similar findings were observed in a fluidised bed where a process is governed by the enhanced reaction kinetics for smaller particles. 10 Wood is mainly made up of cellulose, hemicellulose and lignin. 1 Also, wood is a complex fuel that undergoes dramatic changes as it thermally degrades. When wood is heated, it undergoes pyrolysis process which liberates organic gases and leaves behind carbon-rich charcoal. Besides, it is an exothermic reaction that tends to be self-sustaining once started. 1,11 Indeed, various models have been developed to describe the fast pyrolysis process, where most of these models attempt to predict the heat transfer, mass and momentum transport effects for every single biomass particle undertaking the pyrolysis process inside an EFR. 3 The EFR used in this work is made up of stainless steel which has two inlets with each connected to the biomass particle feeder and the nitrogen gas input respectively. Also, the nitrogen gas is used to transport the pyrolysing biomass particle from the feeder input to the EFR outlet. It is vital to keep in mind that nitrogen gas is an inert gas, thus, it will not react with the biomass particles.

MODEL DESCRIPTION AND ASSUMPTIONS
At this instance, the EFR used by Energy Laboratory, Universiti Sains Malaysia (USM), is illustrated in Figure 2. Fundamentally, the proposed fast pyrolysis model is a type of mass transport model, where the biomass particle is transported by the inert nitrogen gas from the biomass particle injection point up to the outflow of the EFR. 7,12 Recall that when the biomass particles are transported within the EFR, they undergo thermal degradation and produce pyrolysis products such as char, tar and gases. These products are also transported by the inert nitrogen gas throughout the EFR.

Biomass feed
Biomass particle injection Sand (as heat transfer medium) Pyrolysis products outflow (to the cyclone) Nitrogen gas input Also, sand is filled at the bottom of the reactor and the nitrogen gas flows into the reactor from the bottom end of the reactor. The bubbling of the mixture of sand, biomass particle and the nitrogen gas do occur, but the sand's purpose is merely to provide even heat transfer from the EFR environment to the biomass particle. Hence, the process efficiency can be maximised. However, the bubbling mixture of sand and nitrogen gas is neglected. Plus, it is assumed that a perfect and even heat convection occurs within the region of the reactor. Furthermore, the specific heat capacity of the biomass particle is a function of the radius of the particle. Such behaviour is due to the assumption that the reaction changes the particle density, but not the size. 3,4 The size of the particle is assumed to remain constant throughout the entire course of the pyrolysis reaction. 4 In other words, the simulated particle is assumed to have no volume shrinkage during the pyrolysis reaction. In this numerical model, the pyrolysis products are assumed to disintegrate from the portion of the particles shortly after it is produced. 3 However, such a product disintegration does not cause the particle to change its size. Besides, the density and mass fraction of the biomass particle is changed instead.
Additionally, the reaction will occur at the surface of the particle. Then, it continues to react towards the centre of the biomass particle. 3,4 The injected biomass particles do not necessarily have equal sizes. Therefore, it is appropriate to assume that the biomass particle size population is distributed according to the Rosin-Ramler particle size distribution. The Rosin-Ramler probability distribution function will be described later in this article. Also, the biomass feeding injection is supported by the gravitational force. Here, biomass particles are injected along the negative direction of the y-axis (the direction of the gravitational force). At the injection point, the terminal velocity of the biomass particle is v 0 . It is critical to consider the drag force and buoyant force induced by the nitrogen gas continuum that acts upon the biomass particle. Intuitively, the relative velocity of the biomass particle is given by: Here, v(t) is the instantaneous terminal velocity of the biomass particle and U is the phase velocity of the nitrogen gas. 3,4 Principally, the flow velocity of the nitrogen gas will always higher in magnitude compared to the terminal velocity of the biomass particle. The purpose of having the flowing nitrogen gas is to transport the biomass particle and the pyrolysis products to the EFR outlet. The drag force exerted by the nitrogen gas on the biomass particle is: where t is the density of the nitrogen gas, v is the velocity of the biomass particle, C D is the drag coefficient of the nitrogen gas and A is the area of the orthographic projection of the biomass particle on a plane perpendicular to the direction of the trajectory. However, Equation 2 is too general to describe the drag force that impedes the motion of the biomass particle. Theoretically, the specific drag function for the model must be derived using the basic concepts of physics. 13

Derivation of the Drag Function
Formally, the drag function is defined as the drag force per unit mass of the biomass particle. One can use Equation 2 to derive the drag function, δ(t), via some mathematical manoeuvres. There are several assumptions need to be considered to define δ(t). Firstly, the nitrogen gas flow is assumed to be laminar. Thus, the nitrogen gas continuum flows in parallel layers where there is no disruption in between these layers. Furthermore, the external forces exerted by foreign particles upon the biomass particle motion are neglected. Such an assumption is based on the fact that the foreign particles count in the reactor is much smaller than the number of nitrogen gas particles. 14 Secondly, the biomass particles are assumed to be perfectly round with an even surface. 3,4 At the rudimentary level, δ(t) is generally defined as: Theoretically, the Reynolds number for a perfectly spherical particle flowing through a laminar continuum is given by: where t is the density of the continuum, which is the nitrogen gas. Here, d p is the diameter of the particle, v R is the relative terminal velocity of the particle with respect to the nitrogen gas flow and μ is the molecular viscosity of the gas. 3,4 The volume of the biomass particle in terms of its diameter is given as: The area of the orthographic projection of a biomass particle on a plane perpendicular to the direction of the particle trajectory is given by: Hence, by solving for t in Equation 4 and substituting the solved t, Equations 5 and 6 into Equation 3 yields the drag function: Simplifying the form of Equation 7 gives the following: where d is the diameter of the particle, dp is the mean diameter of the biomass particles and a is the spread parameter of the biomass particle size population. 15 It is possible to practically determine the mean diameter of the biomass particle population in the laboratory. Also, the cumulative distribution data of the biomass particle diameter population can be constructed via multiple particle sieving procedure. Here, the average diameter value can be graphically obtained by determining the value of d p when f (d p ) = e −1 ≈ 0.368. To illustrate, Table 1 shows the outcome of multiple sieving technique to screen the biomass size starting from the smallest size. At the same time, one must also record the corresponding mass fraction of the specified particle size. Hence, the mean value of the particle diameter is determined graphically on the cumulative probability density plot of Table 1.

Determination of Biomass Particle's Equation of Dynamics
Vitally, one may introduce a force balance equation so that the biomass particle motion can be analysed further. Such a force balance equation must equate the biomass particle inertia with the forces acting on it. 3,4 There are various types of force which act upon the biomass particle. The primary force is the resultant buoyant force, F B / , experienced by the biomass particle due to the nitrogen gas flow and vice versa. The former type of buoyant force is given by gV t where t is the density of the nitrogen gas and V is the volume of the nitrogen gas displaced by the biomass particle. The similar theoretical concept also applies to the buoyant force acted by the biomass particle upon the nitrogen gas. The latter type of buoyant force is given as gV p t where p t is the density of the biomass particle and V is the volume of the biomass particle displaced by the nitrogen gas. Henceforth, the resultant buoyant force is given by: The buoyant force per unit mass of the biomass particle is given by: The secondary force is the drag force between the biomass particle and the nitrogen gas where it has been mathematically described using the drag function, Finally, the external forces, F ext , acting on the biomass particle are not neglected. These external forces include the impulsive force felt by the biomass particle when it hits the reactor wall. Likewise, the forces due to the collision with other foreign particles are also considered as the external forces. Using Newton's second law of dynamics: Combining the forces acting on the biomass particle gives: Re It is important to note that the particle's density changes with time as the pyrolysis reaction progresses. Henceforward, the particle's density now is a function of time, ] g g (18) Therefore, Equation 17 can be simplified into: If the biomass particle moves ideally without any external forces, f ext must be set to a zero vector. Conveniently, for a biomass particle which is free from external forces, its equation of dynamic is given as: Using numerical methods, Equation 20 can be solved to obtain the velocity function of the biomass particle.

NUMERICAL RESULTS AND DISCUSSION
The dynamic behaviour of a biomass particle is often characterised by its position, momentum and the forces acting on it. The main objective of this section is to verify the discrete phase model of fast pyrolysis process to a test case involving a simple EFR geometry. Also, the change in biomass particle position inside an EFR is determined by the amount of the momentum transfer between the biomass particle and the nitrogen gas continuum. Plus, the momentum transfer between the biomass particle and some other foreign particles may also cause such a change. Figure 3 illustrates the simulated plot of the char mass fraction of the pyrolised biomass particle. The simulated plot of the char mass fraction is based on the reaction kinetic theory of biomass pyrolysis which will be reviewed later in this text. In Figure 3, the char mass fraction is interpreted by the coloured particle tracks. As the biomass particle is transported along the positive axial direction by the nitrogen gas, it underwent a few reactions which are described by the two stages semi-global model of fast pyrolysis reaction of wood. 8,12,16,17 The particle tracking result is given in Figure 3. It shows that as the particle moves in the positive axial direction, the production of char increases. The char is fully produced before the actual path of the biomass particle enters the EFR outflow. 18,19 Figure 4 illustrates the simulated char mass fraction against time.  From the plot in Figure 4, the char mass fraction is monotonically reduced according to the exponential trend. Such a curve trend of mass fraction, W, can be derived from the following rate equation:

Products Formation
where k is the reaction factor and is the order of the reaction. Here, Equation 21 can be reduced to a first-order differential equation by simply setting n = 1 and it has the solution of the form: where W 0 is the initial mass fraction. The value of k for each different types of biomass feedstock can be obtained via the static thermogravimetric analysis (TGA) of the feedstock. Most importantly, the reaction factor k can be theoretically determined by the use of Arrhenius equation: where A is the frequency factor, E a is the reaction activation energy, k B is the Boltzmann constant and T is the reaction temperature. 3,7 In this numerical test, the value of A and E a for wood is obtained from the reference value given by Sinha et. al. which is A = 10 8 s -1 and E a = 125.4 kJ mol -1 . 12 The reaction temperature, T, is chosen to be 600K, which is the optimum value being used by USM Energy Laboratory for a fast pyrolysis process.

Particle Tracks Analysis
The transport phenomenon of the simulated biomass particle is first analysed so that the factors that affect fast pyrolysis process yield can be understood.
At t = 0 s, the initial velocity of the biomass particle is 2 m s -1 heading towards the positive x-direction. This value is equal to the biomass particle injection velocity. At this stage, the biomass particle moves in the direction which is perpendicular to the course of the nitrogen gas flow. Thus, the biomass particle motion is distorted by the continuum flow. The simulation result shows that when t < 0.02 s, the particle bounces several times within the confining reactor walls. At 0.02 s < t < 0.06 s, the biomass particle tries to stream through the nitrogen gas laminar flow. Within this period, the biomass particle requires time to "relax" itself until it attains the same velocity as the nitrogen gas phase velocity. Such a timespan is also known as the biomass particle relaxation time, τ p . Formally, τ p is a typical time for a biomass particle response in a flow, or in other words, it is the typical timescale of the particle's reaction to changes in the nitrogen gas phase velocity. For t > 0.06 s, the velocity magnitude of the biomass particle blows up due to the escape boundary wall. On the other hand, the biomass particle has escaped through the outflow of the EFR.
The result of the above velocity plot is the solution of Equation 19. One can possibly solve Equation 19 using numerical methods to obtain the velocity profile of the biomass particle by considering the three Euclidean basis components of the differential equations. In the three-dimensional Euclidean space, ℝ 3 , the three components of the differential equations are given as follows: Note also that Equations 24 and 26 do not have the term (1−ρ/ρ p )g since the gravitational acceleration in the direction of x-axis and y-axis is zero. In addition, the velocity of the nitrogen gas is divided into three components and can be expressed as: Here, U is a velocity vector field of the nitrogen gas flow. By assuming that the solution of the three differential equation as v x (t), v y (t) and v z (t) and correspondingly, the velocity of the biomass particle can be stated as: v Also, the velocity magnitude of the biomass particle can be easily calculated using: The instantaneous momentum of the biomass particle can be calculated using: Crucially, the biomass particle mass changes over time since they undergo reactions whilst being transported by the nitrogen gas continuum. Theoretically, m(t) is related to the reaction kinetics of the pyrolysis reaction and m(t) has the form of Equation 22. Plus, the reaction kinetics must be consistent with the two stages semi-global model of wood as proposed by this work. 3,4,14

CONCLUSION
In summary, the main goal of this work is to propose a numerical model describing the dynamic transport of the biomass particle inside an entrained flow reactor during the fast pyrolysis process. Based on the discrete phase modelling of the fast pyrolysis process, this work proposed a dynamics equation (given in Equation 20) which describes the dynamics of the biomass particle inside the EFR. To achieve a better and accurate solution of this differential equation, one must include all possible external forces which are acting on the biomass particle inside the EFR.
With the currently available mathematical methods and techniques, the dynamics equation can become a handy tool for predicting the dynamic behaviour of biomass particles. However, there is no specific analytical method exists in an attempt to solve the differential equation since the equation can be in any forms, depending on the biomass particle material characteristics and its surrounding condition. From the numerical result of a test case involving the pyrolysis of wood, it is feasible to simulate the fast pyrolysis process using the discrete phase approach. With this tool, it enables us to estimate the optimum pyrolysis operating conditions that give maximum pyrolysis process yield.