Effects of interphase coupling in one-dimensional particle-laden flow

In Euler–Euler particle-laden flow models, the interphase coupling terms might include, in addition to the drag force, the effects of fluid pressure interaction with the particles and the volumetric displacement of the fluid phase due to the presence of the particles. This paper presents an analytical study on these effects for laminar, one-dimensional, steady, inviscid particle-laden flow. The interphase drag force is represented by a linear Stokes formula that is valid for low particle Reynolds number. It is shown mathematically that a full two-way coupling model, which includes all the aforementioned effects, depends on three parameters: inlet particle volume fraction (concentration), inlet density ratio–between the fluid and particle phases–and inlet velocity ratio. When the inlet particle concentration and density ratio are kept low at 0.01 and 0.001, respectively, both the interphase pressure coupling force and the volumetric displacement of the fluid phase are insignificant. For density ratios of 0.001 and 1 and for the inlet particle concentration values examined, the effect of pressure interaction force is not significant. However, it becomes significant for higher density ratio of 10 at inlet particle concentration of 0.1.


Introduction
Particle-laden flows, in which a large number of small dispersed particles co-flow with liquid or gas, are ubiquitous in both engineering and natural processes. When the size of these particles is smaller than the smallest length scale of the relevant fluid motion, it is valid to approximate them as point-particles [1,2]. This approximation allows the fluid phase to be described as a continuous media in the Eulerian framework. The dispersed particle phase, on the other hand, could be modelled using either the discrete Lagrangian approach [3][4][5] or the continuum Eulerian approach [6][7][8][9]. In the former, the equations of motion for each particle are solved in order to obtain its positions and velocities throughout its time of motion whereas in the latter, the particle phase is assumed to behave as a pseudo-fluid so that its variables are analogous to the hydrodynamic ones . The Lagrangian approach requires much more computational resources to calculate the trajectories and velocities of the massive number of particles encountered in realistic flow situations [10]. This led to extensive research effort in the development of Eulerian models [7,[11][12][13][14][15][16] and their application to different flow problems [7,9,[17][18][19][20][21][22].
In particle-laden flows, the continuum modelling of the fluid phase implies that the two-phases co-exist at every point in the domain while accounting for the fraction of volume occupied by each phase. Since the summation of the volume fractions of the two phases ought to be unity, the value of the volume fraction of either phase is used to formulate the governing equations because it determines how the two phases are coupled to each other. The oneway coupling model is often employed when the particle volume fraction (concentration) is very low compared to that of the fluid phase [6,8,9]. This implies that the volume fraction of the fluid phase remains constant at unity. Furthermore, the effect of the hydrodynamic interactions between the two phases, for example the drag force, is negligible for the fluid phase but significant for the particle phase. Thus, one needs to solve the single phase flow equations for the fluid phase and subsequently, use the fluid flow field as an input parameter to obtain the particle phase flow field.
For moderate particle concentration, in which it is still valid to neglect the particle-particle interactions, a twoway interphase coupling approach is necessary. This is because the particles affect the fluid flow field so it cannot be evaluated independently. There are two differences between this two-way coupling model and the one-way coupling model. First, the hydrodynamic interactions have significant effects on the fluid phase. Second, the fluid phase volume fraction cannot be assumed constant at unity, which means that the volumetric displacement of the fluid by the particles should be accounted for [23]. However, in one class of the two-way coupling models, referred to as the dusty gas model [7], the fluid phase volume fraction is assumed constant at unity. This means that the fluid phase is coupled to the particle phase via the hydrodynamic interactions only. In either case, the twoway coupling models require the solution of the governing equations of the two phases simultaneously.
Owing to the complexity of the fluid flow governing equations and the infinite number of possible flow configurations, it is intuitive to claim that the transition from one-way coupling to two-way coupling is not universal. Even for a specific flow problem, it is not possible to determine quantitatively, prior to establishing the problem solution, the values of inlet or initial particle volume fraction, or other flow parameters, at which this transition takes place. The flow configuration could enrich some regions with particles and deplete other regions [23][24][25], resulting in different levels of interphase coupling across the domain. Thus, understanding the effect of interphase coupling requires solving the governing equations of the corresponding models for a specific flow problem then comparing their theoretical predictions. This approach has been used in previous research by employing the Lagrianian framework to simulate the particle phase using computational approaches [23][24][25]. The aim of this paper is, however, to use the Eulerian framework, for different levels of interphase coupling, in the description of a flow problem that can be solved using analytical methods. This analytical approach enables to elucidate the physics of the problem by understanding the functional relations in the models with different levels of interphase coupling.The rest of the paper is organised as follows. Section 2 presents the problem description and analytical solutions of four Eulerian models with different levels of interphase coupling. Section 3 presents some comparisons between the predictions of these models for a wide range of parameters. Section 4 summarises the paper findings highlighting its importance and limitations. Three appendices are used to explain some mathematical derivations of the analytical solutions.

Problem description and analytical solution
This paper considers a one-dimensional particle-laden flow problem in the streamwise direction with a uniform cross-sectional area. In order to focus on the effects of inertia and interphase coupling terms, the flow is assumed inviscid and the gravity effects are neglected for both phases. Furthermore, it is assumed that the slip velocity between the two phases is small, i.e. low particle Reynolds number, so that the linear Stokes drag relation is valid in describing the hydrodynamic interactions between the two phases. In what follows, these assumptions are applied to the governing equations of different Eulerian models representing various levels of interphase coupling in order to derive their analytical solutions.

One-way coupling model
The one-way coupling (1WC) model is used in the literature when the particle concentration (used interchangeably in this paper with "particle phase volume fraction") is too low to affect the gas phase [6,9]. Thus, the governing equations of the fluid phase are decoupled from that of the particle phase, which means they can be solved independently to obtain the fluid phase flow field. The particle phase, on the other hand, is coupled to the fluid phase so that its governing equations are [6] where is the particle phase volume fraction, p is the particle material density, u pi is the particle phase velocity, u fi is the fluid phase velocity, and is the reciprocal of the particle relaxation time which is constant for low particle Reynolds number and given by [6,9] (1) t p + x j p u pj = 0, where is the dynamic viscosity of the fluid and d is the particle diameter. For one-dimensional steady flow, the continuity equation of the fluid phase implies that its velocity is constant at U f (the subscripts denoting the coordinates are dropped). Moreover, the continuity equation of the particle phase is Thus, the non-conservative form of the particle phase momentum equation in the x direction is In order to write Eq. (5) in the dimensionless form, where the star superscripts denote dimensionless variables, the following provisional scaling is introduced which is substituted in Eq. (5) to give Defining the velocity ratio B as where U p is the particle phase inlet velocity. By choosing the length scale X as the pre-factor of the left-hand side of equation (5) becomes unity. Thus, the dimensionless particle phase momentum equation now reads Because the flow is one-dimensional, the particle phase continuity equation can be written as where 0 is the inlet particle phase concentration which is used in the scaling as the continuity equation in the dimensionless form is thus The non-dimensional particle phase momentum equation (10) can be integrated to obtain a solution for the dimensionless particle phase velocity [9]. Re-writing Eq. (10) as By using the separation of variables method and applying the boundary condition, the dimensionless particle phase velocity is given by the following implicit relation Equation (16) can be inverted to obtain an explicit relation of the velocity in terms of Lambert W function as [9] The particle phase volume fraction can thus be obtained using the continuity equation (13). For fully developed flow, the particle phase velocity is obtained by setting du p ∕dx = 0 in the momentum equation (10) resulting in

Dusty gas model
The dusty gas (DG) model [7] uses a pseudo two-way coupling formulation. The particle phase momentum equation is the same as the one-way coupling model given by Eqs. (1) and (2). However, the particle phase affects the fluid phase via a momentum source term only, which makes the coupling between the phases two-way. Thus, one needs to solve the governing equations for the two phases simultaneously. These governing equations for the fluid phase are [7] where f is the material density of the fluid phase, p is the fluid phase pressure, ij is the Kronecker delta, and ij is the viscous stresses which is neglected in this paper.
It should be noted that the continuity equation (19) of the fluid phase is identical to that of single phase flows, which results in a uniform fluid velocity for the one-dimensional incompressible flow problem under investigation. Thus, the momentum equation of the particle phase is the same as that of the one-way coupling model (5) resulting in the same particle velocity as (17). Moreover, the left-hand side of the fluid phase momentum equation is zero. Thus, the fluid phase momentum equation (20) becomes a function of the pressure only as Using the same scaling for the one-way coupling model, given by Eqs. (6) and (9), and a provisional scaling for the pressure, p = p * P 0 , where P 0 is a reference pressure, the momentum equation (21) in the dimensionless form is which shows that, for incompressible flow is a natural scaling of the fluid phase pressure, because it is divided by the particle phase inlet kinetic energy. Thus, the final dimensionless form of the fluid phase momentum equation is It is possible to substitute (17) into (24) to obtain a differential equation describing the pressure as a function of the dimensionless distance x * . However, it is more informative to use the particle phase velocity as the independent variable in the pressure Eq. (24). Thus, the chain rule is used to express the pressure gradient with respect to the particle phase velocity as By substituting Eqs. (13), (14) and (25) into Eq. (24), it can be shown that the variation of the pressure is constant with respect to the particle velocity as Using the boundary condition shown in Eq. (15) for the particle velocity and the following boundary condition for the pressure the solution is thus given by Note that the fluid pressure obtained here is different from that of the 1WC model because the latter is zero as its fluid flow velocity is constant. For fully developed flow, substitute Eq. (18) into Eq. (28) so that the pressure is

Two-way coupling model
For a full two-way coupling (2WC) between the particle and the fluid phases [8], the continuity equation (1) of the particle phase, used in both the one-way coupling and the dusty gas model, remains intact. On the other hand, the momentum equation of the particle phase is given by [8] The governing equations for the fluid phase are given by [8] One difference between the 2WC model and the DG model is that the volumetric displacement of the fluid due to the presence of the particles is accounted for in the 2WC model. In the continuity equation (31) and the inertia terms on the left-hand side of the momentum equation (32), the fluid density is multiplied by the fluid phase volume fraction, (1 − ) . Furthermore, the two phase are coupled via a pressure interaction term, p ij ∕ x i . In one-dimension, the inviscid 2WC model is Particle phase continuity Fluid phase continuity

Particle phase momentum
Fluid phase momentum Integrating the particle phase continuity equation (34) gives Expanding the left-hand side of the particle phase momentum equation (35) and substituting the continuity equation (37), the particle phase momentum equation can be written as Using the scaling defined previously as Using the scaling given by (12) and (39) and after some manipulation, the fluid phase dimensionless momentum equation is Applying the scaling given by (12) and (39)  In order to solve this 2WC model, the continuity equations of the two phases, (13) and (46), are used to obtain an algebraic relation between u p * and u f * by eliminating * as The first derivative of Eq. (47) is Substituting u f * , given by Eq. (47), into the particle phase momentum equation (40), the fluid phase pressure gradient is expressed as Using the integrated continuity equation of the particle phase, * = 1∕u p * , and Eq. (47), the fluid phase momentum (44) can be written as Substituting du f * ∕dx * , given by Eq. (48), into Eq. (50), the following equation is obtained Substituting dp * ∕dx * given by Eq. (49), into Eq. (51), the following differential equation (see Appendix A) is obtained which is solved for u p * with the boundary condition u p * (0) = 1 using the separation of variables method (Appendix B) resulting in the following algebraic expression where Note that for the limit 0 → 0 , the first two terms on the left-hand side of Eq. (53) cancel each other and 3 = B . Thus, the whole equation is identical to the solution of the DG model given by Eqs. (17) and (28). An explicit solution of the particle phase velocity requires inverting the left-hand of Eq. (53) which does not seem possible.
To derive the solution of the fluid phase pressure, the chain rule is applied to Eq. (49) in order to express the pressure gradient with respect to the particle phase velocity, in a similar fashion to Sect. 2.2, as Substituting du p * ∕dx * , given by Eq. (52), into Eq. (59) gives which after some manipulation (Appendix C) and substituting A from Eq. (45) becomes then solving using the separation of variables method gives where c is the constant of integration, which is evaluated using the boundary condition given by Eq. (27). Thus, the pressure field as a function of the particle phase velocity is Note that the fluid pressure of the 2WC model equals to that of the DG model shown in Eq. (28) multiplied by 1 − CB 2 1 − 0 ∕ u p * − 0 . The fully developed particle phase velocity is obtained by setting its gradients to zero in Eq. (52) resulting in which shows that the fully developed velocity of the 2WC equals to that of the 1WC and DG plus 0 (1 − B) . It is thus expected that the discrepancy between the 2WC and the 1WC/DG increases with 0 . The fully developed pressure is obtained by substituting Eq. (64) into Eq. (63) as which shows that the fully developed pressure of the 2WC equals to that of the DG model multiplied by

Two-way coupling without pressure interaction
In order to obtain an explicit solution for the particle phase velocity, u p * , it is possible to approximate the 2WC by neglecting the pressure interaction term, 0 (dp * ∕dx * ) , dp * du p * = −1 + For fully developed flow, the particle velocity is obtained by setting its gradient to zero in Eq. (69) to obtain which is the same as the full two-way coupling model shown in Eq. (64).
In order to derive the solution for the pressure, substitute the continuity equation (13) and Eq. (47) into the momentum equation of the fluid phase (67) as  45), results in the following differential equation describing the pressure as a derivative of the particle phase velocity which is identical to Eq. (61) of the 2WC model. Thus, the fluid pressure of the 2WCNP is given by Eq. (63), and the fully developed fluid pressure is given by Eq. (65) because the fully developed velocity is the same as that of the 2WC. Moreover, since the fully developed particle phase velocity (73) and the continuity equation (47) are the same as the 2WC model, the fully developed gas phase velocity is the same as Eq. (66).

Results and Discussion
In what follows, the predictions of the dimensionless particle phase velocity and the fluid phase pressure are assessed using the analytical solutions of the models derived in Section 2 and summarised in Table 1. Because the particle phase velocity of the 1WC is identical to that of the DG and the fluid phase pressure of the 1WC is constant at zero (inlet value), the 1WC model is not included in the following analysis. Thus, the focus here is on three models: DG, 2WC and 2WCNP. For the 2WC model, Table 1 shows that both the particle phase velocity and the fluid phase pressure depend on three dimensionless parameters: the inlet velocity ratio, B, the inlet particle concentration, 0 , and the density ratio, C. The velocity obtained using 2WCNP does not depend on the density ratio, C, though its pressure depends on it. Both particle phase velocity and fluid phase pressure of the DG model depend on the velocity ratio, B, only. In order to compare the aforementioned models over an applicable range of parameters, the predictions are obtained for 12 cases using the conditions shown in Table 2. For the density ratio, 0.001, 1 and 10 are used. For the particle concentration, 0.01 and 0.1 are used to represent particle flows of dilute and moderate concentration, respectively. The velocity ratio is kept small at 2 and 0.5, because the underlying assumption of the model is low slip velocity between the two phases. The implicit form of the particle phase velocity for the 2WC is plotted by using the value of the velocity obtained from the 2WCNP model as the independent variable and calculating the dimensionless distance x * as the dependent variable. Figures 1 and 2 show the particle phase velocity and the fluid phase pressure profiles, respectively, for a density ratio C of 0.001. In all cases and for all models, both velocity and pressure change with the distance until reaching the fully developed values. For a low particle concentration of 0.01 (cases 1 and 2), the velocity and pressure predictions of all three models are identical. This behaviour is expected because the DG model is a limiting case, where 0 → 0 , of both 2WC and 2WCNP models.
When the particle concentration is increased to 0.1 (cases 3 and 4), the predictions of the 2WC and 2WCNP models remain identical, suggesting that the additional terms in Eq. (53) have negligible effects, for the low density ratio of 0.001, compared to Eq. (72) or (71), even for high particle concentration of 0.1. For cases 3 and 4, the predictions of both 2WC and 2WCNP deviate from that of the DG model; suggesting that a value of inlet particle concentration of 0.1 is not small enough for the DG model to be accurate. For each case presented in Figs. 1 and 2, it is observed that the pressure trend is consistently opposite to that of the velocity; the fluid pressure decreases when the particle velocity increases and vice versa. This opposite sign of the velocity and pressure gradients is expected for the DG model; because dp * ∕du p * is always negative as shown in Eq. (26). For both 2WC and 2WCNP, Eq. (61), which is identical to Eq. (77), shows that dp * ∕du p * is negative when CB 2 1 − 0 2 ∕ u p * − 0 2 < 1 , which is expected because C is very small at 0.001. Increasing the density ratio to 1, as shown in Fig. 3, does not affect the behaviour of the velocity predictions of the three models, in a similar fashion to Fig. 1. All three models are identical for the low particle concentration of 0.01, and the DG model deviates, for a particle concentration of 0.1, from both the 2WC and 2WCNP. Figure 4 shows the fluid phase pressure profiles for a density ratio of 1. For both values of particle concentration (0.01 and 0.1), there is a significant discrepancy between the predictions of the DG model and that of the 2WC/2WCNP models. When B = 2 (cases 5 and 7), the pressure trends of DG and 2WC/2WCNP are opposite, with the former decreasing monotonically and the latter increasing monotonically. This is because the term CB 2 1 − 0 2 ∕ u p * − 0 2 in Eq. (61), which dictates the pressure sign, becomes greater than unity for the relatively large density ratio, C = 1 . This makes the fluid pressure gradient with respect to the particle velocity positive, i.e. they have the same trend. When B = 0.5 (cases 6 and 8), there is a discrepancy between the 2WC/2WCNP and the DG, though all three models have the same increasing trend which is opposite to that of the particle velocity. This is because the term CB 2 1 − 0 2 ∕ u p * − 0 2 in Eq.  (61) become less than unity, because B 2 is small at 0.025 compared to 16 in cases 5 and 7. Figures 5 and 6 show the velocity and pressure predictions, respectively, using a density ratio C of 10. The velocity predictions of all models are nearly identical, when using a low inlet particle concentration of 0.01, as cases 9 and 10 in Fig. 5 show. In contrast to the previous figures, there is a pronounced discrepancy, for inlet particle concentration of 0.1, in the velocity predictions of the 2WC model compared to that of the 2WCNP model, as cases 11 and 12 in Fig. 5 show. This suggests that the additional terms of Eq. (53) compared to Eqs. (71) or (72) cannot be neglected for this value of density ratio. For an inlet particle concentration of 0.01 and velocity ratio of 2, the pressure predictions of both the 2WC and 2WCNP are identical, and they deviate from that of the DG model, as shown in case 9 of Fig. 6. Here, the discrepancy between the 2WC/2WCNP and the DG is greater than that of case 5 shown in Fig. 4. This is because, for both the 2WC and 2WCNP models, the pressure is a function of the density ratio C, as shown by Eq. (63), whereas the pressure of the DG model does not depend on C, as shown in Eq. (28). In case 10 of Fig. 6, where the velocity ratio is 0.5, the discrepancy between 2WC/2WCNP and DG is also clear. A main difference in the 2WC/2WCNP model predictions here, compared to case 6 of Fig. 4, where the density ratio is 1, is that in the former the pressure decreases whereas in the latter it increases. Again, this is traced to the additional term in Eqs. (61) and (77), CB 2 , which becomes greater than unity for large C, making the pressure gradient with respect to the particle velocity positive. In a similar fashion to the velocity profiles shown in cases 11 and 12 of Fig. 5, there is a pronounced discrepancy between the pressure predictions when using higher inlet particle concentration of 0.1 in cases 11 and 12 of Fig. 6. This is because, despite having identical pressure equations for both the 2WC and 2WCNP models, the pressure is a function of the velocity, which depends on C for the 2WC model and does not depend on C for the 2WCNP model.

Conclusions
This paper presented an analytical study concerning the effects of the interphase coupling on the predictions of onedimensional steady particle-laden flow models in the Eulerian formulation. Four models involving different levels of interphase coupling were analysed: one-way coupling (1WC), dusty gas (DG), two-way coupling with (2WC) and without (2WCNP) interphase pressure interaction. It was shown mathematically that the 2WC model depends on three parameters: the inlet particle volume fraction, inlet density and velocity ratios between the fluid and particle phases. Because the DG model is a special case of the 2WC/2WCNP, where the inlet particle concentration vanishes, its predictions are identical to that of the 2WC/2WCNP when the inlet particle concentration and density ratio are kept low at 0.01 and 0.001, respectively. This confirms that both the interphase pressure coupling force and the volumetric displacement of the fluid phase are insignificant for these values. For density ratios of 0.001 and 1 and for all particle concentrations examined, the effect of pressure coupling force, which was evident by comparing the 2WC and 2WCNP, is not significant. It becomes significant, however, for higher density ratio of 10 at inlet particle concentration of 0.1.The flow problem studied in this paper is fairly simple-realistic particle-laden laminar flows are often multi-dimensional with more complicated forms of interphase interactions, i.e. the drag force is not within the Stokes regime. However, it is necessary to study such simple configuration using analytical methodology as it provides some understanding of more complicated geometries and physics.