Computational Fluid Dynamics-Discrete Element Method (CFD-DEM) Study of Mass-Transfer Mechanisms in Riser Flow.

We report a computational fluid dynamics-discrete element method (CFD-DEM) simulation study on the interplay between mass transfer and a heterogeneous catalyzed chemical reaction in cocurrent gas-particle flows as encountered in risers. Slip velocity, axial gas dispersion, gas bypassing, and particle mixing phenomena have been evaluated under riser flow conditions to study the complex system behavior in detail. The most important factors are found to be directly related to particle cluster formation. Low air-to-solids flux ratios lead to more heterogeneous systems, where the cluster formation is more pronounced and mass transfer more influenced. Falling clusters can be partially circumvented by the gas phase, which therefore does not fully interact with the cluster particles, leading to poor gas-solid contact efficiencies. Cluster gas-solid contact efficiencies are quantified at several gas superficial velocities, reaction rates, and dilution factors in order to gain more insight regarding the influence of clustering phenomena on the performance of riser reactors.

Introduction

Mass- and heat-transfer phenomena under riser flow conditions have been widely investigated during the last decades. The usage of mass- and heat-transfer coefficients is essential in phenomenological and computational models to estimate the performance of chemical processes in bubbling and fast fluidized bed reactors. However, there is still a lack of understanding regarding the applicability of the different correlations that could be employed to estimate mass- and heat-transfer coefficients in fluidized systems. Many mass-transfer correlations for fluidized systems have been proposed to compute the dimensionless mass-transfer coefficient of a particle in a fluidized system.[1,2] Furthermore, Breault reported that Sherwood numbers in riser fluidized systems can differ in several orders of magnitude,[3] so it is evident that there is not a unique and unequivocal equation to compute a mass-transfer coefficient for fluidized systems. The reason for this is thought to be the presence of heterogeneities in the particle distributions that severely influence the flow patterns and consequently the mass transfer to particles.

Flow heterogeneities are especially prevalent in riser flows, which are characterized by a core annulus flow involving a rather dilute region in the core of the riser and a dense particulate phase close to the walls. The dense solids phase can be either formed by a falling solids film (annulus) or can be a more cluster-like flow, which has an intermittent behavior that can be represented by, for example, an intermittency index.[4] When particle clusters are formed, low slip velocities, gas back-mixing, and poor gas–solid contacting are key hydrodynamic phenomena that prevail. Besides heterogeneities, another cause of disagreement between mass-transfer coefficients reported in the literature is caused by different model interpretations and/or definitions of mass-transfer coefficients. Thus, it should be noted that the obtained values of mass-transfer coefficients can be completely meaningful within the context of their respective models, but these data should in principle not be used in other scenarios where (completely) different hypotheses and/or modeling assumptions have been made.

In this paper, we generate insight about the mass transport mechanisms related to cluster formation and the causes of the aforementioned disagreement. For this purpose we apply computational fluid dynamics–discrete element method (CFD-DEM) simulations. CFD-DEM simulations can not only estimate the global process efficiency of a riser reactor but also provide detailed information about cluster-related phenomena. The reason is that in CFD-DEM all particles are explicitly modeled. Therefore, phenomena such as clustering of particles are emergent. This means that closures for mass transfer are needed at the particle level only. Extra mass-transfer resistances related to the presence of particle clusters are a result of the simulations. Thus, by means of CFD-DEM simulations, a detailed analysis is performed regarding the influence of several mass-transfer mechanisms on the performance of a pseudo-two-dimensional (2D) riser reactor.

Frössling-type correlations have been reported to predict accurately the mass-transfer rate of a single particle immersed in a fluid flow.[1,5,6] Convection has a strong impact on the thickness of the momentum and mass boundary layers and for higher Reynolds numbers determines the mass-transfer resistance at the particle level.[7] Because of the thin boundary layer, the particle Sherwood number of an isolated particle at larger Reynolds number will be greater than 2. When particles are surrounded by other particles, such as in packed bed reactors, the mass transfer is influenced significantly compared to isolated particle mass tranfer.[8−10] Gunn’s correlation, which is widely used in the literature,[11] was developed from experiments in a packed bed reactor and a liquid fluidized bed reactor. The Gunn correlation has also been shown to represent well the fluid to particle mass/heat transfer in resolved simulations of dense[12] and more dilute stationary particle arrangements[13,14] as well as solid–liquid systems.[15] Note that the dependence on solids-volume fraction is different if a particle that exchanges mass is surrounded by inert particles that do not exchange mass.[1,16−19] Such “diluted” systems have been employed to not reach saturated gas naphthalene concentrations[1] or to not reach too high conversion rates. Therefore, the presence of other particles significantly influences the mass transfer to individual particles. It has to be noted that, for example, Gunn’s correlation was obtained from experimental measurements in packed bed and liquid fluidized bed systems. However, in riser flows, particles may form heterogeneous structures and the gas–solid distributions may not be as homogeneous as in a fixed bed reactor or a bubbling liquid fluidized bed. It should be emphasized that Gunn’s correlation is valid if we assume that gas and solids distributions are relatively homogeneous. Thus, in our simulations, Gunn’s correlation is employed at a sufficiently small scale where the flow can be regarded as homogeneous, e.g. the CFD-DEM cell size.

For heterogeneous systems, the distinction between the local particle mass-transfer coefficient and the global or overall (i.e., reactor length scale) mass-transfer coefficient should be made. A particle-based mass-transfer coefficient is interpreted as the mass-transfer rate of a single particle that exchanges mass with the surrounding fluid phase, where the driving force is defined by local concentration differences. At the same time, an overall mass-transfer coefficient can well represent the mass-transfer rate at the reactor scale. Here the driving force for mass transfer is globally defined, for example, by means of a gas concentration cup-averaged over the reactor cross section. However, different hydrodynamic models can be found in the literature, and consequently, definitions of global mass-transfer coefficients could also differ. Thus, the measurements of global mass-transfer coefficients in fluidized systems (e.g., assuming a 1D plug flow interpretation model) can be clearly different from Gunn’s correlation predictions (using a global driving force) because the flow is heterogeneous and local particle-scale driving forces are different from the global reactor-scale driving forces.[20,21] Larger (than particle level) scale hydrodynamic resistances can play an important role in mass-transfer phenomena.[22−25]

