On the fast modeling of species transport in fluidized beds using recurrence computational fluid dynamics

Funding information Christian Doppler Association; The Austrian Research Promotion Agency (FFG) Abstract Due to variety of scale dynamics evolved in gas–solid flows, most of its numerical description is limited to expensive short durations. This has made the slow processes therein, such as the chemical species conversion, to be out of an appropriate reach. In this work, an application of the transport-based recurrence computational fluid dynamics (CFD) has been introduced for the fast modeling of passive scalar transport, which is considered as species conversion and heat transfer in fluidized beds. The methodology discloses the recurrent dynamics during a short-term full CFD simulation as Lagrangian shift operations upon which a passive scalar can infinitely be traced. Apart from convecting, a proper approach based on the turbulent kinetic energy of tracked dynamics is introduced for modeling the physical diffusion of the scalar transported. Our outcomes have revealed a subtle chasing to the full CFD species simulation with a speed-up up to 1,600.

The arrival time of a space probe traveling to Saturn can be predicted more accurately than the behavior of a fluidized bed chemical reactor!". 1 Fluidized beds are systems with a bed of granular particles initially resting on a perforated bottom plate. When the inlet fluid is passed upward through the bottom plate, it suspends, or fluidizes, the particles to allow a liquid-like behavior of solids with a high effective contacting process. This contributes to the chief advantage of fluidized beds in providing a rigorous mixing and favorable heat and mass transfer characteristics. With these attributes, fluidized beds have occupied a high-ranking beneficial position in the chemical processing applications, such as the well-known fluidized catalytic cracking of petroleum oil, granulation for powder production, coal carbonization and gasification, coking, coating preparations, and also in nuclear fuel fabrication. However, fluidized beds reactors are challenging to design and scale-up. The dynamics involved composes into complex physical phenomena because of the multiscale nonlinear interactions between the solid and fluid phases. Therewith, the evolution of structures, on the other hand, is strongly dependent on the particle properties such as size, shape, and density. 2 This nature has engrossed Geldart 1 and many other scientists [3][4][5] in order to understand and predict accurately all involved dynamics of fluidized beds. Principally, the fluid flow interacts with the solid particles by the interstitial fluid drag. The particles among themselves induce kinetic, collisional, and frictional stresses and undergo to complicated deformations and influences of pressure gradient, rolling, and gravity forces. The heat is transferred within each single phase and also between them. These (spatial and temporal) small-scale (microscopic) interactions compose to larger (mesoscopic and macroscopic) interactions which take place between clusters of solid (emulsion) and fluid bubbles. The temporal evolution of the (microscopic) particle's collision dynamics lasts about few milliseconds, much shorter than the bubble evolution. While the bubbles in turn, move faster than the fluid slugs in the bed and phenomena as heat transfer and chemical species conversions can take minutes or hours.
Nowadays, computational methods based on fundamental principles for resolving each particle-particle and fluid-particle interactions have allowed to explore many details of the underlying dynamics 6 in fluidized beds. Using the computational fluid dynamics (CFD) techniques or lattice Boltzmann methods, the fluid hydrodynamics (obeying Navier-Stokes equations) is investigated on finer grids than the particles separations disclosing all forces exerted on the particles surface (drag), and deducing improved empirical correlations. 7,8 The solid side contacting forces of collisions are pictured as either an eventdriven hard sphere model or a time step-driven soft sphere approach.
With all relevant active forces on individual grain, each particle's trajectory is computed with the aid of Newton's second law, which is referred as the discrete element method (DEM) introduced by Cundall and Strack. 9 Coupling CFD-DEM has been the most appropriate numerical modeling for the simulation of fluidized beds; however, the huge inherent computational effort makes its application impractical for large particulate systems. 10 It is, therefore, common to investigate fluidized beds in large processing units using averaged equations of motion. 11 Namely, by averaging Navier-Stokes equations over several particle diameters while the particles may still be treated as discrete elements (unresolved CFD-DEM), 12 the computational cost can be reduced, but now the solid-fluid interactions are modeled using empirical closures. Another alternative method is considering the solid particles as a separate continuum, like the fluid phase, and the motion of particles is analogously averaged out, in the so referred two-fluid model (TFM). 13,14 Therein, the solid stresses arising from particleparticle collisions, the transitional dispersion of grains, and rotational speeds are deduced by adopting the kinetic theory of granular flow. 15 Since TFM allows for coarser grids without demanding to track individual particles, it ends to be more suitable for large systems than CFD-DEM; however, it is still forbidding for huge industrial-scale reactors.
With the considerable developments of numerical methods, the coarse-graining models have appeared to upscale both CFD-DEM and TFM to macroscopic sizes. For instance, the parcel-based DEM methods handle in tracking a parcel of several solid constituents which interact with modified material parameters. [16][17][18] Likewise, the coarse-grained TFM leads to simulate larger systems on much coarser grids, 19 but with the requirement of fundamental considering to the subgrid heterogeneities. [19][20][21][22][23][24][25][26][27][28][29] Even though, using these methods, the plant-size reactors may be simulated, but they are bounded by shortterm investigations due to the huge computational power needed.
The time scale simulation remains of small rates in order to capture mathematically all relevant dynamics of collisions; and therefore, the slow processes, which take hours in large fluidized beds, are still inaccessible. Seeking a remedy, some of us have introduced the idea of recurrence CFD (rCFD) for pseudoperiodic flows. [30][31][32][33] Namely, the continuous reappearing (recurrent) structures of such dynamics, for instance, the gas bubbles evolution in bubbling fluidized beds, is considered and analyzed on the base of a short-term full CFD simulation, that is, CFD-DEM or TFM. Using these data fields, the recurrence plots 34 are generated to predict a recurrence path upon which a passive scalar can be transported till infinity. This new methodology has been successfully applied for the fast passive transportation in multiphase flows 30 and heat transfer prediction in fluidized beds, 31,32 basing on CFD-DEM simulations. Therein, the advancing timeextrapolation procedure, by which the passive scalar is propagated on candidate recurrent patterns, is established using three proposed models. The flow-based category 30 where the recurrent flux fields have to be provided in order to resolve either a convection-diffusion equation (Eulerian model) or a stochastic differential equation for a fluid-parcel trajectory (Lagrangian model). In the third model named as the transport-based rCFD model, 32,33 the recurrent fields are shorten to only start-end positions information with no a posteriori need to resolve any equation. In this work, we focus on the test and application of the third rCFD model for the fast passive prediction of species and temperature (considered as a passive scalar) reconstructed transport in a lab-scale bubbling fluidized bed. The TFM is used to perform a short-term full CFD simulation, during which all database is collected and analyzed in terms of the recurrence properties and rCFD process.
The main content, herein, implicates enhancements of the transportbased rCFD methodology in relevance to evaluating the proper physical diffusion approach for the passive scalar transported. Namely, after tracing the scalar, different diffusion approaches have been applied on the base of an approximate local diffusion and the global balance of the mass concentration scalar and enthalpy in the domain.
The methodology is adapted for the reconstruction procedure of species and temperature path, separately, and the advancing timeextrapolation process for the case of a continuous coming-in inlet species on the gas phase. Doing so, the results have shown a successful consistency with the full CFD evolution by consuming very cheap and short runtime computations.
The rest of the paper is arranged as follows. First, the TFM and rCFD methodologies, particularly the transport-based method, are explained in Section 2 and Appendix. Then, details of the simulated fluidized bed with the results of species and heat transport rCFD modeling are presented and discussed in Section 3. Therein, the main feature is about finding the proper approach of the physical diffusion part in the transport-based rCFD methodology. On the other hand, the adapted algorithm of rCFD to transport the solid and gas temperatures as interacted two scalars is also outlined. Finally, relevant results are summarized and conclusions are given in Section 4.
2 | THEORY 2.1 | Two-fluid model density and velocity fields, respectively, in correspondence to the subscript gas (g) or solid (s) phase. In some cases of a cold fluidization is accompanied with nonreactive chemical species conversions, the involved transport dynamics, for each jth specie mass fraction γ j , is described in Eulerian framework, likewise (Equation A14 in the Appendix). The species, therein, are principally diffused obeying the dilute approximation of Fick's law, and roughly consider an identical constant value of the diffusion coefficient D j , on the solid and gas phases. For another particular fluidization, when a participating heat transfer is implied, the temperature on the gas T g and solid T s phases follows, as well, an Eulerian frame which describes the enthalpy transport in the Appendix), assuming constant specific heats c s=g p . The key element, therein, is about measuring the interphase heat transfer coefficient denoting the Nusselt number, Nu, where the solid thermal conductivity κ s is rigorously considered fixed and equivalent to the gas phase one κ g .

