Fluidized bed systems are used in many commercial processes; including mixing particles, particle separation, solids coating, power generation through combustion or gasification and particle drying. A proper understanding of the hydrodynamic behavior of the fluidized bed system is needed for the design and to scale up the new efficient reactor  . Inter-phase forces cause the interaction and coupling between phases. The coupling is achieved through the inter-phase forces such as the drag force. The drag model plays an important role in the two-phase flow modeling   . Different correlations for the drag coefficient are available in various literatures. An interested reader can gain adequate knowledge of the subject by studying the literature review of  . Most researchers used experimental and semi-empirical approaches to find correlations for the drag force. However, computational fluid dynamics has become a useful tool for evaluating the two-phase interaction and hydrodynamic behavior.
New insights have been proposed in recent literature based on the combination of CFD method and empirical correlations. Sun Liyan et al.,  proposed a drag coefficient model for gas-mono-size particle flow to predict the interaction between the gas and solid phases in the tapered bubbling fluidized beds. They predicted distribution of concentration, velocity and moments of particles in the tapered bubbling fluidized bed. Wang Shuai et al.,  proposed a modified cluster structure-dependent drag model. They verified their drag model by CFD simulation and concluded that it was able to capture the axial heterogeneity with the dense bottom and dilute top sections. Hua et al.,  proposed a simple drag model in their work to address the critical role of the particle shape to determine the drag force. This drag model uses the Ganser correlation and Ergun correlation. They calibrated particle sphericity by minimum fluidization velocity and the corresponding voidage. Their results showed that the proposed drag model with the estimation method of particle sphericity is reasonable and convenient for the Eulerian-Eulerian simulation of irregular particles in the fluidized beds. Liu et al.,  presented a new scheme to establish the structural parameters model and they used this model to solve structural parameters based on the available drag models. Their simulation results showed that combination of the structural parameters model with the available structure-based drag model can predict the hydrodynamics for Geldart A and B particles in circulating fluidized beds well.
The approach in which effective correlations for fluidization process were obtained analytically was introduced for the first time by the novel study of Grbavcic et al.,  . By introducing the calculus of variations in their research as an optimization tool, a new approach was created. Grbavcic et al.,  expressed the governing momentum equations for the fluid-particle phase and after replacing the Darcy pressure drop term, they obtained the drag coefficient in terms of voidage. This function was then extremized using the calculus of variations and finally an analytic equation for the fluid velocity in terms of voidage was achieved. Owing to the notability of this achievement, Littman and Morgan,  reviewed the generalities of this approach. They compared relations and results with the experimental data in the form of a new study.
The equation obtained by Grbavcic et al.,  is in excellent agreement with the experimental data for higher bed voidages. But in the voidage range of [0.45, 0.59], the pressure gradients stemming from the experimental method lead to higher drag coefficients in comparison to Grbavcic et al.’s analytical approach that ignored the wall effect. Therefore, deviation is observed between the experimental data and analytical correlations. The present study reveals that such a difference stems from the use of the Darcy pressure drop equation in the process of obtaining the momentum equations. Moreover, neglecting the importance of the linear velocity term in the drag model of the low bed voidage amplifies this deviation.
To fix this problem in this study, the pressure drop equation has been replaced by Ergun’s correlation. For vessels of small diameter, the wall effect becomes important and influences the bed voidage. The frictional pressure drop, proposed by Ergun, is proper  because this correlation consists of a viscous term, and a Reynolds-dependent term to take into account the viscous and inertial effects respectively  .
An acceptable adjustment between the analytical and experimental data is obtained by considering the combination of the first and second order terms of polynomials in the drag equation. The process of obtaining relations in this study includes dividing the fluidization process into two intervals of voidage of [0.45, 0.59] and [0.59, 1] and defining two functions with different criteria for these two intervals. To achieve the required quantities in the first range of voidage, the governing momentum equations for the two-phase of the fluid-solid are written, and the drag term is replaced with the first degree velocity term. Linear velocity term is removed from the drag equation to find the required criteria in the second range of the voidage. In this case, for the first range of voidage, a new equation is introduced but in the second range the equation obtained in the studies of Littman and Morgan,  is written. For equations, verification experimental data is used, and the new assumptions are in good agreement with the experimental data.
2. Equations of Motion
To obtain the required parameters in the fluidization process, the conservation equations of mass and momentum in the vertical transport of fluids and solids with the assumption of no-acceleration motion are used. Also, the motion of fluid and solids is considered to be one-dimensional, laminar and uniform in the steady state condition with no wall shear stress.
2.1. Deriving the Drag Coefficient
The drag coefficient is a dimensionless quantity which is used to calculate the drag force acting on the particles suspended in the fluidized bed. This factor depends on voidage of the particles in the fluidization process.
During the fluidization process, by accelerating the fluid, the voidage of the solids increases and it changes the drag force acting on the solid particles. Therefore, the relationship between the drag coefficient and the voidage of particles is the first important parameter in this research and will be found subsequently.
To find this relationship and for simplification, the one-dimensional conservation equations of mass and momentum for vertical transport of fluid and particles by assuming no acceleration in movement should be expressed and the relation for drag force should be replaced. The conservation equations of mass for two-phases of fluid and solid particles will be generally in form of relations 1) and 2). In these equations, the air moves as the continuous phase and suspended particles in the fluidized bed are described by the Eulerian approach. Accordingly, conservation equations of mass and momentum for the two-phase flow of the particle and air are written in the form of Eulerian-Eulerian, outlined as follows:
Momentum equations representing the motion of the fluid and solid phases are written as relations (3) and (4):
Regardless of the wall shear stress and ignoring the stress of particles on each other, or in other words, assuming the two-way coupling, the momentum equation is reduced and rewritten as follows:
where the average normal pressure is introduced in the form of Equation (6):
In the present study, the drag force is considered a linear function (see Equation (7)). In the research of Grbavcic et al.,  this force is assumed to be a second order polynomial, and the relation they obtained is given in the form of Equation (8).
Davidson  assumed that the drag force in the dense phase with lower voidage is proportional to . It is common to use for the large particles involving higher superficial velocities in the momentum equation.
Current study (7)
Grbavcic et al.  and current study (8)
Now, it is desirable to detect a relation for the drag coefficient from the momentum equation. In order to achieve this aim, the pressure drop equation should be placed in Equation (5). In this study, the pressure drop term is replaced with Ergun’s equation. Ergun’s equation for pressure drop of the solid-fluid phase is introduced in Equation (9). In the research of Grbavcic et al.,  the Darcy pressure drop Equation (10) is used:
Current study (9)
Grbavcic et al.  (10)
where and are constants.
The relation for the drag coefficient in terms of voidage is obtained (Equation (12)) by replacing Ergun’s equation (Equation (9)) in the momentum equation
(Equation (5)) for non-accelerating beds ( ). Similarly, by putting the
Darcy correlation (Equation (10)) in the momentum equation (Equation (5)), the relationship between the drag coefficient and the voidage in the research of Grbavcic et al.,  is obtained:
Current study (11)
Grbavcic et al.  (12)
where superficial velocity is defined as:
After deriving the drag coefficient, it is necessary to extremize this function. In the following sections, the drag coefficient is extremized, the relationship between fluidization velocity and voidage is determined, and the results are discussed.
2.2. Extremizing the Drag Functions
In order to be able to extremize this drag coefficient function (in terms of voidage using the calculus of variations), it is necessary to change it to the shortest path problem. In the fluidization process, the physics of the problem dictates that for different voidages, the drag force acting on the suspended solids should be at their least magnitude. It is physically equivalent to a rope that is attached from both ends to two points. In this situation, its potential energy will stand at the lowest possible level.
To formulate this physical sense mathematically, it is necessary to define the drag coefficient and voidage as dimensionless parameters and then normalize them in the interval of [0, 1]. Doing this on the drag coefficient and voidage leads to emergence of two quantities x and y, in the form of relations (14) and (15).
Since the fluidized bed is addressed between the minimum and maximum fluidity, it is necessary to consider the drag coefficient in the minimum fluidity as the beginning of the fluidization process. This quantity is obtained by considering the voidage coefficient of the minimum fluidity in Equation (11) Equation (12).
Current study (16)
Grbavcic et al.  (17)
Now, after the definition of the normalized quantities, it is necessary to determine the required functional of the calculus of variations. It can be perceived from a mathematical viewpoint that, among all the curves available for the
normalized drag coefficient (y), the area under the curve of is
constant and this relation is considered a restriction for the problem. Then:
The next step is to determine the conditions in which the least possible
lengths of the normalized drag curves occur. This condition is
mathematically defined as Equation (19):
The drag coefficient curve will change the maximum between the first and end point. Therefore, the functional of this problem is defined as Equation (20):
In Equation (20), the parameter λ is called the Lagrange multiplier.
Now, after definition of the functional of the problem, its extremum is obtained using the calculus of variations. Referring to the approach of the calculus of variations of Gelfand and Fomin,  it can be perceived that this problem is isoparametric with one restriction, and its functional must satisfy the Euler-Lagrange differential equation. The Euler-Lagrange equation is defined as the relation (21):
By substituting the relation , which was obtained in Equation (20), into Equation (21) and solving it, the extremum value of y is obtained as:
In the relation , and are constant values and for detecting them, three boundary conditions are required. Relations (23), (24) and (25) are referred to as boundary conditions.
The boundary condition expressed in Equation (25) will be in the form of the following relations by using Equations (14) and (15):
Current study (26)
Grbavcic et al.  (27)
where describes terminal velocity when . Using these boundary conditions, the constant values of , and are obtained:
Current study (30)
Grbavcic et al.  (31)
By substituting the values of y and x into the relation (22), the extremum of the drag coefficient equation in terms of voidage is obtained as Equation (32)  :
This dimensionless drag coefficient is the function of the bed fluidization velocity in the balanced bed situation.
2.3. Determining the Relationship between Fluidization Velocity and Voidage
As mentioned before, another important parameter in the fluidization process is the fluidization velocity for different voidages. To find this parameter, the relations obtained for the drag coefficient in terms of voidage are used. After simplifying relations (11) and (12), the connection between the fluidization velocity and voidage is obtained as Equations (33) and (34):
Current study (33)
Grbavcic et al.  (34)
3. Results and Discussion
In this section, the relations which were mentioned previously are validated with the empirical data. To achieve this purpose, the empirical research data of Grbavcic et al.,  will be accompanied. The test set-up used water as fluid and spherical solid particles with 1.94 in diameter. Also, in this experiment, the initial bed height, bed diameter and the reactor height were 1000 mm, 40 mm and 2.2 m, respectively. According to the result of Grbavcic et al.,  experimental data for are 0.0205, 0.447, and 3.77 and the calculated values of are 0.288, 0.9987, and −0.054, respectively.
The data obtained from the experiment of Grbavcic et al.,  for the drag coefficient in terms of voidage are shown in Figure 1. In this figure, the results of the analytical equations obtained from the present study and the research of Grbavcic et al.,  are illustrated.
As can be seen in Figure 1, in the interval of [0.45, 0.59], the relation obtained in the present study is in excellent agreement with the experimental data. In the voidage range of [0.59, 1], the equation obtained in the investigation of Grbavcic et al.,  predicts experimental data as well. It can be concluded that in the voidage range of [0.45, 0.59], the assumptions of Ergun pressure drop as well as the linear term of superficial velocity in the drag model are more consistent with reality. But in the voidage range of [0.59, 1], the assumptions of Darcy pressure drop and non-linear term of the superficial velocity which are performed by Grbavcic et al.,  in the investigation determine the drag coefficient well.
In the voidage range of [0.49, 0.59], the relation derived in the present study represents a closer approximation of the empirical data and is specified in Figure 2. However, in the voidage range of [0.59, 1], the relation achieved in the research of Grbavcic et al.,  predicts the experimental data in the same way.
In the voidage range of [0.45, 0.59], the measured pressure gradients in the research of Grbavcic et al.,  lead to slightly higher drag coefficients. They concluded that this deviation is due to a small wall friction which was ignored in their considered correlations. For vessels of small diameter, the wall effect becomes important and influences the bed voidage  . The frictional pressure drop proposed by Ergun consists of a viscous term and a Reynolds-dependent term to consider the viscous and inertial effects, respectively  .
Figure 1. Comparison of the drag coefficient equations obtained in current study and by Grbavcic et al.  .
Figure 2. Comparison of the fluidization velocity obtained in current study and by Grbavcic et al.  .
The assumptions of Ergun pressure drop and the linear term of superficial velocity in the drag model are more consistent with reality. However, in the voidage range of [0.59, 1], the assumptions of Darcy pressure drop and non-linear term of superficial velocity, which are used in the investigation of Grbavcic et al.,  determine the drag coefficient perfectly.
 Loha, Ch., Chattopadhyay, H. and Chatterjee, P.K. (2012) Assessment of Drag Models in Simulating Bubbling Fluidized Bed Hydrodynamics. Chemical Engineering Science, 75, 400-407.
 van Wachem, B.G.M., Schouten, J.C., van den Bleek, C.M., Krishna, R. and Sinclair, J.L. (2001) Comparative Analysis of CFD Models of Dense Gas-Solid Systems. AIChE Journal, 47, 1035-1051.
 Gao, J., Lan, X., Fan, Y., Chang, J., Wang, G., Lu, Ch. and Xu, Ch. (2009) CFD Modeling and Validation of the Turbulence Fluidized Bed of FCC Particle. AIChE Journal, 55, 1680-1694.
 Sun, L.Y., Xu, W.G., Liu, G.D., Sun, D., Lu, H.L., Tang, Y.J. and Li, D. (2012) Prediction of Flow Behavior of Particles in a Tapered Bubbling Fluidize Bed Using a Second-Order Moment-Frictional Stresses Model. Chemical Engineering Science 84, 170-181.
 Wang, S., Zhao, G.B., Liu, G.D., Lu, H.L., Zhao, F.X. and Zhang, T.Y. (2014) Hydrodynamics of Gas-Solid Risers Using Cluster Structure-Dependent Drag Model. Powder Technology, 254, 214-227.
 Hua, L., Zhao, H., Li, J., Wang, J. and Zhu, Q. (2015) Eulerian-Eulerian Simulation of Irregular Particles in Dense Gas-Solid Fluidized Beds. Powder Technology, 284, 299-311.
 Liu, W.M., Li, H.Z., Zhu, Q.S. and Zhu, Q.H. (2015) A New Structural Parameters Model Based on Drag Coefficient for Simulation of Circulating Fluidized Beds. Powder Technology, 286, 516-526.
 Grbavcic, Z.B., Garic, R.V., Hadzismajlovic, Dz.E., Jovanovic, S., Vukovic, D.V., Littman, H. and Morgan III, M.H. (1991) Variational Model for Prediction of the Fluid-Particle Interphase Drag Coefficient and Particulate Expansion of Fluidized and Sedimenting Beds. Powder Technology, 68, 199-211.
 Littman, H. and Morgan III, M.H. (2010) Two Sample Applications of a Classical Isoperimetric Problem of the Calculus of Variations to Fluid Mechanical Problems in Fluidization and Spouting. Particuology, 8, 503-506.