Effective lateral dispersion of momentum, heat and mass in bubbling fluidized beds

Abstract

The lateral dispersion of bed material in a bubbling fluidized bed is a key parameter in the prediction of the effective in-bed heat transfer and transport of heterogenous reactants, properties important for the successful design and scale-up of thermal and/or chemical processes. Computational fluid dynamics simulations offer means to investigate such beds in silico and derive effective parameters for reduced-order models. In this work, we use the Eulerian-Eulerian two-fluid model with the kinetic theory of granular flow to perform numerical simulations of solids mixing and heat transfer in bubbling fluidized beds. We extract the lateral solids dispersion coefficient using four different methods: by fitting the transient response of the bed to (1) an ideal heat or (2) mass transfer problem, (3) by extracting the time-averaged heat transfer behavior and (4) through a momentum transfer approach in an analogy with single-phase turbulence. The method (2) fitting against a mass transfer problem is found to produce robust results at a reasonable computational cost when assessed against experiments. Furthermore, the gas inlet boundary condition is shown to have a significant effect on the prediction, indicating a need to account for nozzle characteristics when simulating industrial cases.

Keywords

effective dispersion / heat transfer / mass transfer / mixing / gas-solid fluidized bed

1 Introduction

Bubbling fluidized beds are well-established in a wide range of industries due to their excellent ability to thoroughly mix solids and gases and provide good heat and mass transfer [1]. However, the scale-up of fluidized bed processes, which for solid fuel combustion has yielded units with cross-sectional areas above 100 m², entails the intrinsic challenge of ensuring a relatively uniform lateral distribution of solids reactants (such as catalyzers, emission sorbents, etc.), as well as temperature for optimal performance. Recently, the possibility offered by fluidized bed units to combine generation of heat, power and chemicals in polygeneration schemes has attracted increasing interest [2−9]. One possibility would for example be to use a single fluidized bed divided into different sections [7], where exothermic reactions in one section could drive endothermic reactions in another [10]. This requires a detailed understanding of the heat and mass transfer in the bed, to enable optimal design and operation [11].

Computational fluid dynamics (CFD) simulations represent a well-established tool for detailed investigations of transport phenomena in fluidized beds [10−15]. However, efficient scale-up of existing processes and design of novel ones will also necessitate viable reduced-order models. Models where chemical kinetics and heat and mass transfer are described through effective parameters, to some extent can substitute comprehensive descriptions of the transient motion of solids, and their overwhelming computational cost [16−24]. Based on experimental evidence showing the lateral mixing of solids and gas to follow a diffusion-like pattern at the macroscopic scale (e.g., [25]), in such reduced-order models for fluidized beds, convective and diffusive transport phenomena in the lateral direction are typically lumped together, so that the lateral motion of the bed material can be described using an effective solids lateral dispersion. Effective retrofitting of fluidized bed reactors would thus benefit immensely from a possibility to derive such effective parameters for reduced-order models in silico using detailed CFD simulations.

As the heat and mass transfer in a bubbling fluidized bed is mainly governed by the motion of the bed material, the challenge of understanding and quantifying the heat and mass transfer becomes intimately linked to understanding and quantifying the momentum transfer in the bed. It is thus clear that, in principle, it should be possible to derive an effective solids dispersion coefficient from either a momentum, mass or energy basis, and that one would expect the three approaches to return similar values. CFD simulations offer the possibility to perform such derivations in parallel, without interference between methods, based on the same underlying prediction of the solids mixing. However, to the best of our knowledge, no such comprehensive analysis has been performed on the same numerical data set, although studies employing either one of the three individual approaches (momentum [13], heat [26] or mass [14,15,27−31] transport) are available.

In this work, we perform CFD simulations of two different bubbling fluidized beds. Our primary objective is to assess to what extent the extracted value of the effective lateral solids dispersion depends on whether the analysis is performed using the solids velocity field, the temperature field or a passive tracer concentration field as its basis. We also investigate the influence of the size of the computational domain and the representation of the gas inlet device on the prediction. The accuracy of the numerical results is established by comparisons to experimental measurements from both beds. Finally, we conclude by discussing the differences in accuracy, computational cost and data storage requirements for the investigated dispersion analysis methods.

2 Theory

2.1 The two-fluid model (TFM)

This study makes use of the Eulerian-Eulerian framework in fluidization commonly called the TFM. In the TFM, both phases are treated as interpenetrating continua with their own separate but analogous conservation equations. In this setup the interfacial terms as well as the solid phase stress tensor are unclosed and need to be modeled. Here, these models are provided by the kinetic theory of granular flow (KTGF), as well as an appropriate drag law.

The complete framework of equations used, including the drag law and the submodels for the KTGF, is provided in the Electronic Supplementary Material (ESM). This framework has been used in previous studies to estimate particle dispersion in conjunction with a particle tracking routine [10] and to estimate particle dispersion [32]. Apart from the drag law, the framework is also similar to that used by Oke et al. to estimate particle dispersion [14]. Although the specific combination of closure models used could well influence the result, this has not been investigated further in the present study as the focus lies in the post-processing step of the simulation.

2.2 Heat transfer and solids mixing in fluidized beds