| Recurrence CFD
In previous work, 30 the idea of rCFD had been introduced as an approach that helps to model the long-term processes in massively large systems, such as the chemical species conversions. Its application implies disclosing the recurrence properties of the system, which access the similarity between states and determine its periodicity parameters. Namely, on the base of a short-term full CFD simulation, the evolution in time t of active patterns at different probes in the domain is considered to decide the recurrence period τ rec , which has to exceed several pseudoperiodic periods τ p − p of the system. This last is identified by performing a spectral analysis of the signals (see an example in Figure 4d for the probes taken of the studied fluidized bed). And, hence, the recurrence period should span multiple periods of the corresponding lowest-lying peak frequency (f crit ) of the energy spectra, that is, The recurrence time step Δt rec is estimated from the consideration that a given field φ is not changing too strongly from one time step to the next, inside τ rec . Meaning that, and hÁi denotes the temporal averaging operator. In order to quantify the similarity of flow patterns between states, for example, t and t 0 = t + Δt rec , the recurrence norm is defined as, 30 where φ can be chosen depending on the phenomenon of interest. Basing on the recurrence matrix, a recurrence path can be extracted to convey the system's evolution far beyond the recording time τ rec . The strategy adopted was previously explained in, 30,31 and here we recall it in a simplified flow frames, f, aspect, shown in Figure 1. Namely, basing on R(500, 500), and starting from the end f 500 , we look for the most similar state by jumping backward toward the first half of τ rec , and considering the original flow of the maximum R frame, that is, f 130 ≈ f 500 . This frame will be taken as a base start field for the propagation of the passive scalar. Next, we pick an inter- In the so-called transport based rCFD model, 32 Step (1) in Figure 2). As graphically interpreted, this volume, V inlet tr , can be smaller or bigger than the neighboring-cell volume V bc . In case of V inlet tr = ϵ s=g ku in kA in Δt rec < V bc , the inlet tracers will move during Δt rec and end to those (neighboring-) cells, with a normal physical behavior.
However, when V inlet tr = ϵ s=g ku in kA in Δt rec > V bc , which occurs at significant inlet velocities or long Δt rec , a big portion of the inlet tracers will travel deep inside the domain and they will all end to one cell producing unphysical high accumulation hit of transported concentration. In this regard, a temporal delaying approach on these tracers, by random periods between [0: Δt rec ], is required in order to distribute them homogeneously along its path next to the inlet. Regarding the outlet flux, we respect to those internal tracers that hit the outlet surfaces, where v i inlet is the volume of inlet tracer.
2. Transport the concentration γ following the shift positions of internal tracers. Likewise, for cells of more than one shift operation and, for instance, the cell is hit by an inlet tracer, γ is given as where v i tr is the volume of internal tracer and γ ci is the concentration transported from the cell c i , with the volume fraction ϵ i s=g (following the phase on which the propagation of passive scalar is traced).
3. Fill the holes. Meaning that, those resultant cells which are not hit by any transport information, its value will be interpolated by the surrounding cells concentration. Meaning that, the concentration in the cells where the outlet tracers were located at the beginning of Δt rec , and before leaving the domain.