Breault et al. reported that Sherwood numbers in riser fluidized systems can differ in several orders of magnitude.[3] This spread in reported values is likely due to the presence of different levels of heterogeneity in different experiments. However, another important cause of these differences is that different models and definitions have been employed to report such Sherwood numbers. For instance, Subbarao and Gambhir[26] utilized naphthalene deposition on sand particles at the bottom of a riser as a model experimental system to develop a correlation to compute an overall mass-transfer coefficient. In their work, circulation pattern effects are lumped into global parameters of gas superficial velocity and solids mass flux. Wang and Li[27] used a two-fluid model coupled to a energy-minimization multiscale (EMMS) method to model the local heterogeneous flow structure. The model was validated, using a particle-based Sherwood correlation,[28] with experimental measurements of naphthalene sublimation.[27] Venderbosch et al.[16] employed CO oxidation as a model reaction to determine mass-transfer-limited rates of a diluted (chemically active and inert particles) fluidized system, where platinum-supported FCC particles (catalytically active) were mixed with unsupported FCC particles (inert). Sherwood numbers reported by the previously mentioned authors can easily differ in several orders of magnitude, and all might be correct for their respective models. In Table S4 (in the Supporting Information), an overview of relevant publications related to mass-transfer phenomena on riser flow is provided. It can be noticed that CFD computational strategies, which account for the flow heterogeneities, have become more popular in recent years. However, the computational expenses of these models are still too high to apply them at large industrial scale.

Thus, some attempts have been made to develop correlations that compute overall mass-transfer coefficients that could be easily estimated under riser flow conditions with known operating parameters such as the gas superficial velocity and solids mass flux of a riser.[23,26,29] Other authors suggested that the cluster–bulk mass-transfer coefficient is proportional to the cluster surface renewal rate and that this is defined by their respective size and velocity.[30−33] Many core–annular models have been developed to predict the solids and gas distributions under riser flow conditions[34−37] and to estimate the performance of a riser reactor. However, these are mainly used to model only fully developed sections of a riser and do not model complex flow structures, which are highly dependent on the physical properties of the particles and the reactor geometry. It is known that riser hydrodynamics can be highly affected by the riser diameter, entrance, and exit effects; therefor,e more detailed models are required. Thus, high-precision CFD models are expected to describe riser flow hydrodynamics more accurately and provide a complete picture of the system to facilitate design and upgrading of riser reactors.

EMMS-based computational works[38,39] and filtered closures[40] have been employed to characterize flow heterogeneities. In this way, heterogeneity is modeled in Euler–Euler simulations, in which the resolution scale can be significantly larger than the typical cluster size. For instance, EMMS schemes are based on energy minimization assumptions to model flow heterogeneous structures. Analogous EMMS–mass-transfer methods can be employed to determine interphase mass-transfer coefficients: from gas to particles of a dilute phase, from gas to particles of the dense phase, and from gas to cluster surface.[20,21] Other authors employed EMMS–mass-transfer models to analyze the influence of particle heterogeneity over the riser performance.[41−43] Although some of the EMMS-based computational studies account for clustering phenomena, simplifying assumptions are made. Some of these are (1) clusters have a uniform density, equal to maximum packing fraction[44] and (2) clusters are spherical structures that interact with a lean surrounding phase and are homogeneously dispersed in a control volume.[45] Some of these computational studies have suggested that global scale mass-transfer resistances are related to cluster formation.[20,25,46] However, these findings are not explicitly supported with cluster scale level computations, such as the instantaneous gas–solid contact efficiency of particle clusters.

In riser flows, the definition of a gas–solid contact efficiency has been employed to quantify the deviation of a riser flow from an idealized plug flow, where all catalyst particles are fully exposed to the bulk fluid phase. The contact efficiency is not only related to the exposed area of a single particle to the surrounding gas phase[47] but also related to the particle exposure to a reactant-rich gas phase that tends to circumvent dense particulate regions,[48] which can retain gas pockets of highly depleted reactant. In experimental works,[35,49−52] inefficient contacting was presumed to be due to cluster formation because gas dispersion effects could be discarded.[16,50] In these cases, the quantification of gas–solid contact efficiency was computed from time-averaged properties (solids volume fraction, reactant concentration, and gas velocity) without obtaining instantaneous data that confirmed these hypotheses. Previous authors[25,42] analyzed the influence of the reaction rate on the mass-transfer rate, which is reported to be less affected at lower values of kr.[42] The difficulty of measuring experimental data at the cluster level has prevented the confirmation of these assumptions by a direct method.

CFD-DEM has been shown to be a suitable model to predict riser hydrodynamics and complex clustering phenomena,[53] study particle mixing and deactivation in a catalytic process, or even track real-time particle activity.[54,55] In this paper, different mass-transfer mechanisms are isolated using ozone decomposition as a model reaction. The aim is to quantify their respective contribution to mass-transfer-limited chemical processes and to obtain a close relation between cluster formation and gas–solid contact efficiency. In addition, we quantify the instantaneous gas–solid contact efficiency by relating the ozone gas mass fraction inside the cluster and the bulk ozone gas mass fraction at the cross section where each detected cluster is located. These results illustrate the severe impact that clustering phenomena can have on gas–solid contacting and the estimation of global mass-transfer coefficients. By means of a CFD-DEM, we confirm that clustering phenomena have a more severe impact at higher reaction rates.

Methodology