Lateral heat transfer in dispersive flows can be described by means of an effective thermal conductivity k_eff as shown in the left-hand side of Eq. (1), where Q is the heat flow in the horizontal direction through a certain vertical cross-sectional area A and ∂T/∂x is the temperature gradient in the normal direction to A.

(1)$k_{\mathrm{e}\mathrm{f}\mathrm{f}}=-\frac{Q/A}{\mathrm{\partial }T/\mathrm{\partial }x}=k_{\mathrm{s}\mathrm{t}\mathrm{a}\mathrm{t}\mathrm{i}\mathrm{c}}+k_{\mathrm{m}\mathrm{i}\mathrm{x},\mathrm{g}}+k_{\mathrm{m}\mathrm{i}\mathrm{x},\mathrm{s}}.$

The resulting thermal conductivity can be assumed to be the result from the contribution of three different factors, as depicted in the right-hand side of Eq. (1): a static contribution k_static (dependent on the volumetric fraction and intrinsic properties of the bed material and gas) and contributions induced by the mixing of gas and solids (k_mix,g and k_mix,s), respectively. These latter two can, through Eq. (2), be expressed as a function of their respective dispersion coefficients D_L,i, heat capacity c_p,i, density ρ_i, and volume fraction α_i:

(2)$k_{\mathrm{m}\mathrm{i}\mathrm{x},i}={\alpha }_i{\rho }_i{c}_{\mathrm{p},i}{D}_{\mathrm{L},i}.$

The contribution of the gas mixing to the total effective thermal conductivity is usually negligible compared to the contribution of the solids mixing because of the high solids-to-gas density ratio. Additionally, the ratio $k_{\mathrm{s}\mathrm{t}\mathrm{a}\mathrm{t}\mathrm{i}\mathrm{c}}/k_{\mathrm{m}\mathrm{i}\mathrm{x},\mathrm{s}}\approx D_{\mathrm{t}\mathrm{h},\mathrm{s}}/D_{\mathrm{L},\mathrm{s}}$, with D_th,s being the solids thermal diffusivity. As D_th,s $\le$ 10⁻⁴ m²·s^–1 and typically D_L,s $\gg$ 10⁻⁴ m²·s^–1, it can also be stated that the contribution of the solids mixing is much higher than the contribution of the static term. Thus, the sensible heat transported laterally in fluidized beds is dominated by the solids mixing, and therefore it is reasonable to assume that k_eff ≈ k_mix,s, and thus:

(3)$D_{\mathrm{L},\mathrm{s}}\approx k_{\mathrm{e}\mathrm{f}\mathrm{f}}/{\alpha }_{\mathrm{s}}{\rho }_{\mathrm{s}}{c}_{\mathrm{p},\mathrm{s}}.$

3 Method

Two works using Geldart group B solids are taken as experimental reference: one employed finer material under high fluidization velocity [33]—hereon referred to as high voidage case, HiVoid—and the other using coarser bed material at lower fluidization velocity—hereon referred to as low voidage case, LoVoid. Tab.1 shows a breakdown of the different simulation cases used in this study. All CFD simulations are performed in three-dimensional (3D) as it is well established that the dimensionality plays a dominant role in predictions of the lateral solid mixing, and that two-dimensional simulations therefore are less reliable [34]. The grid resolution used was optimized to 6 and 15 mm for the HiVoid and LoVoid cases respectively, based on the observed convergence of the predicted time-averaged bed height and bubble sizes and trajectories for resolutions of ≤ 30 particle diameters. The time step used was 10⁻⁴ s.

Tab.1 The simulation cases evaluated in the present work

Case	dp/mm	u0/umf	H0/m	Lz/m	Simulation study
					Lx/m	Inlet boundary condition
HiVoid-small-uniform	0.3	33.7	0.2	0.6	0.6	Uniform
HiVoid-large-uniform					1.2	Uniform
LoVoid-small-uniform	0.9	2.3	0.45	1.29	0.99	Uniform
LoVoid-large-uniform					1.995	Uniform
LoVoid-small-nozzle					0.99	Nozzles

The simulation domain is periodic in the z-direction to avoid wall effects, and on the sides in the x-direction a symmetry condition is imposed to minimize wall effects yet prevent particle movement. Two domain sizes are used for both beds to evaluate the influence of the domain size on the simulated solids dispersion. In general, one would expect that the estimated dispersion value is independent of the choice of domain size once the domain is large enough; that is, when the domain is significantly larger than the bubble size [35].

A pressure-outlet boundary condition is applied above the bed. Two different types of inlet boundary conditions for the gas are used: a uniform condition (i.e., the conventional “porous plate” boundary condition) and a more elaborate boundary condition mimicking the effect of industrial nozzles. In the latter boundary condition, approximately 94% of the bottom plane is a no-slip wall, with the remaining area representing separate inlet channels for the gas. These gas inlets are spaced as nozzles in an industrial bed. At the bottom of the gas inlet channels, a pressure-inlet boundary condition drives a flow through the channel so that the average fluidization velocity attained is constant and equal to that specified in the corresponding uniform case (Tab.1). The solid phase velocity is set equal to zero in the inlet channels to prevent bed material from falling out of the bed. This case offers the possibility to study the influence of the velocity boundary condition on the attained solids dispersion.