One step correcting diffusion controlled by the physical global bal-
ance between the coming-in, accumulated, and going-out mass concentration in the total domain, 32 that is, and which can be written as If the left-hand side (LHS) in Equation (8) exceeds the target (right-hand side RHS), the diffusion operator loops over all the internal faces and shifts/swaps rigorously (out) a specific concentration portion from the high-concentration cell, as Here, c i and c i + 1 are the two adjacent cells which share the surface, and f is a constant diffusion factor of magnitudes 1/2 n . In the contrary case, when the system indicates an excessive target the shift direction will be (in) toward the low-concentration cells.
In different than the previous transport-based rCFD version reported in Reference, 32 the aforementioned methodology implicates physically consistent enhancements regarding the inflow and outflow modelings. Namely, the inflow modeling before was by setting a concentration source next to the inlet and convey it following the internal tracers. Therefore, in some situations when the portion of these tracers is insufficient, a discontinuity in the concentration transport can be produced. On the other hand, the outflow was globally aligned to a constant value of the inflow, and was uniformly distributed over the outlet adjacent cells. In this case, the going-out mass concentration is roughly approximated with no local accuracy. All these features are resolved in the current version to make the methodology more consistent and associated totally to the particles dynamics. The functionalities therein besides the incoming experienced diffusion are implemented in the frame of user-defined functions code in the software ANSYS/Fluent.