A CFD-DEM has been utilized to resolve simulations of a pseudo-2D riser bed reactor. The governing equations for the gas phase have been solved using a semi-implicit finite difference technique.[56] For the species field, the gas convection–diffusion equation has been solved as well:The Zehner and Schlünder[57] model is used to evaluate the effective diffusivity:The gas and solids motion is coupled via a sink term that involves the computation of the interphase momentum-transfer coefficient by means of the Beetstra drag correlation.[58] Additionally, gas-phase mass balances are coupled to the mass balance at the particle level by means of a sink term, where the mass-transfer coefficient is computed by means of the Gunn correlation:[8]where Ap and Vp are the surface area and volume of a single particle, respectively; β is the interphase momentum-transfer coefficient, and kmt is the computed particle-based mass-transfer coefficient by Gunn correlation.

Two way coupling is performed by means of a regularized Dirac delta function D( – p), which maps the gas-phase variables from neighboring Cartesian nodes to the particle location to enable evaluation of the drag force and the interphase mass-transfer coefficient. Moreover, changes in particle momentum are fed back from the particle position to the surrounding Eulerian nodes using the same regularized function.[56]

The particle translational and rotational momentum is governed by the Newtonian equations of motion:

The particle collisional forces are deterministically computed, by means of a soft sphere model that was originally proposed by Cundall and Strack[59] and first employed in a gas fluidized system simulation by Tsuji et al.[60] Gas turbulence was assumed to be insignificant compared to the velocity fluctuations due to gas–particle interaction, following other authors’ observations.[61,62] Thus, a subgrid turbulence model was not employed. Additionally, no mass-transfer exchange was taken into account during particle collisions.

A first-order irreversible reaction under riser flow conditions was included in the CFD-DEM to simulate ozone decomposition:Pure ozone was fed from the bottom of the riser. It was assumed that ozone reacts at the particle surface without taking into consideration internal mass-transfer effects throughout the particle volume:

After integration, the particle ozone mass fraction can be expressed as follows:If we compute the characteristic time scale of the exponential term, assuming a 0.85 mm glass bead with a slip velocity of 1 m/s and kr = 10 s–1, we obtainwhich is lower than the computational gas time step of our simulations (5 × 10–5 s). It is noted that the computed characteristic time scale is an upper bound (when kr = 10 s–1) in our set of simulations. Therefore, for all other values of kr considered in this work, this parameter is even lower. Because this reaction is fast compared to the gas residence time in the system, the mass balance at the particle level is well approximated bywhere the Sherwood number was computed by means of the Gunn correlation.[8] As for momentum, the ozone mass fraction is computed for all particles in the system. The gas concentration wA,g is interpolated from the neighboring Cartesian nodes to the particles locations, which are represented as points and act as ozone sink sources. As a result, the particle mass fraction wA,p is computed by means of eq , and the gas mass fraction is updated via two-way coupling.[56]

Thus, the source term (eq ) can be expressed as follows:

Simulation Conditions

Figure shows the simulation domain used. This consisted of a pseudo-2D riser of 1.57 × 0.07 × 0.006 m3. Particles were fixed at the top of the simulation domain to mimic the lateral curved outlet of an experimental unit published elsewhere.[53] In the simulations, a system with a highly diluted ozone content is simulated. Thus, physical properties of air gas are used in the simulations, and corrections for molecular counter diffusion of reaction products are not needed. In addition, the reaction is assumed to be equimolar, and no heat effects are considered. More details about the simulation settings and physical parameters used are listed in Tables and 2. During the whole simulation, Geldart D particles were placed at random x–y coordinates near the bottom of the riser at a height of 2–3 times the particle radius above z = 0. The insertion was accepted only if there was no particle overlap. The simulations were performed under fast fluidization regime at several gas superficial velocities (U = 5.55, 5.95, 6.35, and 6.74 m/s) and a fixed solids mass flux rate of 32 kg/(m2s). The particles were inserted at random positions of the bottom X–Y plane of the simulation domain at a velocity that was nearly zero (0.01 m/s). The particles that reached the lateral outlet left the simulation domain.

Figure 1

Simulation domain.

Table 1

Simulation Conditions

NX	28		dp (mm)	0.85
NY	5		ρs (kg/m3)	2500
NZ	628		μg (kg/m·s)	2.0 × 10–5
X (m)	0.07		T (K)	298
Y (m)	0.006		ep-p	0.96
Z (m)	1.57		et	0.33
ΔtGas (s)	5.0 × 10–5		ep-w	0.86
ΔtDEM (s)	5.0 × 10–6		μfr	0.15
Gs (kg/(m2·s))	5.0 × 10–5		kn (N/m)	1600
U (m/s)	5.16–6.74		P	1 atm

Table 2

Characteristic Dimensionless Numbersa

U (m/s)	⟨φs⟩	Re	Pep	St
5.16	0.0567	1858	274	28
5.55	0.0410	1998	295	31
5.95	0.0228	2142	316	33
6.35	0.0147	2286	337	35
6.74	0.0094	2426	358	37

The Damköhler number for all simulations presented in this work ranges between 0.014 (U = 6.74 m/s and kr = 100 s–1) and 13.5 (U = 5.95 m/s and kr = 1000 s–1). These parameters have been computed assuming an effective molecular diffusivity of 1.6 × 10–5 m2/s and a characteristic hydraulic diameter equal to the depth of the system (6 mm).

At the top, front, back, and right walls, no-slip boundary conditions were applied. With a prescribed inflow axial velocity equal to U, gas was supplied at the bottom of the domain. The left side wall (x = 0) was subdivided into two regions: for the top-left outflow region of 0.07 m, the pressure P0 is described and Neumann conditions were applied for the species field. Below this region, no-slip and no-flux boundary conditions were applied for the gas momentum and ozone mass fraction respectively, as Figure illustrates.

Gas–Solid Contact Efficiency