The bed is initialized as static and made to fluidize. After fluidization is achieved and observed to be pseudo-steady (as judged from observing the pressure drop over the bed), the evaluation is started. The start of the evaluation consists of applying a marker function to half of the domain, marking the bed solids in that region. In addition, two source terms are added to the energy equation in the outer parts of the domain. The source terms maintain a constant bed temperature in those regions to provide a set-up analogous to heat transfer through a slab with fixed temperatures at the ends. The setup is shown in Fig.1, which depicts the marked bed regions and the frozen temperature fields. Note that the solution to the energy equation has no influence on the bed motion, and also the aforementioned marker function provides no backreaction to the bed.

Fig.1 The domain with (a) imposed initial and (b) boundary conditions for methods 2 and 1, respectively.

Full size|PPT slide

Solution data are exported continuously during the simulation. For most cases the sampling interval was 0.05 s—this sampling frequency was motivated by a power spectral density (PSD) analysis, as presented in the ESM. The exception was the LoVoid-large-uniform case where additional data was exported on two planes at different elevations in the domain. On those two planes, data was exported at a frequency of 400 kHz. After data collection, the time-resolved data needs to be post-processed further to yield an effective solids dispersion value. In the present work, four different methods for this post-processing were tried.

3.1 Method 1—Fitting a heat conduction value from a transient response (CondTrans)

A straightforward way to quantify the solids lateral dispersion is to calculate it via the effective conduction of the bed (cf. Section 2.2). The effective conduction can be obtained by comparing the transient heat transfer response in the bed to an analytical system in which the conduction is known (a slab with fixed temperatures at each end, see Fig.1(b)) and the temperature profile, T^*, can be obtained analytically. For the 3D-CFD temperature field to be compared to the one-dimensional (1D) analytical solution, averaging is required. The average is twofold, both along the z-axis and the y-axis, the latter average being a weighted mean with the solid volume fraction as weight. Details on the calculation of T^* and the averaging procedure are provided in the ESM. Given both the analytical temperature profile and the simulated profile, the effective conduction value, k_eff, for the simulated system is then obtained as the value minimizing the loss function:

(4)$L\equiv \left|T^{\ast }\left(x,t;k_{\mathrm{e}\mathrm{f}\mathrm{f}}\right)-\overline{T}\left(x,t\right)\right|,$

where T^* represents the analytical temperature profile and $\overline{T}$ the averaged profile from the simulated system. Having obtained the value for the effective lateral conduction, Eq. (3) can be used to calculate the output of method 1, i.e., the bed-averaged D_L-value that produces the closest resemblance to the transient analytical 1D heat conduction problem (the solution of which is given in the ESM).

3.2 Method 2—Fitting a solids dispersion value from a transient response (DispTrans)

An alternative and more direct route could be to quantify the average lateral dispersion coefficient of the bed through a fitting method where the analytical system instead is initially segregated, such that half of the domain is marked with an indicator function. The concentration profile across the domain for such a system is gradually flattened over time through dispersive mixing following the diffusion equation. The analytical solution to the ideal transient 1D diffusion problem is presented in the ESM. The output of this method is thus the bed-averaged D_L-value that produces the closest resemblance to this analytical solution.

This fitting method was well described by Oke et al. [14], and before them Liu and Chen [32] had used essentially the same procedure, albeit with tracer particles instead of a marker function. Given the ideal solution, the procedure is analogous to that of method 1 in that a loss function is formulated and the simulation data averaged first along the z-axis and then volume-fraction weighted. In essence, the major difference between methods 1 and 2 is the problem setup for the ideal system. In method 1 the initial condition is a constant value across the domain and the boundary values are fixed. In method 2 the initial condition is a step function, and the boundary values are free to slide.

3.3 Method 3—Getting a heat conduction value from the average heat flux at pseudo-steady-state (CondSteady)

A third estimate of the solids lateral dispersion can be gained from the pseudo-steady-state heat conduction behavior of the simulated beds. It is known from Fourier’s law that the temperature profile in a slab with fixed end temperatures should tend to be a straight line at steady-state. A conduction value can then be obtained algebraically via $k_{\mathrm{e}\mathrm{f}\mathrm{f}}=-q_x\mathrm{\Delta }x/\mathrm{\Delta }T$, and recast into a dispersion via Eq. (2). Here, q_x is extracted from the simulation, and Δx and ΔT are known from the geometry and boundary condition specification. While the evaluation procedure for this method might seem simple, the pseudo-steady-state might not be easily recognized in a simulated fluidized bed: due to the chaotic motion of the bed and the high thermal mass of the sand, rather large fluctuations are to be expected. Thus, for this method, the quantity of interest is the time-averaged heat flux. As with any average, sufficient data needs to be sampled so as to get a stable value. Further, over sufficient long time, a pseudo-steady temperature profile with constant gradient, i.e., constant heat flux, should be expected as the hot and cold regions become thoroughly mixed. Related to this, an estimation of when the system has reached a pseudo-steady-state is of high importance since analysis of data from the time before the pseudo-steady-state would overestimate the heat flux and thereby also the solids dispersion derived from it.

3.4 Method 4—Getting a solids dispersion from the velocity fluctuations (DispFluct)