| RESULTS AND DISCUSSION
We consider the simulation of a lab-scale (Geldart B particle) bubbling fluidized bed, displayed in Figure 3a, within the framework of the following separated cases. First, introduce a mixture of two species on the gas phase by injecting a red-color species into the lateral gas inlet and a blue-color species into the bottom gas inlet. Second, set an initial mixture of a red-color (small rectangular region) and blue-color (elsewhere) species, on the solid phase ( Figure 3c). Third, consider the enthalpy transport with an initial hot rectangular region of the solid phase, that is, T s = 500 K (Figure 3d), and a continuous cold air, that is, T g = 300 K, passing within all air inlets.
To do so, the mathematical models described in Section 2.1 and Appendix are numerically resolved in Cartesian coordinates r = (x, y, z), For details about the numerical methods, algorithms and solver, the reader is referred to. 36 The consistent time step Δt is defined by an appropriate value Δt = 5 × 10 −4 s, which satisfies the Courant number Co = u g Δt/Δr ≤ 0.1, at the operating conditions. We firstly run ). They have been utilized, a priori, to decide the recurrence properties, that is, Δt rec and τ rec , given by Equation (1) and (2), respectively.
Following the mean-square values and temporal derivatives, ≈0:077s , for p, w s , and ϵ s , respectively. Using the aid of different similar analyzing points, not shown here, the consistent value is decided as Δt rec = 0.008 s, that resolves adequately the global pseudoperiodic characteristics. All simulation details and recurrence properties are summarized in Table 1, where τ rec = τ γ g rec , is identified for the period used in the gaseous species modeling, while τ γ s rec and τ T rec are defined for those periods in the solid species and temperature transport, respectively. In different than the case of continuous coming-in gaseous species, the solid species and heat transfer problems are set by an initial value, as a small rectangular region, for the solid mass fraction and solid temperature. In the solid species problem, this initial mass fraction (unity) will be diffused in time till the reach of the uniform concentration.
Hence, in order to avoid this final attainment, τ γ s rec is decided shorter than in the gaseous species. Likewise, in the temperature transport problem, the heat will be transfered inside the solid phase itself and toward the gas phase. This last in turn, is continuously going-out through the outlet, and thus τ T rec is selected shorter in order to avoid the uniform cold extent.  Figures 5a,b,  deep fluctuated pressure, the solid flux ϵ s u s , reveals a coherent relevance to the uniform flow periodicity (note the signal of w s in Figure 4b). Namely, it highly discloses the fingerprints of rising bubbles passage, with high effective solid mixing through the wake effect particles, emulsion drift particles, and the bubble eruptions particles. 37 The gas flow, however, evolves as dense flow in rich solid areas, visible bubbles, and throughflow which bypasses through the bubbles. 38 Hence ϵ g u g includes high fluctuated dynamics that absentees the recurrence between states and can be less appropriate for the periodicity disclosure.
The volume fraction, ϵ s/g , in turn, stands to be a strong and directly sensitive field to the bubbles identification, where its recurrence norm value is identical for the gas and solid phases,