To quantify the gas–solid contacting, it is necessary to provide a parameter that captures this effect, e.g., a contact efficiency. Otherwise, by assuming a homogeneous system, e.g., ideal plug flow, overestimation of the conversion rate in heterogeneous systems is likely to result. The contact efficiency could be determined by assuming a riser as a steady-state plug flow reactor:[24,25,48,50,63]with a solution:where φs is the averaged solids volume fraction in a slice of thickness Δz, U the gas superficial velocity, Kov the apparent volumetric reaction rate constant, γpf a gas–solid contact efficiency, and Da the Damköhler number defined asThe gas–solid contact efficiency γpf could be computed at different axial increments of Δz, by quantifying the contact efficiency as the ratio between the apparent conversion rate and the conversion obtained when a 1D plug flow model is assumed, ignoring any heterogeneity.

In this paper, we perform CFD-DEM simulations to quantify the instantaneous cluster-level gas–solid contact efficiency, which is the ratio between the gas ozone mass fraction inside a cluster region and an average gas mass fraction that will be precisely defined below.

As previously reported,[64] clusters are defined as connected regions with local solids fractions exceeding 0.2 everywhere that have a minimum (projected) area of 60 mm2 and a dense core with at least one grid cell with φs > 0.4. The minimum area requirement limits the amount of noise in our measurements that would be caused by the frequent appearance and disappearance of small clusters. The area of 60 mm2 corresponds to an equivalent circle diameter of 8 mm. The detection of clusters was performed by postprocessing simulation data by means of a Matlab script. In Figure , we show a snapshot of some particle clusters obtained from a typical CFD-DEM simulation. The ozone mass fraction of red-colored cells are averaged in order to compute the ozone mass fraction of a cluster as follows:where wA, is the ozone mass fraction in cell n that is part of the cluster under consideration, φs, the solids volume fraction in that cell, and N the total number of cells that are occupied by that particular cluster.

Figure 2

Snapshot of particle clusters from CFD-DEM.

This average concentration inside the cluster is compared towhere A, is the cross-sectional-averaged ozone mass fraction of slice k, where n is the number of cells of the kth slice that are occupied by the cluster under consideration. The cross-sectional averaged ozone mass fraction, A,, is computed by excluding cells that are identified as part of a cluster. These cross-sectional-averages are thus weighted by the number of cells that the cluster occupies at each cell row. Thus, if a cluster consist of 3 cells at the kth row and 1 cell at the (k – 1)th row, the bulk gas mass fraction (averaged value of those cells that are not occupied by a cluster at that particular row of cells) of cell row k is weighted 3 times over the bulk gas mass fraction of cell row (k – 1), producing a unique value of the bulk gas ozone mass fraction.

The contact efficiency is defined as the ratio between the gas ozone mass fraction inside the cluster, wA,cluster, and the average cross-sectional gas ozone mass fraction of the bulk gas, A:

Efficiency values very close to 1 can then be interpreted as highly efficient contacting, where all catalyst particles are fully exposed to the bulk gas concentrations and gas diffusion through the particle cluster is much faster than the intrinsic reaction rate. Conversely, numbers very close to zero indicate poor gas–solid contacting due to either diffusional limitations or gas bypassing around the clusters. It is worth mentioning that the γpf term can be interpreted as the ratio of external surface area of the catalyst that is exposed to the gas phase, while γcl represents the ratio of external surface area of cluster particles that are exposed to the bulk gas phase. However, these two terms are not comparable because they belong to two different interpretations. The gas–solid contact efficiency, γcl, is a parameter that quantifies the instantaneous gas bypassing around particle clusters and differs from γpf, which measures the deviation from a steady-state 1D ideal plug flow model. Thus, if the dominant mass-transfer resistances are found to be at the particle level, it is expected that γpf ≪ 1, while γcl ≈ 1. The cluster-level contact efficiency can then be employed to identify well the level at which the mass-transfer resistance lies.

Results and Discussion

In this section, we first show some results of mass-transfer coefficients at different operating conditions and values of kr. In this case, a 1D plug flow model is assumed to compute mass-transfer coefficients from time-averaged ozone gas mass fraction profiles using the CFD-DEM generated axial solids distribution. The aim is to show that low Sherwood numbers are obtained when these assumptions are made for riser flows.

In the following subsections, one of our main objectives is to identify and quantify the influence of different mass-transfer mechanisms on the performance of a riser reactor when clusters are present. The influence of reduced slip velocity, axial gas dispersion, and gas bypassing are evaluated. In addition, the influences of the reaction rate and cluster phenomena on the riser performance are quantified and analyzed.

Concentration Profiles: 1D Plug Flow Model

In Figure , time-averaged ozone gas mass fraction profiles at several gas superficial velocities are shown. It can be observed that at the bottom region of the riser, higher conversion rates are found. This is because a dense bottom region exists. It has to be noted that in riser flows there is a trade-off between catalyst holdup and cluster formation. Although high gas superficial velocities can lead to more homogeneous systems (less clustering), both the solids inventory and the gas-phase residence time drop. This gives a lower conversion rate at higher superficial velocities, as we can see in Figure .

Figure 3

Axial profiles of time-averaged ozone gas mass fraction at kr = 100 s–1.

To compute a global mass-transfer coefficient for each one of these cases, a plug flow model can be assumed (see eq ).[24,25,65] The values of Kov were solved through linear regressions above heights of z = 0.2 m, where a constant decaying trend of the ozone gas mass fraction profiles was obtained. The values are provided in Table .

Table 3

Mass-Transfer Coefficientsa

U (m/s)	Kov (s–1)	γpf = Kov/kr	kmt·av (s–1)	kmt (m/s)
5.16	35.48	0.35	56.0	0.0079
5.55	42.69	0.43	74.5	0.0106
5.95	65.09	0.65	186.5	0.0264
6.35	81.18	0.81	431.6	0.0611
6.74	85.75	0.86	601.8	0.0853

U influence at kr = 100 s–1.