A fourth method to assess solids dispersion was presented by Jiradilok, Gidaspow and Breault [13], who proposed to make an estimate of the solids radial and axial dispersion values based on the flow alone without introducing any analogy to heat or mass transfer. They therefore compute the turbulent dispersion (the dispersion due to oscillations of particle clusters in the bed) based on the fluctuating solids velocity. The procedure in which the turbulent dispersion coefficient, D_T, is calculated from the sampled velocity fluctuations is explained in the ESM.

This method yields local dispersion values at each point in the domain, and to arrive at a value that is comparable to the previous global assessments, these local values must be averaged. This averaging is performed along a plane in the bed, as suggested by Jiradilok et al. [13]. Two planes were evaluated, one at y = 0.3525 m and one at y = 0.1125 m. It should be noted that the auto-correlation is a time-averaged quantity, thus once again sufficient time-resolved data needs to be obtained to capture the relevant phenomena in the bed. Note also that since the auto-correlation is evaluated in terms of a specific lag τ, the time series available to form the average with should in fact be shortened by as much as the maximum τ evaluated, leading to even longer simulation times required.

An advantage with method 4, in relation to the other methods, is that it does not require fitting of a dispersion value by comparison with an ideal reference problem. For this reason, method 4 can be applied to quantify the mixing in all directions simultaneously. Another point worth noting is that it is not obvious that the values returned by this method should match those of the other three methods. The dispersion coefficient returned by this method is based on the deviations from the mean flow in each cell. While this makes the value truer to the concept of dispersion, it also means that on a larger scale any contribution to the lateral mixing by the mean flow will be neglected. The other three methods have had all lateral transport mechanisms lumped together, i.e., both the “true dispersive” motion and the convection by the mean flow. Thus, the values outputted by this method are not meant for comparison with those from the other methods used in this work, but rather as data offering a deeper analysis of the solids flow. The interested reader is directed to Askarishahi et al. [36] for a comprehensive discussion on the differences between solids dispersion and solids diffusion.

3.5 Experiments

Experimental data on the lateral solids dispersion for the bed simulated in the HiVoid cases are available in Sette et al. [33]. These experiments used an indirect tracking method for bed material using magnetic separation and are thus based on a mass-transfer approach. For the LoVoid cases, new experiments were carried out to determine the lateral dispersion of solids in fluidized beds under conditions relevant for industrial operation using a heat-transfer approach. The parametric study included the effect of fluidization velocity, solids particle size and pressure drop over the distributor. The results from these experiments are used in the present work for model validation. A description of the experimental setup and methodology is provided in the ESM; for a complete account of the study the reader should refer to Martinez Castilla et al. [37].

4 Results and discussion

4.1 General characterization of the simulated flow

To provide some visual understanding of the main bed characteristics, instantaneous voidage contours of the simulated beds are presented first in Fig.2. The smaller particles and much higher fluidization velocity in the HiVoid cases (Tab.1) yield the bed characteristics as illustrated in Fig.2(a), exhibiting a more aerated bed volume than for larger particles and lower fluidization velocity presented in Fig.2(b) and exemplifying LoVoid cases. Also worth noting is the much more structured character of the bed used in the LoVoid-small-nozzle, which is depicted in Fig.2(c). Here, a nozzle boundary condition is employed and, as can be expected for a shallow bed like this, the nozzles seem to lock the bubble paths into a self-reinforcing pattern [38].

Fig.2 A visual comparison of instantaneous volume fraction fields in the different beds studied. (a) The HiVoid case with smaller sand corresponds to the solids used in the experiments of Sette et al. [33], and (b, c) the LoVoid cases with larger sand corresponds to the solids used in the experiments of Martinez Castilla et al. [37]. The LoVoid-small-nozzle case (c) employs “nozzle” inlet boundary condition instead of the uniform “porous plate” inlet boundary condition used in the other cases (a, b).

Full size|PPT slide

4.2 Comparison with experiments

Fig.3 shows the collected results for the HiVoid cases (higher fluidization velocity, higher voidage). The figure shows data for methods 2 and 3, but not for methods 1 and 4 since these proved to be unreliable for reasons discussed in the followings. The deviation between the solids dispersion coefficient values obtained from simulations and the experimental one is between 5% and 56%, where the cases with a larger simulation domain yield a better agreement to experiments. Hence, both methods 2 and 3 can give a reasonable estimate given a large enough domain. Method 3 seems less sensitive to domain effects compared to method 2.

Fig.3 Experimental [33] and simulated (current work) values of the solids lateral dispersion coefficient for the HiVoid cases, extracted using methods 2 (DispTrans) and 3 (CondSteady).

Full size|PPT slide

The corresponding comparison for the LoVoid (lower fluidization velocity) cases is shown in Fig.4. Worth noting is that, regardless of the domain size, all the cases which employed a uniform (porous plate) inlet boundary condition yield a higher dispersion coefficient than the experimentally-derived value (by a factor between roughly 3 and 4). In contrast, the simulations using nozzles yield dispersion coefficient values almost fully within the experimental uncertainty range. This observation is consistent with the fact that these experimental results were obtained with a low-pressure-drop gas distributor.

Fig.4 Experimental [37] and simulated (current work) values for the denser bed (LoVoid cases). The experimental error bars indicate lower and upper bounds for the method used (more information is provided in the ESM).

Full size|PPT slide

4.3 Evaluation of methods