| Reconstruction of species transport on gas phase
In this context, we continue the rCFD application by assessing the enhanced methodology of transport-based rCFD model (Figure 2) for the reproduction of gaseous species transport all along τ γ g rec . In other words, the tracer-based shift information of N tr = 10 6 tracers, linked to the gas dynamics, and stored for 400 flow frames are retrieved in sequential order, that is, the recurrence path is {f 1 , f 2 ,…,

| Global mass species conservation
Following the global conservation condition of integral mass concentration over Δt rec , and which is given in Equation (7), the face-swap concentration Δγ (Equation 9) between two cells, is passed in or out, as explained in Step (5). Note that this criterion is physically based and comprehensive to correct the mass balance deviation in localcareless way. Talking in global sense, the diffusion amount needed can notably vary between flow frames, in relevance to the species lifetime that is essentially associated to the gas dynamics (i.e., ϵ g , ϵ g u g ). Therefore, instead of using a constant diffusion factor f, its value is defined in an approximate dynamical way. Namely, given an upper, 1/2 3 , and lower, 1/2 7 , bound, the value of the diffusion factor is increasing/decreasing if the RHS of Equation (7), tends to be bigger/ smaller than the LHS, from one frame to the next. Following that attitude, the diffusion procedure is applied standalone in only one loop step for each frame, and the pertinent outcomes have revealed a feasible chasing to the actual full CFD species simulation. For example, by considering a visual instantaneous comparison (at t = 6.2 s) between the full CFD results and rCFD, as shown in Figures 6a, one can observe the good tracing using very short runtime; however, an excess of the concentration transported is sustained (Figure 6a, middle). This excess is clearly pronounced in the evolution of the global mass quantity, that is, m γ = P cells i = 1 ρ g ϵ i g γ ci V cell computed all along τ γ g rec , and represented in Figure 6b (bottom). Therein, the rCFD results indicate a higher concentration than it is out to be in the mass fraction species. In a detailed review, as well, the high concentration is remarked on γ histogram chart generated for the same instant in On the basis of that, one can conclude that applying one loop diffusion standalone cannot return the proper diffusive consequences even though it is a very fast procedure from practical point of view.
Another approach is sought in the same line by doing multiple diffusion loops till the convergence between the two sides in Equation (7).
For the sake of quickness and simplicity, the diffusion factor is chosen as a constant value, that is, f = 1/2 9 , which is proportional to the molecular diffusion coefficient, that is, f ≈ D j Δt/A. It has to be certainly smaller than the lower bound in the dynamical procedure above, where A here is supposed to be the face area through which the concentration swapping takes place. This factor contributes into the molecular diffusive part, and by repeating it several times in consistency with the mass balance converging, the proper amount of diffusion in each particular frame will be recovered. Each frame requires a different number of diffusion loops to eventually make the approach somewhat expensive (see the performance Figure 16a in Section 3.5).
If we explore the same mere optical comparison in Figure 6, the new outcomes, referred as confined molecular diffusion, give a subtle capturing to the actual full CFD evolution in instant (Figures 6a, right) and global (Figure 6b, bottom) terms. They deliver a nice reproduction to the species path, but in such locally rough transport. Looking closely, the procedure relatively dissipates the concentration in some regions, for instance, outside of the bed, more than what is out to do. While, in other parts close to the source, the operator does not dissipate sufficiently. Note the histogram, particularly, in Figure 6b (top). Finally, these globally controlled models are feasible in terms of the main purpose of rCFD and saving considerable computational effort, but in spite of that a better calibration using a local diffusion is somehow mandatory.

| Local diffusion
In order to derive, accurately, all compensative diffusion effects, we firstly execute a one globally controlled diffusion loop characterized by a dynamical determination of the diffusion factor (see Section 3.1.1). This step satisfies an initial correction to the methodology's error in a featureless global way. Then we follow it by a kind of a fast and local physical diffusion approach. To that end, the oscillated (random) walk of species parcel along its path, from the beginning to the end Δt rec , is assumed as a turbulent diffusion and computed on the base of the turbulent kinetic energy of the tracked tracers. Taking the online temporal sampling of tracers velocity fluctuations u 0 p , that is, velocity variance, during Δt rec , the turbulent kinetic energy results as The sample variance of each velocity component, as a scalar hφ 0 2 i, is computed in online running way. Namely, we adopt a fast and cheap algorithm for evaluating the running variance directly at the same arrive moment without the need to save data for a second pass. 39 The method implies calculating of two variants, that is, M and S, for each individual tracer, and which are being updated each particle time step along Δt rec . Giving initial values, M 1 = φ 1 and S 1 = 0, the subsequent i moment for these quantities is computed as the consequential turbulent diffusion is estimated by a mixing length model, 40 that is, passing in γ ci and out γ ci + 1 . Here, Δr ci,ci + 1 is the vector distance from c i to c i + 1 and n A is the surface-normal vector on the same direction. Probing different constants C in Equation (12), with a special near-wall treatment, that is, Δ = 0:5V 1=3 cell , the tests are demonstrated in Figures 7 and 8, together with the outcomes of the previous approaches. In global-scale investigations, the suitable concentration modeling has suggested a constant value of C = 0.2, that best fits the full CFD evolution (Figure 7a), Therewith, the new evolution of the global mass shows identical behavior to the confined molecular diffusion by performing only two diffusion loops (one dynamical + one physical), to eventually be a very economic procedure (see the performance Figure 16a in Section 3.5).
On the other hand, the instant reproduction of species structures is enhanced using the new strategy, as can carefully be

| Reconstruction of species transport on solid phase
In this problem, and identical to the reproduction procedure of gas-  Figure 9 for the field (ϵ s γ). It can be noticed an acceptable matching between the rCFD results and full CFD counterparts. Moreover, the mid-width profiles, averaged over 200 flow frames, for hϵ s γi, and extracted inside the bed, at posi- 0, y, z = {0.01, 0.02, 0.03, 0.04, 0.05, 0.06, 0.07, 0.08, 0.09} m), have revealed a satisfied rCFD modeling, as represented in