When one assumes that there is external mass-transfer resistance only at the particle level, a global resistance analysis can be used to decompose an overall mass-transfer coefficient for each simulation aswhere av is the specific particle surface area av = 6/dp.

The computed overall mass-transfer coefficients are listed in Table for different values of the superficial gas velocity. It can be seen that the overall mass-transfer coefficient increases with increasing gas superficial velocity. It can be noticed that at U values exceeding 5.95 m/s, the order of magnitude of the mass-transfer rates (kmt·av) is similar to that of the reaction rate (kr = 100 s–1). It is clear that the hydrodynamic resistances play an important role even at high superficial velocities. At higher superficial velocities, the mass-transfer rates increase. This is consistent with an assumption of external mass-transfer limitations at the particle level. Note, however, that with increasing superficial velocity also the size and amount of particle clusters change. In our previous study, it was shown that the formation of clusters is highly influenced by the operating conditions, as well as cluster-related properties such as size and aspect ratio.[53] Therefore, (part of) the dependency of Kov versus superficial velocity might actually be indirect, i.e., due to changing characteristics of clusters.

In Figure , the axial profiles of the time-averaged ozone gas mass fraction profiles at different reaction rates are shown. In Table , overall mass-transfer coefficients obtained from these profiles at different reaction rates are shown. As expected, the conversion rate increases at higher values of kr.

Figure 4

Axial profiles of time-averaged ozone gas mass fraction at U = 5.95 m/s.

Table 4

Mass-Transfer Coefficienta

kr (s–1)	Kov (s–1)	γpf = Kov/kr	kmt·av (s–1)	kmt (m/s)
10	8.04	0.80	41.2	0.0058
50	36.17	0.72	130.9	0.0185
100	65.09	0.65	186.5	0.0264
200	115.94	0.58	275.9	0.0391
400	196.67	0.49	386.9	0.0548
1000	367.58	0.37	581.2	0.0823

Kinetic constant influence at U = 5.95 m/s.

The “plug flow contact efficiency”, γpf, decreases at higher reaction rates. As expected, the kinetic resistance decreases with increasing reaction rates, while the hydrodynamic resistance becomes more dominant. However, the dependence of kmtav on the reaction rate, kr, is inconsistent with the assumption of external mass-transfer limitations at the particle level. If kmtav represented the external mass-transfer limitation at the particle level, one would expect it to remain constant when only kr is changed. One can draw a more general conclusion, namely, that the resistances in series analysis of eq where kmtav is an external mass-transfer resistance at whatever level is not valid here.

There can be several mass-transfer mechanisms that cause the under-performance of the riser reactors. Falling clusters close to the walls may reduce local slip velocities, leading to lower particle-based Reynolds numbers and thus decreasing local mass-transfer coefficients (through Gunn’s correlation). Or maybe the bulk gas stream may circumvent dense particle regions (clusters), leaving the system without contacting with all particles. Gas back-mixing effects can also impede the reactor performance. All these hypotheses are assessed in the next sections.

Slip Velocity

We will first focus on the particle-level mass-transfer resistance. The particle-level mass-transfer coefficient depends on the slip velocity of the particle with respect to the surrounding gas. For larger slip velocities, the mass boundary layers around particles are thinner, which leads to increased mass-transfer rates. This dependence is captured by mass-transfer correlations such as the one reported by Gunn,[8] where increased slip velocities lead to larger particle Reynolds numbers and consequently larger interphase mass-transfer coefficients. As Helland et al.[62] suggested, a falling cluster exerts a local reaction force (via two-way coupling) on the gas phase, decelerating the gas motion, such that it could even follow the cluster trajectory. This phenomenon leads to a local drop of the slip velocity and therefore a lower local mass-transfer coefficient.

In order to analyze the influence of the slip velocity on the computed local mass-transfer coefficient (via the Gunn correlation), the particle-averaged mass-transfer coefficient is computed in each computational cell. In this way, we can evaluate whether the drop in the computed particle mass-transfer coefficient is the main source of lowered riser reactor performance.

In Figure , it can be seen that the particle-averaged mass-transfer coefficient close to the walls is significantly lower than that in the core of the riser. This confirms that cluster formation leads to lower particle-level mass-transfer rates because of a drop in the particle-based Reynolds number.

Figure 5

Particle-averaged mass-transfer coefficient from a CFD-DEM simulation of 40 s. U = 5.95 m/s.

In Figure , the probability distribution function (pdf) of the instantaneous particle Reynolds number at different gas superficial velocities is shown. It can be seen that these profiles describe bimodal data distributions. There is a high peak at rather small Reynolds numbers (5–20), while at high Reynolds numbers a peak appears which grows and shifts to the right with increasing superficial gas velocity. These profiles are consistent with a “slow-moving” solid phase, which can be characterized by particles that are immersed in dense areas where the local slip velocities are low, and a “fast-moving” solid phase that could be characterized by particles located in dilute areas where the slip velocities are relatively large.

Figure 6

Probability density distribution of instantaneous particle-based mass-transfer coefficient.

This pattern can also be observed in Figure , where the probability density distribution of the particle-based mass-transfer coefficient at the same operating conditions are plotted. We see in Figure that the pdf values of the instantaneous particle mass-transfer coefficient describe bimodal data distributions as well. It should be noted that the Gunn correlation has a strong dependence on the particle-based Reynolds number and the solids volume fraction. At similar values of the particle-based Reynolds number, dilute areas acquire lower mass-transfer coefficients than dense regions. At higher gas superficial velocities, clusters are less likely to form and slip velocities increase (see Figure ). It is then expected that the occurrence probability of the dense phase decreases as well. Therefore, looking at Figure , we can state that dense particle regions are characterized by slow motion and low local mass-transfer coefficients, while the dilute solid phase is characterized by high slip velocities and relatively high mass-transfer coefficients.

Figure 7

Probability density distribution of instantaneous particle-based mass-transfer coefficient.