To assess the different methods of estimating D_L, the complete set of results is summarized in Fig.5. The most striking observation in this comparison is that method 1 is consistently lower than the other methods—for the lighter bed even drastically so—suggesting that there is some systematic difference at play. Another interesting observation is that the spread between the methods seem to decrease when using a larger simulation domain. Further, the uniform (porous plate) boundary condition cases yield higher values than the case with nozzles.

Fig.5 A collected comparison of the different analysis methods tried for all the available data sets. The solid black lines indicate the experimentally obtained DL for the HiVoid case and the experimentally obtained uncertainty interval for DL for the LoVoid cases (Fig. 4).

Full size|PPT slide

Owing to the long simulation times and large data storage requirements of method 4 (DispFluct), only one case was analyzed using this method. From the current data, this prediction seems to follow the overall trend albeit being lower compared to the other methods. Further investigation and application to different fluidization conditions would be beneficial.

4.3.1 Method 1

The consistently lower values obtained from method 1 can be explained by further examining the simulated temperature profiles. As discussed in Section 2.2, the value of k_eff obtained is mainly characterizing the solids lateral motion, and thus constitutes convective fluxes. However, at any given time, the net flux across the boundary between a fixed-temperature region and the inner bed (Fig.1(b)) must be one-directional. Hence, the average temperature of the boundary felt by the inner region of the bed will always differ somewhat from the desired fixed boundary temperature, since only incoming particle fluxes (to the inner bed) contribute to maintaining this temperature, whereas outgoing particle fluxes do not implement any external boundary condition. To this background, method 1 will always tend to underestimate k_eff, and hence also D_L. A more elaborate analysis of the temperature profiles attained is presented in Fig. S3 (cf. ESM).

4.3.2 Method 2

The boundary condition for this method is that no mass is transferred through the wall, a simulation condition which accurately resembles the experimental reference system. The ESM contains a comparison of the analytical to the simulated concentration profiles over time which exemplifies the good agreement between the theoretical setup and experimental data (Fig. S4, cf. ESM). Note that the simulation time needed with this method is relatively short since the valuable information is gathered before the system reaches steady-state, rather than after.

Overall method 2 would seem to be the most promising method, owing to well-described boundary conditions, simplicity of numerical implementation, reduced simulation time, and, nonetheless, the good agreement with experimental values.

4.3.3 Method 3

While the results from method 3 seem to follow those of method 2 in most cases, there are some points to be made about its reliability still. While the conduction value can be deduced through a simple algebraic relation, it is based on averaging the time-varying heat flux through the domain. As with the formation of any average, this method requires sufficient data to provide a stable value. An illustration of the slow convergence of the rolling mean of the heat flux is provided in Fig. S5 (cf. ESM). The straightforward interpretation of method 3 is a desirable trait and, based on the results in Fig.3 and Fig.4, it does seem to approximate the experiments reasonably well given a suitably large domain. The biggest downside of this method is the long simulation times needed to reduce the uncertainty of the average heat flux value to desirable levels.

4.3.4 Method 4

For this method, as with method 3, the same kind of complications with respect to the large data requirement for statistical robustness of the time-average value formation arises. Fig.6 shows the calculated time-averaged turbulent dispersion coefficients in the x, y, and z directions, over two planes at different elevations in the dense bed. It should be noted that none of the planes was situated in the splash zone of the bed, where solids ejected by bubble eruptions might possibly increase the lateral dispersion [34].

Fig.6 The computed turbulent dispersion DT in each direction (averages along two planes at different heights in the bed).

Full size|PPT slide

As seen previously from Fig.5, the calculated value for the turbulent dispersion is lower than the lateral dispersion values from methods 2 and 3. This is in line with the fact that the turbulent dispersion represents only one part of the lumped dispersion coefficient assessed by the other methods. Thus, the value obtained with method 4 tells about the breakdown of the lumped mixing into convective and true dispersive mixing.

Fig.6 shows the expected result that the dispersion is about one order of magnitude higher in the axial direction (y) than in the lateral directions (x, z)—which exhibit similar dispersion coefficients although with some deviation likely related to the asymmetry of the simulation domain.

4.4 Requirements on the computational domain

The mixing in a bubbling fluidized bed is governed by the convective motion of counter-current vortices around the bubble paths in the bed [38]. The smallest tolerable domain size is therefore dictated by the inherent length scales of the specific bed to be simulated. To ensure that the mixing is not influenced by the domain size, the domain needs to be large enough not to limit the sizes of the aforementioned vortices. In studies where the influence of bed width has been investigated, the lateral dispersion has been observed to increase with bed width until it eventually reached a stationary value [14,32].

Moreover, care should be taken to model the bottom inlet geometry well, since in a shallow bed the bubble paths tend to be fixed above the fluidization nozzles [39,40], thus establishing additional restrictions on the size of the vortices. In other words, the minimum domain size might well be determined mostly by the nozzle spacing for the type of beds investigated here.