| Reconstruction of heat transport
Different than the previous offline rCFD methodology, which predicts the transport of standalone-phase species, the heat transfer, herein, sustains an interaction between the gas and solid phase. In order to capture this interaction, the values of heat transfer coefficient Nu are computed during the online tracking procedure, that is, τ T rec , beforehand, and stored as an auxiliary online data. In that stage, both the gas and solid dynamics are simultaneously tracked by injecting two groups of N tr = 5 × 10 5 , massless tracers, to eventually produce two shift information data sets. Afterward, the reconstructed transportation of T g and T s scalars with the selfsame full CFD initial conditions is commenced in the framework of the following modified rCFD algorithm: 1. Transport T g following the inlet gas tracers shift positions, and employing a two-halfs temporal integration step for the transportation. Namely, we use the analytical solution (complete mixing) of temperature evolution in the bulk granular material 41 to describe the interacting transport of T g from c 0 (start cell) to c 1 (end cell) during Δt rec , as follows: where T c0 = ϵ g,c0 ρ g c g p d 2 s 6ϵ s,c0 κ g Nu c0 and T c1 = ϵ g,c1 ρ g c g p d 2 s 6ϵ s,c1 κ g Nu c1 ð15Þ are the time scales of T g change. Here, T g,c0 = T g,inlet = 300K for all inlet tracers. If we name the exchanging inlet temperature as T inlet g,c1 = T Δtrec g,c1 , the resultant T g transported into the target cell c 1 , in case of multiple hits, is given by During this transportation, the quantity of heat transfered due to the interaction with T s is also summed up as 2. Transport T g following the shift positions of internal gas tracers, and employing the same method before. The exchanging tempera- is evaluated from Equation (14) and the resultant T g convected into c 1 follows: Likewise, the amount of heat transfered Q g,c1 is summed up in c 1 using Equation (17).
3. Fill the holes of T g and Q g fields. 4. Consider (store) T g of the start cell outlet gas tracers T outlet g,c , which leave the domain during Δt rec .