In this subsection, we have confirmed that cluster formation leads to a lower particle-level mass-transfer coefficient. As can be seen from Figure , the mean particle mass-transfer coefficient of a dense solid phase ranges between 0.22 and 0.27 m/s, while the mean value of the same property for the dilute solid phase ranges between 0.43 and 0.49 m/s.

When comparing the particle-level mass-transfer coefficient measured here with the values reported in Tables and 4, we clearly see that the actual particle-level mass-transfer coefficients are much larger than those in the tables. The values reported in the tables are obtained by assuming that the dominant mass-transfer resistance is at the particle level. The disagreement shows that this assumption is incorrect. We conclude that, while the presence of clusters significantly influences the particle-level mass-transfer coefficients, the dominant mass-transfer resistance is not at the particle level.

Axial Gas Dispersion Coefficient

Axial dispersion might lead to a lower apparent Kov when experiments are interpreted using a plug flow model without axial dispersion. Some of the measured low Sherwood numbers might be explained by this type of “misinterpretation”. Therefore, the influence of the gas axial dispersion is evaluated in this subsection. In order to analyze whether gas dispersion effects play a major role in these deviations, a 1D convection dispersion equation can be employed to compute the apparent reaction rate in the riser reactor. Given the time-averaged axial ozone mass fraction profile, the influence of the axial gas dispersion coefficient can be determined. Changes in Dax are evaluated in order to analyze the deviations between the attained mass fraction profiles and those obtained by assuming a steady-state axially dispersed plug flow model.For Dax ≪ U2/(Kov⟨φs⟩) ,it has solutions λ1 ≈ −(Kov ⟨φs⟩) /U, λ2≈ U/Dax. Because λ2 ≫ −λ1, λ2 corresponds to a shorter length scale. In fact, this second solution (with positive exponential factor λ2) influences only the mass fraction near the exit of the column. That is, the contribution of exp(λ2z) will be significant only near the exit and is sensitive to the outlet boundary condition. Away from the exit, only a single-exponent solution, namely exp(λ1z), is relevant; therefore, wA,g is expected to decay exponentially in this region. This means that λ1z can be locally estimated using λ1Δz = ln⟨wA,⟩/⟨wA,⟩, such that Kov can be obtained from the characteristic eq asThe axial gas dispersion coefficient has been determined in CFD-DEM by injecting a pure ozone gas pulse of 0.01 s over steady-state simulations in an inert environment (no chemical reaction). The gas velocity fluctuations are expected to be larger in simulations with a higher degree of clustering (low U). In Figure , the obtained gas residence time, E(t), for U = 5.16 m/s is plotted. The axial dispersion coefficient was computed by means of eq :[66]where tm = ∫0∞t·E(t) dt is the mean gas residence time, L the riser length, ϵbed the bed porosity which amounted to 0.943, and U the gas superficial velocity.

Figure 8

Gas residence time distribution at U = 5.16 m/s.

The mean gas residence time was around 0.307 ± 0.11 s, and the axial dispersion coefficient amounted to Dax = 0.527 m2/s. By feeding this input parameter into the previous convection–dispersion equation, we can obtain the order of magnitude of the deviation in Kov (kr = 1000 s–1) when axial gas dispersion effects are accounted for in the interpretation model. In Figure , it can be seen that the change due to axial dispersion is small. Therefore, discarding gas axial dispersion in an interpretation model is not a major cause of overestimation of the conversion rate. For the case that dispersion has only a limited influence on the determined mass-transfer coefficient, we find that the relative contribution is approximately (Dax⟨φs⟩Kov)/U2. For the measured dispersion coefficient this gives a 13% deviation (with ⟨φs⟩ = 0.0535). Therefore, in this case, gas back-mixing is not a major cause for (apparent) mass-transfer limitations on riser reactor performance.

Figure 9

Gas axial dispersion influence (kr = 1000 s–1).

Gas–Solid Contact Efficiency

Gas–solid contact efficiency is related to gas bypassing:[49] some of the reactant will have an intimate contact with the catalyst particles, and the rest may leave the system chemically unchanged because of a very poor exposure to the particulate phase. In Figure , we present an illustrative snapshot when gas bypassing occurs.

Figure 10

Left: Ozone mass fraction at the particle surface. Right: Gas velocity field.

On the left-hand side of Figure , particles are colored according to their respective ozone mass fraction. These values can be assumed as the ozone mass fraction at the surface of each particle. On the right-hand side of the figure, the gas velocity field is shown. It can be seen that the gas flows at high velocities in the core of the pseudo 2D riser. Red-colored particles (rich in ozone content) are mainly encountered in regions exposed to the main bulk stream, where the velocities are higher. We can see that cluster regions are mostly composed of blue-colored particles that possess low ozone content presumably because of higher gas residence times or trapped gas pockets inside the clusters. Actually, Ouyang et al.[51] suggested that falling particle clusters could capture and retain gas, and these observations confirm this suggestion.

The gas–solid contact efficiencies have been computed for several values of the reaction constant, kr, at the same operating conditions (U = 5.95 m/s) and for different gas superficial velocities at a fixed kr = 100 s–1 to analyze the influence of cluster characteristics on the gas–solid contacting.

Reaction Rate Effect

In Figure , the pdf values of the cluster contact efficiency of all these simulation cases are shown. It can be noticed that at the lowest kinetic constant kr = 10 s–1, CFD-DEM predicts the major part of occurrences have contact efficiency values ranging between 0.8 and 1. Therefore, assuming ideal plug flow would in this case be a reasonable assumption if we want to estimate the riser reactor performance. However, when the reaction rate is increased, larger errors result. For instance, if a significant improvement is made on a catalyst by increasing its activity with corresponding change in kr from 10 to 1000 s–1, a plug flow model assumption can lead to larger overestimations because the gas–solid contact efficiency would be much lower than at low reaction rates (see Figure ).

Figure 11