It is also worth to highlight that the different methods discussed here could be helped or even hampered to different amounts by extension of the domain side lengths, L_x and L_z. An extension of L_x means for the fitting methods (1 and 2) that minor kinks in the curves due to velocity fluctuations will have less impact on the results, which is beneficial. However, for method 3 the effects are more mixed. A longer mixing path means that the profiles will evolve toward pseudo-steady-state under a longer period of time, which means that the time until the data collection can start for method 3 increases. On the other hand, a longer mixing path would make the material entering the hot or cold zones more homogeneous, possibly lessening the fluctuations somewhat. However, a further downside is the fact that for method 3 the mixing length ∆x is directly included in the algebraic estimation of k_eff, acting as a multiplier on any fluctuations. The effects on method 4 should be minimal, since no reference system is introduced and all that is needed is realistic flow patterns which should be provided already by the smallest tolerable bed size as discussed above.

The consequences of an extension of L_z are clearer. Aside from method 4 (for the reasons just discussed), all methods should benefit from an extended z-dimension. The reason for this is simply that, when extending L_z, more material is being mixed, without taking longer (physical) time to mix. For the fitting methods 1 and 2 this is akin to running multiple parallel mixing experiments and ensemble averaging the results. Such a procedure would lead to smoother curves and more true profiles since a larger amount of different flow interactions will take place within the simulated experiment. For method 3, more material being mixed should likely make the average mass flow across the hot and cold planes approximated better. The only drawback is the increase in computational cost.

4.5 Simulation time

While the simulation time needed is not something that needs to be decided beforehand to set up the simulation, it is still good to be able to assess how long a given simulation will have to run; especially since 3D multiphase simulations like these can tend to take prohibitively long time. However, just as with the domain requirements, the time needed is different for each method examined. To help explain this more clearly a diagram is provided in ESM. The main conclusion is that methods 1 and 2, which only require part of a transient for the fitting procedure to converge, are significantly faster than methods 3 and 4 that depend on time-averaged quantities. In particular, method 3 requires that the bed is in a statistically steady-state before the data collection can even begin.

4.6 General discussion

While there are many different factors that influence how accurate a numerical prediction of the lateral bed dispersion will be, most factors can be related to one of the following three categories.

The fidelity of the underlying TFM simulation. For example, Oke et al. [14] evaluated the dispersion using submodels with different onset of frictional stress, and in their analysis an earlier onset of frictional stress led to smaller bubble sizes in the bed and through that also reduced lateral dispersion. Farzaneh et al. [41] evaluated three different models for the frictional stress and found that, based on what frictional stress model that was used, the general flow pattern changed significantly, and the range and localization of the local horizontal dispersion coefficient changed as well. Furthermore, it is well known that the drag law used has a significant influence and could well represent one of the most important factors in an accurate TFM implementation for dense fluidized bed simulations [42]. In conclusion, the influence of submodels and inlet geometries on lateral dispersion is something that needs further attention, especially in 3D.

The choice of analysis method. As the results in Fig.5 show, even seemingly similar methods can give different results. If a method based on fitting simulation results to an analytical solution is to be used, care has to be taken to ensure that the simulated setup really can be compared to the analytical one. Method 1 provides a clear example of a case where not all parts of the simulated problem were compatible with the analytical problem. The choice of method also influences the choice of data collection since different methods utilize different parts of the generated data. Lastly, the various methods differ in their inherent assumptions, and depending on the intended use for the estimated dispersion coefficient and the available experimental data, a method which is more costly but has fewer or less severe inherent assumptions could still be preferred over a cheaper one.

The approach to data collection (i.e., what part of the domain that is recorded and with what frequency). While this aspect is mostly limited by practical issues like disk size, post-processing time and degree of automation, it should still be considered. Sufficient data to make unbiased spatial averages must be exported at a rate that will resolve the diffusion process well. However, what is enough will change on a case-to-case basis and vary for different methods. As an example, for the LoVoid-large-uniform case, the PSD showed that a sampling frequency of 20 Hz would capture most of the observed variation. The fitting of effective conduction and dispersion coefficients in methods 1 and 2, on the other hand, is rather robust and works well with even an order of magnitude lower frequency in the stored data, as not more is needed to resolve the transition to the steady-state in the corresponding analytical problem well enough. Moreover, for method 3, complete reconstruction of the signal is not necessary since the quantity of interest is the mean. In other words, given that the sampling frequency is somewhat higher than the lower frequencies in the bed and that the total sampling time at pseudo-steady-state is long enough, capturing the smaller fluctuations is not necessary to obtain the same average value over time. Method 4, on the other hand, relies on characterizing local velocity fluctuations, indicating that a higher frequency is motivated as also the smaller fluctuations contribute. For methods 1 and 2, the duration of the transient is also affected by the length of the domain, meaning that a larger domain could have a lower sampling frequency and still be well resolved. In further studies, a more thorough statistical investigation could be made to quantify the effects of data collection more in depth.

5 Conclusions

Four different methods have been assessed for the numerical evaluation of the lateral dispersion coefficient in a bubbling fluidized bed. Methods are based on: fitting the transient response of an ideal heat or mass transfer problem to an analytical solution, time-averaging a heat transfer signal, and using information from the velocity field in an analogy with turbulent dispersion.

Based on comparisons with experiments, we conclude that the TFM can indeed be used to provide a quantification of the solids lateral mixing in shallow fluidized beds. However, care must be taken when selecting a method for postprocessing the simulated data, as the investigated methods may differ substantially in their predictions. Given the results from this study, fitting a simulated concentration profile to an analytical one has the greatest potential. Among the other methods, only the calculation of a dispersion coefficient through evaluation of the time-averaged heat flux was able to provide a reasonable estimation albeit at the additional cost of needing much longer simulation times.