F I G U R E 9
Comparative instantaneous pictures between the actual full computational fluid dynamics (CFD) red-color solid specie (ϵ s γ) (left) and recurrence CFD (rCFD) counterparts (right), at the moments, t = 3.16 s (a) and t = 4.6 s (b). rCFD results are obtained adopting the local diffusion approach (C = 0.8) and using N tr = 10 6 [Color figure can be viewed at wileyonlinelibrary.com] 5. One step correcting T g controlled by the physical global balance between the coming-in, accumulated, exchanged and going-out gas enthalpy through the total domain, that is, which can be written as Similar to Step (5) in Section 2.2, the swapped/shifted portion of temperatures reads, where f is dynamically decided, similar to Section 3.1.1.
6. One-step physical diffusion modeling that compensates the heat transfered by conduction in the gas phase as where D g,ci = κ g = ρ g c g p and D t g,ci = C g Δk 1=2 g =Pr t . Here, D g,ci is the molecular thermal diffusivity and D t g,c i is the mixing length turbulent diffusion which is computed on the base of the turbulent kinetic energy of gas tracers k g with an empirical constant C g = 0.2. Pr t = 0.4 is the turbulent Prandtl number. Afterward, we iterate one-pass over the solid phase, transporting T s following the solid tracers and adding the quantities of Q g , as follows: 1. Transport T s following the shift positions of internal solid tracers.
The exchanging temperature will be T Δtrec s,c1 = T s,c0 + and the resultant T s convected into c 1 is given by which can be written as Identically to the gas phase, the specific portion is swapped in or out toward T s cells and f is dynamically decided.
5. One-step physical diffusion modeling that compensates the heat transfered by conduction in the solid phase as  Figure 12 and z = 0.3 in Figure 13).

| Time extrapolation of species transport on gas phase
The main concept of rCFD is tracing the passive scalar on stitched recurrent moments prolonged far beyond τ rec , in a fast and low-cost modeling. Using the robust and randomly based approximate, depicted in Figure 1, we start from f 1 on the recurrence matrices ( Figure 5) and move with arbitrary intervals and jumps, to generate in due course, a sequence of, for example, 2,000 recurrent flow frames.
Sweeping over this path, the gaseous passive propagation is predicted till 22 s on the base of τ γ g rec = 3:2s database and time step Δt rec = 16Δt. For a comparative perspective, the concentration is traced using sequences based on R ϵs ð Þ and R ϵsus ð Þ recurrence matrices, and applying the favorable (one step dynamical + one step physical physical) diffusion approach in Step (5). While an accurate description of species trajectories might not be of significant interest in a long-term transfer, but its distribution on the mean picture yes. To that end, the temporal mean patterns extracted at the mid-width (y, z) plane, are drawn in