Probability density distribution of gas–solid contact efficiency at different kr.

In Figure , the cluster-averaged contact efficiency for each simulation is plotted, where the error flags represent the confidence intervals of the 68.2% of the cluster contact efficiency data. It can be seen that it drops at higher values of the kinetic constant as previously stated by other authors.[16,42]

Figure 12

Cluster-averaged contact efficiency.

Influence of Gas Superficial Velocity

In a previous study[53] it was shown that complex clustering phenomena can be well predicted by means of CFD-DEM. In risers, the total cluster population increases at low gas superficial velocities.[64] Larger populations of falling clusters close to the walls can retain gas pockets of highly depleted reactant,[51] leading to inefficient gas–solid contacting. At higher gas superficial velocities, the system becomes more dilute and the particle shielding effect does not become that influential, as Figure reveals (see U = 6.74 m/s line). It is noticed then that clustering and consequently operating conditions play a major role in the performance of a riser reactor (see Figure ).

Figure 13

Probability density distribution of cluster contact efficiency at several gas superficial velocities at kr = 100 s–1..

These results show that clustering phenomena are a major cause of inefficient contacting. From Figure , it can be noticed that the cluster-averaged contact efficiency significantly increases at higher gas superficial velocities. Thus, the measurement of global mass-transfer coefficients requires an accurate estimate of cluster-related properties. This seems to be the cause of so much disagreement between global Sherwood number data. In systems where clustering and particle shielding phenomena are very pronounced or in systems in which the reaction rate is very high, the global Sherwood number will tend to zero.

Figure 14

Cluster-averaged contact efficiency at several gas superficial velocities at kr = 100 s–1..

Influence of Dilution Ratio

In this subsection, we present gas–solid contact efficiency results of CFD-DEM simulation at different dilution ratios of active particles (number of active over total number of particles). Diluted fluidized systems have been employed in the past to measure mass-transfer coefficients. Active spheres can be mixed with inert ones to experimentally measure mass-transfer coefficients. In CFD-DEM, all particles are numbered and tracked. By means of a simple algorithm, a fixed number of particles could be labeled as active or inert. Each particle label was permanent for the whole simulation, and the particles were assumed to be homogeneously mixed in the system.

It should be noted that the Gunn correlation was also utilized in these simulations to compute the particle-based Sherwood number. Mass-transfer correlations for dilute particle systems and Gunn correlation differ in the asymptotic behavior at low Reynolds number (2·ϵ/τ and 2, respectively). These differences in the diffusional contribution of the Sh number were negligible for this set of simulations. The reasons are as follows: First, the mass transfer at the particle level is not limiting, especially for particles inside clusters. Second, for particles outside clusters, the higher Reynolds number contribution to the particle-level Sherwood number is relevant.

In Figure , the probability density distributions of cluster gas–solid contact efficiency at different dilution ratios are shown. We can see that the gas–solid contact efficiency is higher at increasing dilution ratios. At a fixed catalyst activity (kr = 100 s–1), lower dilution rates (more active spheres) will lead to more severe particle shielding effects when clusters are formed.[16]

Figure 15

Contact efficiency pdf at U = 5.55 m/s and kr = 100 s–1..

As expected, Figure shows that the cluster-averaged contact efficiency drops at decreasing dilution ratio.

Figure 16

Cluster-averaged contact efficiency.

It can be seen that an increase of the dilution ratio effect is comparable to an increased catalyst activity (see Figure ). The gas inside the cluster becomes more depleted of reactant; consequently, the average gas reactant concentration is lower, leading to poorer gas–solid contact efficiencies.

In diluted systems, if the active particles are homogeneously mixed, the performance of a system where a fraction, φactive, of the particles are active and kr = 100 s–1 is expected to be similar to that of one with kr = φactive·100 s–1 where all particles are active. To prove this statement, we ran simulations at U = 5.55 m/s at equivalent kr values, where all the particle are active (see Table ).

Table 5

Mass-Transfer Coefficient, Plug Flow Modela

simulation	kr (s–1)	% active particles	U (m/s)
1	100	10	5.55
2	10	100	5.55
3	100	50	5.55
4	50	100	5.55
5	100	90	5.55
6	90	100	5.55

Kinetic constant influence at U = 5.95 m/s.

In Figure , the contact efficiency pdf of simulations 1–6 are plotted (see Table ). If we compare pdf profiles of simulations 1 and 2, it can be noticed that the effect of the particle dilution ratio is equivalent to the effect of the catalyst activity. The same trend is shown for the remaining simulation pairs. Thus, we confirm previous authors’ observations,[18] namely, that mass-transfer coefficients obtained from diluted systems should not be comparable to undiluted fluidized systems.

Figure 17

Contact efficiency pdf at U = 5.55 m/s.

Mass Transport Inside Clusters

The situation of mass-transfer resistance inside a cluster of reacting particles is qualitatively analogous to internal mass-transfer limitation inside a catalytic porous particle. If the analogy also holds beyond a qualitative similarity, a type of Thiele modulus could be applicable to determine the reaction effectiveness inside clusters and hence permit the development of a cluster-based mass-transfer model. In this case, a correlation between cluster size and gas–solid contact efficiency should be obtained, regardless of the gas-phase velocity. In this section, we will see that this is quite challenging because of the large data scattering that such a correlation shows.

The scatter plot of cluster contact efficiency, γcl, against the equivalent cluster diameter () in Figure shows no clear correlation between the two quantities. If the gas–solid contact efficiency was assumed to depend only on the cluster size, a trend should be visible. Moreover, the same master curve would be expected at different gas superficial velocities. However, it is observed that there is no clear correlation, especially in denser systems (i.e., at lower U values) where clustering phenomena are more intense.

Figure 18

Gas–solid contact efficiency versus equivalent cluster diameter of 500 random clusters at kr = 100 s–1: (a) U = 5.16 m/s, (b) 5.55 m/s, (c) 5.95 m/s, (d) 6.35 m/s, and (e) 6.74 m/s.