Simulations showed that the domain size had a significant impact on the simulated solids dispersion values. The current work also highlights the significant improvements in the agreement to experimental data when using a boundary condition mimicking the effect of nozzles as compared to the conventional porous plate (uniform) boundary condition. This observation should be all the more relevant for industrial conditions. Further work could investigate the influence of using different submodels in the TFM formulation, as well as extend the validation to beds fluidized at different conditions.

参考文献

[1] YangW C. Handbook of Fluidization and Fluid-Particle Systems. Taylor & Francis, 2003

[2] Kraussler M , Binder M , Schindler P , Hofbauer H . Hydrogen production within a polygeneration concept based on dual fluidized bed biomass steam gasification. Biomass and Bioenergy, 2018, 111: 320–329

[3] Heinze C , May J , Peters J , Ströhle J , Epple B . Techno-economic assessment of polygeneration based on fluidized bed gasification. Fuel, 2019, 250: 285–291

[4] Eschenbacher A , Jensen P A , Henriksen U B , Ahrenfeldt J , Jensen C D , Li C , Enemark-Rasmussen K , Duus J Ø , Mentzel U V , Jensen A D . Catalytic upgrading of tars generated in a 100 kWth low temperature circulating fluidized bed gasifier for production of liquid bio-fuels in a polygeneration scheme. Energy Conversion and Management, 2020, 207: 112538

[5] Salman C A , Naqvi M , Thorin E , Yan J . Impact of retrofitting existing combined heat and power plant with polygeneration of biomethane: a comparative techno-economic analysis of integrating different gasifiers. Energy Conversion and Management, 2017, 152: 250–265

[6] Salman C A , Omer C B . Process modelling and simulation of waste gasification-based flexible polygeneration facilities for power, heat and biofuels production. Energies, 2020, 13(16): 4264

[7] Gadsbøll R Ø , Clausen L R , Thomsen T P , Ahrenfeldt J , Henriksen U B . Flexible TwoStage biomass gasifier designs for polygeneration operation. Energy, 2019, 166: 939–950

[8] Yu P , Luo Z , Wang Q , Fang M . Life cycle assessment of transformation from a sub-critical power plant into a polygeneration plant. Energy Conversion and Management, 2019, 198: 111801

[9] Naqvi M , Dahlquist E , Yan J , Naqvi S R , Nizami A S , Salman C A , Danish M , Farooq U , Rehan M , Khan Z . . Polygeneration system integrated with small non-wood pulp mills for substitute natural gas production. Applied Energy, 2018, 224: 636–646

[10] Zhao K , Thunman H , Pallarès D , Ström H . Control of the solids retention time by multi-staging a fluidized bed reactor. Fuel Processing Technology, 2017, 167: 171–182

[11] Salehi M S , Askarishahi M , Radl S . Quantification of solid mixing in bubbling fluidized beds via two-fluid model simulations. Industrial & Engineering Chemistry Research, 2020, 59(22): 10606–10621

[12] Diba M F , Karim M R , Naser J . Numerical modelling of a bubbling fluidized bed combustion: a simplified approach. Fuel, 2020, 277: 118170

[13] Jiradilok V , Gidaspow D , Breault R W . Computation of gas and solid dispersion coefficients in turbulent risers and bubbling beds. Chemical Engineering Science, 2007, 62(13): 3397–3409

[14] Oke O , Lettieri P , Salatino P , Solimene R , Mazzei L . Numerical simulations of lateral solid mixing in gas-fluidized beds. Chemical Engineering Science, 2014, 120: 117–129

[15] Hernández-Jiménez H , Sánchez-Prieto J , Cano-Pleite E , Soria-Verdugo A . Lateral solids meso-mixing in pseudo-2D fluidized beds by means of TFM simulations. Powder Technology, 2018, 334: 183–191

[16] Yu M , Miller D C , Biegler L T . Dynamic reduced order models for simulating bubbling fluidized bed adsorbers. Industrial & Engineering Chemistry Research, 2015, 54(27): 6959–6974

[17] Wang H , Li Z , Li Y , Cai N . Reduced-order model for CaO carbonation kinetics measured using micro-fluidized bed thermogravimetric analysis. Chemical Engineering Science, 2021, 229: 116039

[18] Yuan T , Cizmas P G , O’Brien T . A reduced-order model for a bubbling fluidized bed based on proper orthogonal decomposition. Computers & Chemical Engineering, 2005, 30(2): 243–259

[19] Li C , Dai Z , Sun Z , Wang F . Modeling of an opposed multiburner gasifier with a reduced-order model. Industrial & Engineering Chemistry Research, 2013, 52(16): 5825–5834

[20] Saastamoinen J J . Simplified model for calculation of devolatilization in fluidized beds. Fuel, 2006, 85(17-18): 2388–2395

[21] Kaushal P , Abedi J . A simplified model for biomass pyrolysis in a fluidized bed reactor. Journal of Industrial and Engineering Chemistry, 2010, 16(5): 748–755

[22] Gómez-Barea A , Leckner B . Modeling of biomass gasification in fluidized bed. Progress in Energy and Combustion Science, 2010, 36(4): 444–509