| Impacts and efficiency
At the crossroad between modeling accuracy and the storage needed for the empirical framework of transport-based rCFD, a Lagrangian coarse graining procedure is noteworthy. 32 Namely, by reducing the total number of tracers N tr , up to which the internal seeding volume is V tr ≥ V cell , the size of recurrence data frames becomes meaningfully smaller. And hence, the entire transport runtime becomes of very short extent. In order to test this impact, the reconstruction initiative of gaseous species transport is repeated on the same 400 frames using a total number of tracers, N tr = 70,000, and adopting the favorable local diffusion approach. By doing so, we redo the selfsame aspects, previously outlined in Figure 8, and show them in Figure 15.
From a precise look, the coarsening stratagem implies a bold transport for concentration that needs to be smoothed in-between the structures ( Figure 15a). This is obviously returned to the less number of tracers that could convey more information on ϵ g from the neighboring hit cells. In general, all procedures therein are conservative to γ mass due to the proper local diffusion applied. Nonetheless, the spatial coarse-grained transport which transfers the inertial energy following the tracers, as a large negative straining, has its influence on γ mapping. In spite of that, an attainable tracing of gaseous species transport can be obtained using N tr = 70,000, and relatively gives analogous results on the mean fields (Figure 15b), spending even cheaper cost (note the performance Figure 16a,b).
Going through the rCFD performance, the computational time consumed for all current applications is graphically outlined in Figure 16. We firstly carry out an initial simulation of 3 s, needed to hit the pseudoperiodic flow, by the use of 12 CPUs and persisting over 7 hr CPU time per processor. This is followed by a short-term full CFD simulation lasting about $3 s, thereby the 400 flow frames of recurrence database are cached. For cases of solid species transport and heat transfer, the database is shorten to 200 flow frames with a lasting runtime $1.5 s. After establishing these requisites, the efficiency of the passive transportation rCFD is depicted by a comparison of one CPU runtime between rCFD and the long-term full CFD counterpart. Note that the transport-based rCFD methodology is also working, in parallel mode, on the same number of CPUs. Experiencing different approaches, as displayed in Figure 16a, the most feasible rCFD indicates a reduction of CPU time from 7 hr to 28 s, in the framework of gaseous species reconstruction and N tr = 10 6 tracers.
This dramatic lowering in cost implicates a speed-up ratio of 900, which is further raised to 1,575, using N tr = 70,000 tracers. For the case of applying a confined molecular diffusion, the computational runtime is slightly increased to 180 s, however, it is still located in the limits of a very reliable application. Regarding the solid species reconstruction, the feasible rCFD application implies a reduction from 3.5 hr to 9 s with a 1,400 speed-up ratio, while the temperature reconstruction reduces 8 hr CPU time to 18 s with 1,600 speed-up ratio.
If we extend the gaseous species prediction much longer time, for instance, till 22 s, with a time-extrapolating rCFD, the computational efficiency is rather higher. It drastically saves an amount of 37.3 hr − 148 s = 37.29 hr (see Figure 16d) for N tr = 10 6 , while in the coarsegraining case the cost becomes almost negligible. F I G U R E 1 6 Performance of recurrence computational fluid dynamics (rCFD) modeling for the reconstruction procedure of (a) the gaseous species using different assumptions and inputs, (b) the solid species using N tr = 10 6 and local diffusion approach, and (c) the temperature using N tr = 5 × 10 5 on each phase and the local diffusion approach. (d) The performance of the long-term time-extrapolation rCFD modeling for gaseous species [Color figure can be viewed at wileyonlinelibrary.com] possible fast and active interactions between chemical species. Therewith, the physical properties can be changed, and more consequences like the trigger of new solutions and heat source have to be considered. In this article, we have taken a step forward preparation to access the modeling of interacted species transportation with chemical reactions. For this regard, further constituents on reaction kinetics as production/consumption mass phase have to be tracked leading in consequence, to higher storage amount of recurrence database. This will be our primary concern in the future.

| CONCLUDING REMARKS
In this work, we have employed the transport-based rCFD modeling for the fast and low-cost prediction of passive species and heat transport in bubbling fluidized beds. To do so, a typical short-term full CFD simulation has been performed using the framework of TFM, to compose afterward into the backbone recurrence database. In the conceptual transport-based rCFD methodology, massless particles are used to chase the gas and solid dynamics during each recurrence period, and store the trajectories as start-end positions. 32 Herein, we have implied a physically consistent modeling for the inlet and outlet flows to be completely based on the tracers themselves. Doing so, the timeextrapolating rCFD modeling for passive transportation is applied in a more proper aspect. After the establishment of species and temperature convection, the foremost appropriate diffusion is found to be in performing one globally based diffusion step plus a kind of local physical diffusion concern. Similarly to the mixing length assumption, 40 the local physical diffusion has been approximated on the base of the turbulent kinetic energy of tracers, sampled all along the recurrence period. As a consequence, the rCFD outcomes of species and temperature reconstruction path have revealed a very reliable agreement with the full CFD (TFM) evolution by consuming a very low computational cost, and speed-up gaining up to 1,600. For large-scale simulations which commonly use coarse grids, the efficiency of rCFD can be even improved since the small scales, limiting Δt rec , are not resolved but modeled in the context of filtered TFM. 19,28,29 Considering this fast and very cheap capacity, the transport-based rCFD method can constitute a strong helpful tool that allows to access and resolve many queries in the huge industrial applications.