At higher gas superficial velocities (U = 6.35 and 6.74 m/s), the systems are rather dilute. Here, clusters are less likely to interact with each other, and this might be the cause of less data scattering (although this is still quite large, as is evident from Figure ).

Thus, it is worthwhile to show the causes of this scattering and why cluster-based mass-transfer models should not merely depend on the equivalent cluster diameter.

In Figure , we present a snapshot sequence of two clusters in a relatively dilute region of the riser domain. The gas velocity vector field is superimposed on the porosity field. It can be observed that in the first snapshot, the gas stream does penetrate into the cluster wake of the cluster located at the left. Although this riser section is quite dilute, we can observe (if we follow the sequence) how the gas passes through the smallest cluster located at the top as well. Van der Ham et al.[29] suggested that the increase in the gas–solid contact efficiency could be due to the breakup of cluster structures. Although cluster formation leads to poor gas–solid contacting, we also see that the gas can pass through the cluster structure, causing only a change or orientation of the cluster shape without destroying it. We observe that the cluster structure is quite dynamic and can adopt different shapes and orientations in time that can be more susceptible to gas permeation.

Figure 19

Porosity field with the gas velocity vector field superimposed in a dilute region of the riser domain. U = 5.16 m/s. Snapshot time frame is 0.01 s.

In Figure , another sequence in a denser region of the riser is shown. In this figure, it can be observed more clearly how the gas stream accelerates because of the high cluster content at the bottom of the riser. Denser regions will lead to not only enhanced bypassing but also the formation of gas streams with large velocities that can eventually pass through some of the clusters. This phenomenon causes cluster particles found upstream to be more easily accessible to the gas phase and able to experience a more efficient gas–particle contact.

Figure 20

Porosity field with the gas velocity vector field superimposed in a dense region of the riser domain. U = 5.16 m/s. Snapshot time frame is 0.01 s.

In general, we observe in our simulations a rather chaotic behavior of particle clusters. They not only form, grow, break up, or merge, but also they can adopt different shapes, densities, aspect ratios, and orientations. All these phenomena have an effect on the gas–solid efficiency of the cluster itself or/and other neighboring clusters that are found upstream. Although it seems clear that clustering phenomena enhance gas bypassing and poor gas–solid contacting, these phenomena feature such a broad variety of structures that it remains very challenging to develop closures for cluster-based mass-transfer models.

Conclusions

In this work we have performed CFD-DEM simulations in order to generate more insight about mass-transfer mechanisms that take place under riser flow conditions. The instantaneous cluster-level contact efficiency between the gas phase and cluster particles has been computed at several reaction rates and gas superficial velocities. This work explicitly confirms and corroborates suggestions made by other authors,[16,20,46,51,63,67] namely, that particle clusters have a large influence on the gas–solid contact efficiency and on global mass-transfer phenomena. We clearly showed that for the system studied here the increased mass-transfer resistance is due to the presence of particle clusters and not due to axial dispersion effects, or changes of the particle-level mass-transfer coefficient. Moreover, we showed that a decreasing gas superficial velocity leads to lower γcl values. At lower U values, the fact that clusters are larger and a less intense convective mass transfer exists inside these particle structures seem to be the main causes of obtaining such a pattern. In addition, increasing reaction rate was shown to decrease γcl, thus increasing the influence of hydrodynamic resistances at cluster level, as other authors suggested.[21,42] Diluted fluidized systems were found to lead to higher gas–solid contacting rates.[16,18] For the system studied here it was proved that the dilution rate effect is equivalent to reaction rate effects, as Venderbosch et al. suggested.[16] Therefore, the effect of dilution by inactive particles can be easily understood in terms of an equivalent decrease of the reaction rate.

Although in the literature there is general agreement about the relevance of clusters, it is less clear whether the mass-transfer resistance lies in the external mass transfer to the cluster surface[23−25] or whether clusters can be assumed as large porous spheres where only diffusional transport takes place.[20,42] In this work, we have shown that cluster contact efficiency does not correlate well with the cluster size because there is large data scattering. Thus, clusters cannot be assumed to be large porous particles, where effective molecular diffusivity is the only mass transport mechanism. Besides, convective mass transfer can play an important role when high γcl values are attained. Convective mass exchange between dilute and dense phases exists, and it could be enhanced by the formation of gas jets that pass through the cluster structures.

Particle clusters are transient entities that show a broad variety of shapes, sizes, and orientations, as suggested by other authors.[16,64,68] The large amount of properties that characterize particle clusters, altogether with their location, local density, and proximity of high-velocity gas streams can cause interactions of a very diverse nature. The scattering pattern of gas–solid contact efficiency data suggests that convective mass transfer inside clusters can be enhanced or limited by all cluster properties previously mentioned, obtaining quite unpredictable behavior if we analyze only a single parameter, for example, cluster size.

Although we have shown that clusters enhance gas bypassing and result in poor gas–solid contacting, we find it challenging to develop mass-transfer closures at the cluster scale level. It seems very hard to capture the influence of a cluster using simple parameters such as cluster size. The main reason is that the cluster contact efficiency is very much influenced by convection through the cluster and that this convection depends on the configuration of clusters downstream. Our tentative conclusion is that accurate coarse-graining of the influence of particle clusters is difficult and that in fact CFD-DEM is the tool to predict the performance best. Related to the particle-based closures used in CFD-DEM, we point out the following behavior: Most mass transfer seems to take place at the boundaries of clusters, where flow can still partly penetrate the clusters. At these locations, the solids volume fractions quickly change. However, the particle-level mass-transfer correlations used in CFD-DEM were developed for (locally) homogeneous systems. This raises the question whether the particle-level correlations are accurate enough. Therefore, we recommend performing direct numerical simulations of freely evolving clusters to validate local particle-level Sherwood correlations of heterogeneous particle structures.