[23] Kaushal P , Tyagi R . Advanced simulation of biomass gasification in a fluidized bed reactor using ASPEN PLUS. Renewable Energy, 2017, 101: 629–636

[24] Nikku M , Bordbar H , Myöhänen K , Hyppänen T . Effects of heterogeneous flow on carbon conversion in gas-solid circulating fluidized beds. Fuel, 2020, 280: 118623

[25] Sternéus J , Johnsson F , Leckner B . Characteristics of gas mixing in a circulating fluidised bed. Powder Technology, 2002, 126(1): 28–41

[26] Chirone R , Miccio F , Scala F . On the relevance of axial and transversal fuel segregation during the FB combustion of a biomass. Energy & Fuels, 2004, 18(4): 1108–1117

[27] OkeOLettieriPSalatinoPMazzeiL. CFD simulations of lateral solid mixing in fluidized beds. CFB-11: Proceedings of the 11th International Conference on Fluidized Bed Technology, 2014, 287–292

[28] Oke O , Lettieri P , Salatino P , Solimene R , Mazzei L . Eulerian modeling of lateral solid mixing in gas-fluidized suspensions. Procedia Engineering, 2015, 102: 1491–1499

[29] Vandewalle L A , Francia V , Van Geem K M , Marin G B , Coppens M O . Solids lateral mixing and compartmentalization in dynamically structured gas-solid fluidized beds. Chemical Engineering Journal, 2022, 430: 133063

[30] Luo G , Cheng L , Ma Z , Li L , Li Z , Wang P , Li L , Rong H . MP-PIC simulation on solid dispersion in a 350 MW CFB boiler. Industrial & Engineering Chemistry Research, 2022, 61(45): 16857–16868

[31] Luo G , Cheng L , Kang Q , Zhang Q , Li K , Guo Q , Li W . Investigation on solid dispersion in a CFB dense region with bluetooth tracking method and MP-PIC simulation. Powder Technology, 2023, 429: 118892

[32] Liu D , Chen X . Lateral solids dispersion coefficient in large-scale fluidized beds. Combustion and Flame, 2010, 157(11): 2116–2124

[33] Sette E , Pallarès D , Johnsson F . Experimental evaluation of lateral mixing of bulk solids in a fluid-dynamically down-scaled bubbling fluidized bed. Powder Technology, 2014, 263: 74–80

[34] Oke O , Van Wachem B , Mazzei L . Lateral solid mixing in gas-fluidized beds: CFD and DEM studies. Chemical Engineering Research & Design, 2016, 114: 148–161

[35] BakshiA. Multiscale continuum simulations of fluidization: bubbles, mixing dynamics and reactor scaling. Dissertation of PhD degree. Boston: Massachusetts Institute of Technology, 2017

[36] Askarishahi M , Salehi M S , Molaei Dehkordi A . Numerical investigation on the solid flow pattern in bubbling gas-solid fluidized beds: effects of particle size and time averaging. Powder Technology, 2014, 264: 466–476

[37] Martinez Castilla G , Larsson A , Lundberg L , Johnsson F , Pallarès D . A novel experimental method for determining lateral mixing of solids in fluidized beds—quantification of the splash-zone contribution. Powder Technology, 2020, 370: 96–103

[38] PallarèsDDíezPJohnssonF. Experimental analysis of fuel mixing patterns in a fluidized bed. The 12th International Conference on Fluidization-New Horizons in Fluidization Engineering, 2007

[39] Norouzi H R , Mostoufi N , Mansourpour Z , Sotudeh-Gharebagh R , Chaouki J . Characterization of solids mixing patterns in bubbling fluidized beds. Chemical Engineering Research & Design, 2011, 89(6): 817–826

[40] Olsson J , Pallarès D , Johnsson F . Lateral fuel dispersion in a large-scale bubbling fluidized bed. Chemical Engineering Science, 2012, 74: 148–159

[41] Farzaneh M , Almstedt A E , Johnsson F , Pallarès D , Sasic S . The crucial role of frictional stress models for simulation of bubbling fluidized beds. Powder Technology, 2015, 270: 68–82

[42] Ostermeier P , DeYoung S , Vandersickel A , Gleis S , Spliethoff H . Comprehensive investigation and comparison of TFM, DenseDPM and CFD-DEM for dense fluidized beds. Chemical Engineering Science, 2019, 196: 291–309

Competing interests

The authors declare that they have no competing interests.

Acknowledgements

This work was financially supported by the Swedish Energy Agency through the Swedish Centre for Biomass Gasification (SFC, project number P34721-3). The computations were enabled by resources provided by the Swedish National Infrastructure for Computing (SNIC) at Chalmers Centre for Computational Science and Engineering (C3SE) partially funded by the Swedish Research Council through grant agreement No. 2018-05973.

Electronic Supplementary Material

Supplementary material is available in the online version of this article at https://doi.org/10.1007/s11705-024-2503-4 and is accessible for authorized users.

Funding Note

Open access funding provided by Chalmers University of Technology.

Open Access

This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit https://creativecommons.org/licenses/by/4.0/.

版权

2024 The Author(s) 2024. This article is published with open access at link.springer.com and journal.hep.com.cn