# Three-dimensional unsteady numerical simulation of a 150 kW full-loop chemical looping combustion pilot with biomass as fuel: A hydrodynamic investigation

Liyan Sun, Enrica Masi, Olivier Simonin, Øyvind Langørgen, Inge Saanum, Nils Erland L. Haugen

► To cite this version:

Liyan Sun, Enrica Masi, Olivier Simonin, Øyvind Langørgen, Inge Saanum, et al. Three-dimensional unsteady numerical simulation of a 150 kW full-loop chemical looping combustion pilot with biomass as fuel: A hydrodynamic investigation. Chemical Engineering Science, 2022, 260, pp.117835. 10.1016/j.ces.2022.117835. hal-04533880

HAL Id: hal-04533880

https://ut3-toulouseinp.hal.science/hal-04533880

Submitted on 5 Apr 2024

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
# Three-dimensional unsteady numerical simulation of a 150 kW full-loop chemical looping combustion pilot with biomass as fuel: A hydrodynamic investigation

Liyan Sun, Enrica Masi, Olivier Simonin, Øyvind Langørgen, Inge Saanum, Nils Erland L. Haugen

► To cite this version:

Liyan Sun, Enrica Masi, Olivier Simonin, Øyvind Langørgen, Inge Saanum, et al. Three-dimensional unsteady numerical simulation of a 150 kW full-loop chemical looping combustion pilot with biomass as fuel: A hydrodynamic investigation. Chemical Engineering Science, 2022, 260, pp.117835. 10.1016/j.ces.2022.117835. hal-04533880

HAL Id: hal-04533880

https://ut3-toulouseinp.hal.science/hal-04533880v1

Submitted on 5 Apr 2024

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
Three-dimensional unsteady numerical simulation of a 150 kW<sub>th</sub>
full-loop chemical looping combustion pilot with biomass as fuel: a
hydrodynamic investigation.

Liyan Sun<sup>a</sup>, Enrica Masi<sup>a,*</sup>, Olivier Simonin<sup>a</sup>, Øyvind Langørgen<sup>b</sup>, Inge Saanum<sup>b</sup>, Nils Erland
L. Haugen<sup>b</sup>

**Abstract**

A hydrodynamic model for a full Chemical Looping Combustion (CLC) unit was estab-
lished, and simulations performed using the code NEPTUNE_CFD, which is based on an
Euler-Euler approach. The unit is a 150 kW<sub>th</sub> pilot constructed at SINTEF Energy Research.
Three-dimensional unsteady numerical simulations were carried out for studying the local and
instantaneous behavior inside the system, and its effect on the mean quantities relevant to the
process. Solid volume fraction, mass flow rate and phase velocities were computed and analyzed.
Comparison with experimental results showed that the pressure was globally well predicted.
Two collision models were also investigated. The agitation between neighboring particles was
found to be rather uncorrelated; for this reason, the two collision models led to almost the
same results. This work represents a hydrodynamic assessment of CLC using biomass as fuel.
It allows to provide insight in the flow within the system, with fairly moderate computational
costs.

Keywords: Chemical Looping Combustion; Fluidized Bed Reactor; CFD; Gas-particles flows;

Eulerian approach.

1. Introduction

Chemical looping combustion (CLC) is a novel technology for controlling the CO₂ emission
from combustion processes by separating CO₂ from the combustion products with a very low
energy penalty (Lyngfelt et al. (2001)). It is viewed as an economic method for CO₂ capturing
due to the inherent CO₂ separation. A CLC unit mainly comprises a fuel reactor (FR), an
air reactor (AR), cyclones, loop seals and connecting devices. At pilot scale, CLC reactors are
mostly designed as two circulating fluidized bed (CFB) reactors (e.g., the 100 kWₜₕ unit at
Chalmers (Linderholm et al. (2016)), the 120 kWₜₕ unit in Vienna (Pröll et al. (2009)), and the
1 MWₜₕ unit in Darmstadt (Ströhle et al. (2014))). The oxygen is transported between reactors
by an oxygen carrier (OC) that is a solid phase exchanging mass, momentum and heat with the
gas phase. The oxygen carrier reduced in the fuel reactor is regenerated in the air reactor. The
fuel conversion takes place in the fuel reactor where N₂ will not mix with the gaseous products
(Mattisson et al. (2018)). This makes the CO₂ separation possible by condensation of water,

<footnote><sup>a</sup>Institut de Mécanique des Fluides de Toulouse (IMFT), Université de Toulouse, CNRS, Toulouse, France</footnote>

<footnote><sup>b</sup>Department of Thermal Energy, SINTEF Energy Research, Trondheim, Norway</footnote>

<footnote>*Corresponding author</footnote>

<footnote>Email address: enrica.masi@imft.fr (Enrica Masi)</footnote>
Three-dimensional unsteady numerical simulation of a 150 kW<sub>th</sub>
full-loop chemical looping combustion pilot with biomass as fuel: a
hydrodynamic investigation.

Liyan Sun<sup>a</sup>, Enrica Masi<sup>a,*</sup>, Olivier Simonin<sup>a</sup>, Øyvind Langørgen<sup>b</sup>, Inge Saanum<sup>b</sup>, Nils Erland
L. Haugen<sup>b</sup>

**Abstract**

A hydrodynamic model for a full Chemical Looping Combustion (CLC) unit was estab-
lished, and simulations performed using the code NEPTUNE_CFD, which is based on an
Euler-Euler approach. The unit is a 150 kW<sub>th</sub> pilot constructed at SINTEF Energy Research.
Three-dimensional unsteady numerical simulations were carried out for studying the local and
instantaneous behavior inside the system, and its effect on the mean quantities relevant to the
process. Solid volume fraction, mass flow rate and phase velocities were computed and analyzed.
Comparison with experimental results showed that the pressure was globally well predicted.
Two collision models were also investigated. The agitation between neighboring particles was
found to be rather uncorrelated; for this reason, the two collision models led to almost the
same results. This work represents a hydrodynamic assessment of CLC using biomass as fuel.
It allows to provide insight in the flow within the system, with fairly moderate computational
costs.

Keywords: Chemical Looping Combustion; Fluidized Bed Reactor; CFD; Gas-particles flows;

Eulerian approach.

1. Introduction

Chemical looping combustion (CLC) is a novel technology for controlling the CO₂ emission
from combustion processes by separating CO₂ from the combustion products with a very low
energy penalty (Lyngfelt et al. (2001)). It is viewed as an economic method for CO₂ capturing
due to the inherent CO₂ separation. A CLC unit mainly comprises a fuel reactor (FR), an
air reactor (AR), cyclones, loop seals and connecting devices. At pilot scale, CLC reactors are
mostly designed as two circulating fluidized bed (CFB) reactors (e.g., the 100 kWₜ₡ unit at
Chalmers (Linderholm et al. (2016)), the 120 kWₜ₡ unit in Vienna (Pröll et al. (2009)), and the
1 MWₜ₡ unit in Darmstadt (Ströhle et al. (2014))). The oxygen is transported between reactors
by an oxygen carrier (OC) that is a solid phase exchanging mass, momentum and heat with the
gas phase. The oxygen carrier reduced in the fuel reactor is regenerated in the air reactor. The
fuel conversion takes place in the fuel reactor where N₂ will not mix with the gaseous products
(Mattisson et al. (2018)). This makes the CO₂ separation possible by condensation of water,

<footnote><sup>a</sup>Institut de Mécanique des Fluides de Toulouse (IMFT), Université de Toulouse, CNRS, Toulouse, France</footnote>

<footnote><sup>b</sup>Department of Thermal Energy, SINTEF Energy Research, Trondheim, Norway</footnote>

<footnote>*Corresponding author</footnote>

<footnote>Email address: enrica.masi@imft.fr (Enrica Masi)</footnote>
without resorting to any additional separation methods. This inherent feature of the process is the main strength of the CLC technology.

Since it was first proposed, the CLC process using gaseous fuel has been widely developed and studied in different forms and by different approaches (Li et al. (2017)), including computational fluid dynamics (Wang et al. (2014); Hamidouche et al. (2019)). For solid fuels, the process is more challenging because it involves many additional mechanisms, including pyrolysis, gasification, unsteady feeding of solid fuels and separation of partially converted fuel particles and oxygen carriers (Lyngfelt (2014)). To improve fuel conversion and ensure efficiency in capturing CO₂, the CLC system may indeed require some changes to accommodate the solid fuel, as the addition of an external carbon stripper (Markström et al. (2013); Abad et al. (2020)), or the realization of new reactor designs (Berguerand and Lyngfelt (2008); Kim et al. (2013); Penthor et al. (2016); Haus et al. (2020)). CLC can also resort to the use of CLOU (chemical-looping with oxygen uncoupling) materials, which improve conversion efficiency thanks to their ability to release oxygen in gas phase (Pérez-Vega et al. (2020)). Several experimental works have been carried out to characterize the behavior of the solid-fuel CLC process (Leion et al. (2008); Ströhle et al. (2015)) and inherent reactions (Cao et al. (2006); Siriwardane et al. (2009); Abad et al. (2011)). The power of existing CLC units range from 500 W<sub>th</sub> to 3 MW<sub>th</sub> (Lyngfelt and Linderholm (2017)). Such units use different solid fuels, as coal (Abad et al. (2015)) or biomass (Shen et al. (2009)), together with different oxygen carriers, as ilmenite (Thon et al. (2014)), hematite (Ma et al. (2018)), or manganese ore (Pérez-Astray et al. (2020)), for example. The complexity of the solid-fueled CLC concept makes its modeling and design a real challenge.

With the continuous development of supercomputers and their increasing performance, numerical approaches are becoming more and more powerful for studying industrial applications at the conception stage or processes that need retrofitting. Nowadays, the Computational Fluid Dynamics (CFD) benefits from High Performance Computing (HPC) systems based on massively parallel architectures, which make its use possible even at industrial scale. In particular, unsteady numerical simulations have the advantage to give access to the local and instantaneous fields inside the system. This feature makes the unsteady numerical approach very useful to reproduce the characteristics of a process or salient parts of it, providing complementary information to the experimental research. For this reason, numerical studies have increased significantly in recent years, including on CLC technology. An overview of numerical works on CLC is provided by a recent review of Shao et al. (2021). Only some of them concern solid-fueled CLC systems. The 2D numerical simulation, which has successfully been used over the years for reproducing crucial parts of the CLC (see, e.g., Mahalatkar et al. (2011)) or the entire loop (Su et al. (2015)), is leaving the place to the 3D numerical simulation (see May et al. (2018) as an example) which is more representative of the complex structures of the flow due to the three dimensional nature of its behavior (turbulent conditions, loss of symmetry close to injections, etc.). Concerning full-loop solid-fueled reactive CLC, three-dimensional studies have become available in the literature (Parker (2014); Reinking et al. (2019)) but still few of them assess the numerical results compared to experimental measurements (Chen et al. (2019)).

It is well known that hydrodynamics strongly affects the reactive predictions. The reason is that characteristic reaction times are usually very large compared to those of the flow evolution or momentum transfers in this type of process (in absence of CLOU materials). This is why cold flow models are frequently employed to characterize the CLC behavior, before moving on to ultimate reactive conditions. The increased complexity of the flow when working with solid fuels makes hydrodynamic investigations even more useful in a first stage. Examples of numerical studies based on full-loop cold-flow CLC are the recent works of Wang et al. (2020a) and Wang et al. (2020b), who studied the hydrodynamics of CLC units conceived
without resorting to any additional separation methods. This inherent feature of the process is the main strength of the CLC technology.

Since it was first proposed, the CLC process using gaseous fuel has been widely developed and studied in different forms and by different approaches (Li et al. (2017)), including computational fluid dynamics (Wang et al. (2014); Hamidouche et al. (2019)). For solid fuels, the process is more challenging because it involves many additional mechanisms, including pyrolysis, gasification, unsteady feeding of solid fuels and separation of partially converted fuel particles and oxygen carriers (Lyngfelt (2014)). To improve fuel conversion and ensure efficiency in capturing CO₂, the CLC system may indeed require some changes to accommodate the solid fuel, as the addition of an external carbon stripper (Markström et al. (2013); Abad et al. (2020)), or the realization of new reactor designs (Berguerand and Lyngfelt (2008); Kim et al. (2013); Penthor et al. (2016); Haus et al. (2020)). CLC can also resort to the use of CLOU (chemical-looping with oxygen uncoupling) materials, which improve conversion efficiency thanks to their ability to release oxygen in gas phase (Pérez-Vega et al. (2020)). Several experimental works have been carried out to characterize the behavior of the solid-fuel CLC process (Leion et al. (2008); Ströhle et al. (2015)) and inherent reactions (Cao et al. (2006); Siriwardane et al. (2009); Abad et al. (2011)). The power of existing CLC units range from 500 W<sub>th</sub> to 3 MW<sub>th</sub> (Lyngfelt and Linderholm (2017)). Such units use different solid fuels, as coal (Abad et al. (2015)) or biomass (Shen et al. (2009)), together with different oxygen carriers, as ilmenite (Thon et al. (2014)), hematite (Ma et al. (2018)), or manganese ore (Pérez-Astray et al. (2020)), for example. The complexity of the solid-fueled CLC concept makes its modeling and design a real challenge.

With the continuous development of supercomputers and their increasing performance, numerical approaches are becoming more and more powerful for studying industrial applications at the conception stage or processes that need retrofitting. Nowadays, the Computational Fluid Dynamics (CFD) benefits from High Performance Computing (HPC) systems based on massively parallel architectures, which make its use possible even at industrial scale. In particular, unsteady numerical simulations have the advantage to give access to the local and instantaneous fields inside the system. This feature makes the unsteady numerical approach very useful to reproduce the characteristics of a process or salient parts of it, providing complementary information to the experimental research. For this reason, numerical studies have increased significantly in recent years, including on CLC technology. An overview of numerical works on CLC is provided by a recent review of Shao et al. (2021). Only some of them concern solid-fueled CLC systems. The 2D numerical simulation, which has successfully been used over the years for reproducing crucial parts of the CLC (see, e.g., Mahalatkar et al. (2011)) or the entire loop (Su et al. (2015)), is leaving the place to the 3D numerical simulation (see May et al. (2018) as an example) which is more representative of the complex structures of the flow due to the three dimensional nature of its behavior (turbulent conditions, loss of symmetry close to injections, etc.). Concerning full-loop solid-fueled reactive CLC, three-dimensional studies have become available in the literature (Parker (2014); Reinking et al. (2019)) but still few of them assess the numerical results compared to experimental measurements (Chen et al. (2019)).

It is well known that hydrodynamics strongly affects the reactive predictions. The reason is that characteristic reaction times are usually very large compared to those of the flow evolution or momentum transfers in this type of process (in absence of CLOU materials). This is why cold flow models are frequently employed to characterize the CLC behavior, before moving on to ultimate reactive conditions. The increased complexity of the flow when working with solid fuels makes hydrodynamic investigations even more useful in a first stage. Examples of numerical studies based on full-loop cold-flow CLC are the recent works of Wang et al. (2020a) and Wang et al. (2020b), who studied the hydrodynamics of CLC units conceived
to be used with coal as fuel, with inherent separation or gasification systems. Their studies focused, respectively, on the investigation of the separation (coal from OC) particle efficiency by the high-flux carbon stripper integrated into the process, and on the characterization of the flow in gasifier and reduction reactors, depending on the operating conditions. Results were validated by comparing with previous experimental investigations from the same laboratory. Another example of cold-flow CLC investigation is the three-dimensional numerical simulation carried out by Shao et al. (2020), who analyzed a novel two-stage air reactor and its response under different operating conditions on the whole CLC behavior. Results about the pressure predictions were validated against experimental data.

In the present work, we also explore the hydrodynamics of a full-loop solid-fueled CLC, but comparing pressure predictions with experimental measurements from a hot instead of a cold experimental unit. The hot unit is a 150 $kW_{th}$ pilot operating at SINTEF Energy Research (Trondheim, Norway). The 3D unsteady numerical simulations are performed using an Euler-Euler approach. The latter is considered for its efficiency and low computational costs. In fact, in an Euler-Euler approach, most of the efforts are spent on the development and validation of the modeling to account for additional physical effects (as, for example, particle rotation with friction (Goniva et al. (2012)) or triboelectric charging (Kolehmainen et al. (2018); Montilla et al. (2020))), as well as on the numerical implementation in industrial codes. Depending on the particle characteristics and dimensions at industrial scale, a filtered formulation or heterogeneity models may also be needed (Schneiderbauer and Pirker (2014)). But in the end, the result is an approach that is not excessively time consuming, unless to solve for a distribution of particle sizes, and that has the numerical advantage of treating the separate sets of phase equations analogously, which allows a strong, potentially implicit, coupling between the phases, implying a true mathematical convergence rate with respect to the mesh size and time step.

The alternative Euler-Lagrange particle approaches, such as the Discrete Element Method (CFD-DEM) (Cundall and Strack (1979); Tsuji et al. (1993)), have modeling advantages (easier implementation of additional physical aspects such as polydispersion, particle rotation, particle-particle friction, non-spherical shape, etc.) and numerical advantages (for example non-diffusive Lagrangian numerical schemes leading to less sensitivity to the mesh size), but they are hugely expensive in terms of computational costs already at pilot scale, and definitely unusable at industrial scale. Emerging alternative methods are the Euler-Lagrange approaches using parcels instead of particles (Pirker et al. (2010)), directly accounting for collisions between parcels and at the wall (an overview is given in the review of Di Renzo et al. (2021)), or modeling collisions on a continuum basis (see, e.g., Snider (2001); Cloete et al. (2012)). Unlike a CFD-DEM approach where particles are tracked individually, Euler-Lagrange approaches using parcels need additional assumptions to model mechanisms acting on the particles, which are not taken directly into account (such as solids contacts). These approaches are promising and offer an affordable alternative to the Eulerian models, for example in polydisperse flows when applications require to account for several particle sizes. Depending on the model, they can also provide accuracy improvement, especially in regimes with fine clusters and large-scale crossing (Cloete et al. (2012)). However, efforts still have to be made to reach a degree of maturity equivalent to the Euler-Euler methods. For this reason, an Euler-Euler approach still remains the most reliable and competitive in dense (or moderately dense) regimes; this is why it was considered for this work.

The numerical simulations performed in this study use a non-reactive isothermal model. The model considers fuel and OC conversion by accounting for additional gas injection due to the products from the full fuel conversion and redox reactions. The OC flow behavior inside the reactors and the effect of the coupling of the two reactors are analyzed to improve understanding
to be used with coal as fuel, with inherent separation or gasification systems. Their studies focused, respectively, on the investigation of the separation (coal from OC) particle efficiency by the high-flux carbon stripper integrated into the process, and on the characterization of the flow in gasifier and reduction reactors, depending on the operating conditions. Results were validated by comparing with previous experimental investigations from the same laboratory. Another example of cold-flow CLC investigation is the three-dimensional numerical simulation carried out by Shao et al. (2020), who analyzed a novel two-stage air reactor and its response under different operating conditions on the whole CLC behavior. Results about the pressure predictions were validated against experimental data.

In the present work, we also explore the hydrodynamics of a full-loop solid-fueled CLC, but comparing pressure predictions with experimental measurements from a hot instead of a cold experimental unit. The hot unit is a 150 $kW_{th}$ pilot operating at SINTEF Energy Research (Trondheim, Norway). The 3D unsteady numerical simulations are performed using an Euler-Euler approach. The latter is considered for its efficiency and low computational costs. In fact, in an Euler-Euler approach, most of the efforts are spent on the development and validation of the modeling to account for additional physical effects (as, for example, particle rotation with friction (Goniva et al. (2012)) or triboelectric charging (Kolehmainen et al. (2018); Montilla et al. (2020))), as well as on the numerical implementation in industrial codes. Depending on the particle characteristics and dimensions at industrial scale, a filtered formulation or heterogeneity models may also be needed (Schneiderbauer and Pirker (2014)). But in the end, the result is an approach that is not excessively time consuming, unless to solve for a distribution of particle sizes, and that has the numerical advantage of treating the separate sets of phase equations analogously, which allows a strong, potentially implicit, coupling between the phases, implying a true mathematical convergence rate with respect to the mesh size and time step.

The alternative Euler-Lagrange particle approaches, such as the Discrete Element Method (CFD-DEM) (Cundall and Strack (1979); Tsuji et al. (1993)), have modeling advantages (easier implementation of additional physical aspects such as polydispersion, particle rotation, particle-particle friction, non-spherical shape, etc.) and numerical advantages (for example non-diffusive Lagrangian numerical schemes leading to less sensitivity to the mesh size), but they are hugely expensive in terms of computational costs already at pilot scale, and definitely unusable at industrial scale. Emerging alternative methods are the Euler-Lagrange approaches using parcels instead of particles (Pirker et al. (2010)), directly accounting for collisions between parcels and at the wall (an overview is given in the review of Di Renzo et al. (2021)), or modeling collisions on a continuum basis (see, e.g., Snider (2001); Cloete et al. (2012)). Unlike a CFD-DEM approach where particles are tracked individually, Euler-Lagrange approaches using parcels need additional assumptions to model mechanisms acting on the particles, which are not taken directly into account (such as solids contacts). These approaches are promising and offer an affordable alternative to the Eulerian models, for example in polydisperse flows when applications require to account for several particle sizes. Depending on the model, they can also provide accuracy improvement, especially in regimes with fine clusters and large-scale crossing (Cloete et al. (2012)). However, efforts still have to be made to reach a degree of maturity equivalent to the Euler-Euler methods. For this reason, an Euler-Euler approach still remains the most reliable and competitive in dense (or moderately dense) regimes; this is why it was considered for this work.

The numerical simulations performed in this study use a non-reactive isothermal model. The model considers fuel and OC conversion by accounting for additional gas injection due to the products from the full fuel conversion and redox reactions. The OC flow behavior inside the reactors and the effect of the coupling of the two reactors are analyzed to improve understanding
of the CLC system. The results obtained from these numerical simulations should help in the design and operation of CLC units.

## 2. Numerical approach and mathematical modeling

Unsteady 3D numerical simulations of the CLC unit are carried out using the N-Euler approach for gas-solid turbulent flows implemented in NEPTUNE_CFD by IMFT (Institut de Mécanique des fluides de Toulouse), in collaboration with EDF (Electricité de France) R&D (Hamidouche et al. (2018), Neau et al. (2020)). NEPTUNE_CFD is a multiphase CFD code developed in the framework of the NEPTUNE project, financially supported by EDF, CEA (Commissariat à l'Énergie Atomique), IRSN (Institut de Radioprotection et de Sûreté Nucléaire) and Framatome. The code solves the coupled partial differential equations by a finite-volume approach using an adaptive time step determined by a CFL (Courant-Friedrichs-Lewy) criterion for each phase. The solver is based on a cell-center type finite volume method and an elliptic fractional time-step method. The latter relies on an alpha-pressure cycle (alpha stands for phase volume fraction), which is an iterative method to ensure mass and energy conservation. First, at the beginning of each substep, the velocity is predicted for each phase without accounting for volume fraction and pressure variations in time, while accounting for inter-phase coupling by using an implicit formulation during the sub-step iterations. Then, mass and energy equations are integrated enforcing conservativity, and the velocity is corrected by accounting for volume fraction and pressure variations. Then, the pressure is computed by solving an elliptic equation and the velocities are corrected with respect to the pressure time increment. Convergence criteria of the alpha-pressure cycling is based on the condition of volume conservation of the mixture (EDF R&D (2017)). The code is characterized by a calculation of co-localized gradients with reconstruction methods and a distributed-memory parallelism by domain decomposition (MPI parallelization). It uses unstructured meshes with all cell types and connections. Further details about the numerical code can be found in Neau et al. (2020).

In the present work, the multiphase N-Euler approach implemented in NEPTUNE_CFD has been used to model the evolution of both the gas and solid phase under isothermal conditions. In this section, the corresponding mathematical modeling is presented. More details about the approach may be found in the work of Simonin (2000).

In the current study, a non-reactive isothermal hydrodynamic investigation is carried out. On this basis, the mass balance equations are written as follows:

$$ \frac{\partial(\alpha_g \rho_g)}{\partial t} + \frac{\partial(\alpha_g \rho_g U_{g,j})}{\partial x_j} = 0, \qquad (1) $$

$$ \frac{\partial(\alpha_s \rho_s)}{\partial t} + \frac{\partial(\alpha_s \rho_s U_{s,j})}{\partial x_j} = 0, \qquad (2) $$

where $\rho$, $\alpha$ and $U$ are mean density, volume fraction and velocity, respectively. The subscripts $g$ represents the gas phase while $s$ stands for the solid phase. Since reactions are not taken into consideration in the current work, source terms related to the mass transfer are set to zero. The momentum equations are given by:

$$ \alpha_g \rho_g \left( \frac{\partial U_{g,i}}{\partial t} + U_{g,j} \frac{\partial U_{g,i}}{\partial x_j} \right) = -\alpha_g \frac{\partial P_g}{\partial x_i} + \alpha_g \rho_g g_i + I_{s \to g,i} + \frac{\partial \sum_{g,ij}}{\partial x_j}, \qquad (3) $$

$$ \alpha_s \rho_s \left( \frac{\partial U_{s,i}}{\partial t} + U_{s,j} \frac{\partial U_{s,i}}{\partial x_j} \right) = -\alpha_s \frac{\partial P_s}{\partial x_i} + \alpha_s \rho_s g_i + I_{g \to s,i} + \frac{\partial \sum_{s,ij}}{\partial x_j}. \qquad (4) $$
of the CLC system. The results obtained from these numerical simulations should help in the design and operation of CLC units.

## 2. Numerical approach and mathematical modeling

Unsteady 3D numerical simulations of the CLC unit are carried out using the N-Euler approach for gas-solid turbulent flows implemented in NEPTUNE_CFD by IMFT (Institut de Mécanique des fluides de Toulouse), in collaboration with EDF (Electricité de France) R&D (Hamidouche et al. (2018), Neau et al. (2020)). NEPTUNE_CFD is a multiphase CFD code developed in the framework of the NEPTUNE project, financially supported by EDF, CEA (Commissariat à l'Énergie Atomique), IRSN (Institut de Radioprotection et de Sûreté Nucléaire) and Framatome. The code solves the coupled partial differential equations by a finite-volume approach using an adaptive time step determined by a CFL (Courant-Friedrichs-Lewy) criterion for each phase. The solver is based on a cell-center type finite volume method and an elliptic fractional time-step method. The latter relies on an alpha-pressure cycle (alpha stands for phase volume fraction), which is an iterative method to ensure mass and energy conservation. First, at the beginning of each substep, the velocity is predicted for each phase without accounting for volume fraction and pressure variations in time, while accounting for inter-phase coupling by using an implicit formulation during the sub-step iterations. Then, mass and energy equations are integrated enforcing conservativity, and the velocity is corrected by accounting for volume fraction and pressure variations. Then, the pressure is computed by solving an elliptic equation and the velocities are corrected with respect to the pressure time increment. Convergence criteria of the alpha-pressure cycling is based on the condition of volume conservation of the mixture (EDF R&D (2017)). The code is characterized by a calculation of co-localized gradients with reconstruction methods and a distributed-memory parallelism by domain decomposition (MPI parallelization). It uses unstructured meshes with all cell types and connections. Further details about the numerical code can be found in Neau et al. (2020).

In the present work, the multiphase N-Euler approach implemented in NEPTUNE_CFD has been used to model the evolution of both the gas and solid phase under isothermal conditions. In this section, the corresponding mathematical modeling is presented. More details about the approach may be found in the work of Simonin (2000).

In the current study, a non-reactive isothermal hydrodynamic investigation is carried out. On this basis, the mass balance equations are written as follows:

$$ \frac{\partial(\alpha_g \rho_g)}{\partial t} + \frac{\partial(\alpha_g \rho_g U_{g,j})}{\partial x_j} = 0, \qquad (1) $$

$$ \frac{\partial(\alpha_s \rho_s)}{\partial t} + \frac{\partial(\alpha_s \rho_s U_{s,j})}{\partial x_j} = 0, \qquad (2) $$

where $\rho$, $\alpha$ and $U$ are mean density, volume fraction and velocity, respectively. The subscripts $g$ represents the gas phase while $s$ stands for the solid phase. Since reactions are not taken into consideration in the current work, source terms related to the mass transfer are set to zero. The momentum equations are given by:

$$ \alpha_g \rho_g \left( \frac{\partial U_{g,i}}{\partial t} + U_{g,j} \frac{\partial U_{g,i}}{\partial x_j} \right) = -\alpha_g \frac{\partial P_g}{\partial x_i} + \alpha_g \rho_g g_i + I_{s \to g,i} + \frac{\partial \sum_{g,ij}}{\partial x_j}, \qquad (3) $$

$$ \alpha_s \rho_s \left( \frac{\partial U_{s,i}}{\partial t} + U_{s,j} \frac{\partial U_{s,i}}{\partial x_j} \right) = -\alpha_s \frac{\partial P_s}{\partial x_i} + \alpha_s \rho_s g_i + I_{g \to s,i} + \frac{\partial \sum_{s,ij}}{\partial x_j}. \qquad (4) $$
In the above equations, $P_g$ is the gas pressure and $I_{g→s} (= -I_{s→g})$ is the mean gas to solid interphase momentum transfer after subtracting the mean gas pressure gradient contribution (Archimedes' force). It will be detailed later. $\sum_{,ij}$ are stress tensors defined as:

$$ \sum_{g,ij} = -\alpha_g \rho_g \langle u''_{g,i} u''_{g,j} \rangle_g + \Theta_{g,ij} = -\alpha_g \rho_g R_{g,ij} + \Theta_{g,ij}, \quad (5) $$

$$ \sum_{s,ij} = -\alpha_s \rho_s \langle u''_{s,i} u''_{s,j} \rangle_s + \Theta_{s,ij} + \phi_{s,ij} = -\alpha_s \rho_s R_{s,ij} + \Theta_{s,ij} + \phi_{s,ij}, \quad (6) $$

where $u''_{,i} = u_{,i} - U_{,i}$.

For the gas phase, $R_{g,ij}$ and $\Theta_{g,ij}$ represent the turbulent-Reynolds and viscous stress tensors. They are written as:

$$ R_{g,ij} = -\nu_g^t \left( \frac{\partial U_{g,i}}{\partial x_j} + \frac{\partial U_{g,j}}{\partial x_i} \right) + \frac{2}{3} \delta_{ij} \left( k + \nu_g^t \frac{\partial U_{g,m}}{\partial x_m} \right), \quad (7) $$

$$ \Theta_{g,ij} = \alpha_g \mu_g \left( \frac{\partial U_{g,i}}{\partial x_j} + \frac{\partial U_{g,j}}{\partial x_i} - \frac{2}{3} \frac{\partial U_{g,m}}{\partial x_m} \delta_{ij} \right), \quad (8) $$

where $\delta_{ij}$ is the Kronecker delta. $k$ and $\mu_g$ are turbulent kinetic energy and laminar dynamic viscosity, respectively. $\nu_g^t$ is the turbulent kinematic viscosity written as (Vermorel et al. (2003)):

$$ \nu_g^t = \frac{2}{3} k \tau_g^t \left[ 1 + C_{12} \frac{\alpha_s \rho_s \tau_{gs}^t}{\alpha_g \rho_g \tau_{gs}^F} \left( 1 - \frac{q_{gs}}{2k} \right) \right]^{-1}, \quad (9) $$

where the constant $C_{12}=0.34$. The quantity $q_{gs}$ is the fluid-particle velocity covariance and it will be presented later with the solid phase. $\tau_{gs}^t$ and $\tau_{gs}^F$ are timescales related to the interaction between the gas and the solid phases. The eddy-particle interaction time is the characteristic time for the gas turbulence seen by the particles (Simonin et al. (1993)):

$$ \tau_{gs}^t = \frac{\tau_g^t}{\sigma_k} \left( 1 + C_\beta \frac{V_{r,i} V_{r,j}}{\frac{2}{3} k} \right)^{-1/2}, \quad (10) $$

and $\tau_{gs}^F$ is the mean particle relaxation time (detailed later). The fluid turbulent timescale is defined as $\tau_g^t = C_\mu \frac{3k}{2\epsilon}$. A $k-\epsilon$ model is adopted for closing the above equations. According to this model (Vermorel et al. (2003)), the transport equations for the gas turbulent kinetic energy and dissipation rate are written as:

$$ \alpha_g \rho_g \left( \frac{\partial k}{\partial t} + U_{g,j} \frac{\partial k}{\partial x_j} \right) = \frac{\partial}{\partial x_j} \left( \alpha_g \rho_g \frac{\nu_g^t}{\sigma_k} \frac{\partial k}{\partial x_j} \right) - \alpha_g \rho_g R_{g,ij} \frac{\partial U_{g,i}}{\partial x_j} - \alpha_g \rho_g \epsilon + \Pi_{s \to g}^k, \quad (11) $$

$$ \alpha_g \rho_g \left( \frac{\partial \epsilon}{\partial t} + U_{g,j} \frac{\partial \epsilon}{\partial x_j} \right) = \frac{\partial}{\partial x_j} \left( \alpha_g \rho_g \frac{\nu_g^t}{\sigma_\epsilon} \frac{\partial \epsilon}{\partial x_j} \right) - \alpha_g \rho_g C_{\epsilon 1} \frac{\epsilon}{k} R_{g,ij} \frac{\partial U_{g,i}}{\partial x_j} - \alpha_g \rho_g C_{\epsilon 2} \frac{\epsilon^2}{k} + \Pi_{s \to g}^\epsilon, \quad (12) $$

where $\Pi_{s \to g}^k$ and $\Pi_{s \to g}^\epsilon$ account for the effect of the solid phase on the gas turbulence. Assuming that particle size is comparable or less than the Kolmogorov length scale, the interphase terms are given by

$$ \Pi_{s \to g}^k = \frac{\alpha_s \rho_s}{\tau_{gs}^F} (-2k + q_{gs} + V_{d,i} V_{r,i}), \quad (13) $$

and

$$ \Pi_{s \to g}^\epsilon = C_{\epsilon 3} \frac{\epsilon}{k} \Pi_{s \to g}^k, \quad (14) $$
In the above equations, $P_g$ is the gas pressure and $I_{g→s} (= -I_{s→g})$ is the mean gas to solid interphase momentum transfer after subtracting the mean gas pressure gradient contribution (Archimedes' force). It will be detailed later. $\sum_{,ij}$ are stress tensors defined as:

$$ \sum_{g,ij} = -\alpha_g \rho_g \langle u''_{g,i} u''_{g,j} \rangle_g + \Theta_{g,ij} = -\alpha_g \rho_g R_{g,ij} + \Theta_{g,ij}, \quad (5) $$

$$ \sum_{s,ij} = -\alpha_s \rho_s \langle u''_{s,i} u''_{s,j} \rangle_s + \Theta_{s,ij} + \phi_{s,ij} = -\alpha_s \rho_s R_{s,ij} + \Theta_{s,ij} + \phi_{s,ij}, \quad (6) $$

where $u''_{,i} = u_{,i} - U_{,i}$.

For the gas phase, $R_{g,ij}$ and $\Theta_{g,ij}$ represent the turbulent-Reynolds and viscous stress tensors. They are written as:

$$ R_{g,ij} = -\nu_g^t \left( \frac{\partial U_{g,i}}{\partial x_j} + \frac{\partial U_{g,j}}{\partial x_i} \right) + \frac{2}{3} \delta_{ij} \left( k + \nu_g^t \frac{\partial U_{g,m}}{\partial x_m} \right), \quad (7) $$

$$ \Theta_{g,ij} = \alpha_g \mu_g \left( \frac{\partial U_{g,i}}{\partial x_j} + \frac{\partial U_{g,j}}{\partial x_i} - \frac{2}{3} \frac{\partial U_{g,m}}{\partial x_m} \delta_{ij} \right), \quad (8) $$

where $\delta_{ij}$ is the Kronecker delta. $k$ and $\mu_g$ are turbulent kinetic energy and laminar dynamic viscosity, respectively. $\nu_g^t$ is the turbulent kinematic viscosity written as (Vermorel et al. (2003)):

$$ \nu_g^t = \frac{2}{3} k \tau_g^t \left[ 1 + C_{12} \frac{\alpha_s \rho_s \tau_{gs}^t}{\alpha_g \rho_g \tau_{gs}^F} \left( 1 - \frac{q_{gs}}{2k} \right) \right]^{-1}, \quad (9) $$

where the constant $C_{12}=0.34$. The quantity $q_{gs}$ is the fluid-particle velocity covariance and it will be presented later with the solid phase. $\tau_{gs}^t$ and $\tau_{gs}^F$ are timescales related to the interaction between the gas and the solid phases. The eddy-particle interaction time is the characteristic time for the gas turbulence seen by the particles (Simonin et al. (1993)):

$$ \tau_{gs}^t = \frac{\tau_g^t}{\sigma_k} \left( 1 + C_\beta \frac{V_{r,i} V_{r,j}}{\frac{2}{3} k} \right)^{-1/2}, \quad (10) $$

and $\tau_{gs}^F$ is the mean particle relaxation time (detailed later). The fluid turbulent timescale is defined as $\tau_g^t = C_\mu \frac{3k}{2\epsilon}$. A $k-\epsilon$ model is adopted for closing the above equations. According to this model (Vermorel et al. (2003)), the transport equations for the gas turbulent kinetic energy and dissipation rate are written as:

$$ \alpha_g \rho_g \left( \frac{\partial k}{\partial t} + U_{g,j} \frac{\partial k}{\partial x_j} \right) = \frac{\partial}{\partial x_j} \left( \alpha_g \rho_g \frac{\nu_g^t}{\sigma_k} \frac{\partial k}{\partial x_j} \right) - \alpha_g \rho_g R_{g,ij} \frac{\partial U_{g,i}}{\partial x_j} - \alpha_g \rho_g \epsilon + \Pi_{s \to g}^k, \quad (11) $$

$$ \alpha_g \rho_g \left( \frac{\partial \epsilon}{\partial t} + U_{g,j} \frac{\partial \epsilon}{\partial x_j} \right) = \frac{\partial}{\partial x_j} \left( \alpha_g \rho_g \frac{\nu_g^t}{\sigma_\epsilon} \frac{\partial \epsilon}{\partial x_j} \right) - \alpha_g \rho_g C_{\epsilon 1} \frac{\epsilon}{k} R_{g,ij} \frac{\partial U_{g,i}}{\partial x_j} - \alpha_g \rho_g C_{\epsilon 2} \frac{\epsilon^2}{k} + \Pi_{s \to g}^\epsilon, \quad (12) $$

where $\Pi_{s \to g}^k$ and $\Pi_{s \to g}^\epsilon$ account for the effect of the solid phase on the gas turbulence. Assuming that particle size is comparable or less than the Kolmogorov length scale, the interphase terms are given by

$$ \Pi_{s \to g}^k = \frac{\alpha_s \rho_s}{\tau_{gs}^F} (-2k + q_{gs} + V_{d,i} V_{r,i}), \quad (13) $$

and

$$ \Pi_{s \to g}^\epsilon = C_{\epsilon 3} \frac{\epsilon}{k} \Pi_{s \to g}^k, \quad (14) $$
173 using the relative velocity

$$V_{r,i} = (U_{s,i} - U_{g,i}) - V_{d,i}, \qquad (15)$$

174 and accounting for the turbulent drift velocity (Simonin et al. (1993))

$$V_{d,i} = -D_{gs}^t \left( \frac{1}{\alpha_s} \frac{\partial \alpha_s}{\partial x_i} - \frac{1}{\alpha_g} \frac{\partial \alpha_g}{\partial x_i} \right). \qquad (16)$$

175 The constants involved in the $k-\epsilon$ model are $C_\mu =0.09$, $C_{\epsilon1} =1.44$, $C_{\epsilon2} =1.92$, $C_{\epsilon3} =1.2$, $\sigma_k =1.0$ and $\sigma_\epsilon =1.3$.

176 For the solid phase, the effective stress tensor (Eq. 6) comprises a kinetic part, $R_{s,ij}$, which is dominant in dilute flow, and a collisional part, $\Theta_{s,ij}$, which is dominant in dense flow. 177 Also the frictional part, $\phi_{s,ij}$, contributes to the effective stress tensor in zones with very high 178 concentration and long solid-solid contact. The kinetic and collisional contributions of the 179 effective particle stress tensor are written, respectively, as (Boelle et al. (1995); Gobin et al. 180 (2003); Jenkins and Richman (1986); Simonin (2000)) 181 182

$$R_{s,ij} = -\nu_s^{kin} \left( \frac{\partial U_{s,i}}{\partial x_j} + \frac{\partial U_{s,j}}{\partial x_i} \right) + \frac{2}{3} \delta_{ij} \left( q_s^2 + \nu_s^{kin} \frac{\partial U_{s,m}}{\partial x_m} \right), \qquad (17)$$

183

$$\Theta_{s,ij} = - \left[ \frac{2}{3} \alpha_s \rho_s q_s^2 2 \alpha_s g_0 (1 + e_c) - \Lambda_s \frac{\partial U_{s,m}}{\partial x_m} \right] \delta_{ij} + \alpha_s \rho_s \nu_s^{col} \left( \frac{\partial U_{s,i}}{\partial x_j} + \frac{\partial U_{s,j}}{\partial x_i} - \frac{2}{3} \frac{\partial U_{s,m}}{\partial x_m} \delta_{ij} \right). \qquad (18)$$

184 In the above equations, $e_c$ is the normal restitution coefficient, $g_0$ is the radial distribution 185 function, $\nu_s^{kin}$ is the particle kinetic viscosity, and $\nu_s^{col}$ represents the particle collisional viscosity:

186

$$\nu_s^{kin} = \left[ \nu_{gs}^t + \frac{\tau_{gs}^F}{2} \frac{2}{3} q_s^2 (1 + \alpha_s g_0 \Phi_c) \right] \left( 1 + \frac{\tau_{gs}^F \sigma_c}{2 \tau_s^c} \right)^{-1} \qquad (19)$$

187

$$\nu_s^{col} = \frac{4}{5} \alpha_s g_0 (1 + e_c) \left( \nu_s^{kin} + d_s \sqrt{\frac{2 q_s^2}{3 \pi}} \right). \qquad (20)$$

188 $\Lambda_s$, in Eq. 18, is defined as:

$$\Lambda_s = \alpha_s \rho_s \frac{4}{3} \alpha_s g_0 (1 + e_c) d_s \sqrt{\frac{2 q_s^2}{3 \pi}} \qquad (21)$$

189 The particle fluctuant kinetic energy, $q_s^2$, is defined as $q_s^2 = \langle u_{s,i}'' u_{s,i}'' \rangle_s / 2$. The transport equation 190 of $q_s^2$ is

$$\alpha_s \rho_s \left( \frac{\partial q_s^2}{\partial t} + U_{s,j} \frac{\partial q_s^2}{\partial x_j} \right) = \frac{\partial}{\partial x_j} \left( \alpha_s \rho_s \kappa_s^{eff} \frac{\partial q_s^2}{\partial x_j} \right) + \sum_{s,ij} \frac{\partial U_{s,i}}{\partial x_j} - \alpha_s \rho_s \epsilon_s + \Pi_{qs}, \qquad (22)$$

191 where $\kappa_s^{eff}$ is the particle effective diffusivity coefficient, $\kappa_s^{eff} = \kappa_s^{kin} + \kappa_s^{col}$, formed by the 192 following contributions:

193

$$\kappa_s^{kin} = \left[ \frac{1}{3} \tau_{gs}^t q_{gs} + \frac{5}{9} \tau_{gs}^F \frac{2}{3} q_s^2 (1 + \alpha_s g_0 \varphi_c) \right] \left( 1 + \frac{5}{9} \tau_{gs}^F \frac{\xi_c}{\tau_s^c} \right)^{-1} \qquad (23)$$

$$\kappa_s^{col} = \alpha_s g_0 (1 + e_c) \left[ \frac{6}{5} \kappa_s^{kin} + \frac{4}{3} d_s \sqrt{\frac{2 q_s^2}{3 \pi}} \right]. \qquad (24)$$
173 using the relative velocity

$$V_{r,i} = (U_{s,i} - U_{g,i}) - V_{d,i}, \qquad (15)$$

174 and accounting for the turbulent drift velocity (Simonin et al. (1993))

$$V_{d,i} = -D_{gs}^t \left( \frac{1}{\alpha_s} \frac{\partial \alpha_s}{\partial x_i} - \frac{1}{\alpha_g} \frac{\partial \alpha_g}{\partial x_i} \right). \qquad (16)$$

175 The constants involved in the $k-\epsilon$ model are $C_\mu =0.09$, $C_{\epsilon1} =1.44$, $C_{\epsilon2} =1.92$, $C_{\epsilon3} =1.2$, $\sigma_k =1.0$ and $\sigma_\epsilon =1.3$.

176 For the solid phase, the effective stress tensor (Eq. 6) comprises a kinetic part, $R_{s,ij}$, which is dominant in dilute flow, and a collisional part, $\Theta_{s,ij}$, which is dominant in dense flow. 177 Also the frictional part, $\phi_{s,ij}$, contributes to the effective stress tensor in zones with very high 178 concentration and long solid-solid contact. The kinetic and collisional contributions of the 179 effective particle stress tensor are written, respectively, as (Boelle et al. (1995); Gobin et al. 180 (2003); Jenkins and Richman (1986); Simonin (2000)) 181 182

$$R_{s,ij} = -\nu_s^{kin} \left( \frac{\partial U_{s,i}}{\partial x_j} + \frac{\partial U_{s,j}}{\partial x_i} \right) + \frac{2}{3} \delta_{ij} \left( q_s^2 + \nu_s^{kin} \frac{\partial U_{s,m}}{\partial x_m} \right), \qquad (17)$$

183

$$\Theta_{s,ij} = - \left[ \frac{2}{3} \alpha_s \rho_s q_s^2 2 \alpha_s g_0 (1 + e_c) - \Lambda_s \frac{\partial U_{s,m}}{\partial x_m} \right] \delta_{ij} + \alpha_s \rho_s \nu_s^{col} \left( \frac{\partial U_{s,i}}{\partial x_j} + \frac{\partial U_{s,j}}{\partial x_i} - \frac{2}{3} \frac{\partial U_{s,m}}{\partial x_m} \delta_{ij} \right). \qquad (18)$$

184 In the above equations, $e_c$ is the normal restitution coefficient, $g_0$ is the radial distribution 185 function, $\nu_s^{kin}$ is the particle kinetic viscosity, and $\nu_s^{col}$ represents the particle collisional viscosity:

186

$$\nu_s^{kin} = \left[ \nu_{gs}^t + \frac{\tau_{gs}^F}{2} \frac{2}{3} q_s^2 (1 + \alpha_s g_0 \Phi_c) \right] \left( 1 + \frac{\tau_{gs}^F \sigma_c}{2 \tau_s^c} \right)^{-1} \qquad (19)$$

187

$$\nu_s^{col} = \frac{4}{5} \alpha_s g_0 (1 + e_c) \left( \nu_s^{kin} + d_s \sqrt{\frac{2 q_s^2}{3 \pi}} \right). \qquad (20)$$

188 $\Lambda_s$, in Eq. 18, is defined as:

$$\Lambda_s = \alpha_s \rho_s \frac{4}{3} \alpha_s g_0 (1 + e_c) d_s \sqrt{\frac{2 q_s^2}{3 \pi}} \qquad (21)$$

189 The particle fluctuant kinetic energy, $q_s^2$, is defined as $q_s^2 = \langle u_{s,i}'' u_{s,i}'' \rangle_s / 2$. The transport equation 190 of $q_s^2$ is

$$\alpha_s \rho_s \left( \frac{\partial q_s^2}{\partial t} + U_{s,j} \frac{\partial q_s^2}{\partial x_j} \right) = \frac{\partial}{\partial x_j} \left( \alpha_s \rho_s \kappa_s^{eff} \frac{\partial q_s^2}{\partial x_j} \right) + \sum_{s,ij} \frac{\partial U_{s,i}}{\partial x_j} - \alpha_s \rho_s \varepsilon_s + \Pi_{qs}, \qquad (22)$$

191 where $\kappa_s^{eff}$ is the particle effective diffusivity coefficient, $\kappa_s^{eff} = \kappa_s^{kin} + \kappa_s^{col}$, formed by the 192 following contributions:

$$\kappa_s^{kin} = \left[ \frac{1}{3} \tau_{gs}^t q_{gs} + \frac{5}{9} \tau_{gs}^F \frac{2}{3} q_s^2 (1 + \alpha_s g_0 \varphi_c) \right] \left( 1 + \frac{5}{9} \tau_{gs}^F \frac{\xi_c}{\tau_s^c} \right)^{-1} \qquad (23)$$

193

$$\kappa_s^{col} = \alpha_s g_0 (1 + e_c) \left[ \frac{6}{5} \kappa_s^{kin} + \frac{4}{3} d_s \sqrt{\frac{2 q_s^2}{3 \pi}} \right]. \qquad (24)$$
$$ \varphi_c = \frac{3}{5} (1 + e_c)^2 (2e_c - 1) \qquad (25) $$

and

$$ \xi_c = \frac{(1 + e_c)(49 - 33e_c)}{100} \qquad (26) $$

$\varepsilon_s$, in Eq. 22, is the particle kinetic energy dissipation rate due to the inelastic collisions (Simonin et al. (2002)):

$$ \varepsilon_s = \frac{1}{3}(1 - e_c^2) \frac{\delta q_s^2}{\tau_s^c}. \qquad (27) $$

where $\delta q_s^2$ represents the uncorrelated part of the random particle kinetic energy, also named granular temperature ($\Theta_s = 2/3 \delta q_s^2$) (Fox (2014)). $\Pi_{qs}$, in Eq. (22), is the interphase turbulent kinetic energy transfer rate and is written as

$$ \Pi_{qs} = -\alpha_s \rho_s \frac{1}{\tau_{gs}^F} (2q_s^2 - q_{gs}), \qquad (28) $$

where $q_{gs} = \langle u''_{g,i} u''_{s,i} \rangle_s$ is the fluid-particle velocity covariance, which is solved by the following transport equation (Simonin (2000))

$$ \alpha_s \rho_s \left( \frac{\partial q_{gs}}{\partial t} + \partial U_{s,j} \frac{\partial q_{gs}}{\partial x_j} \right) = \frac{\partial}{\partial x_j} \left( \alpha_s \rho_s \frac{\nu_{gs}^t}{\sigma_k} \frac{\partial q_{gs}}{\partial x_j} \right) - \alpha_s \rho_s \varepsilon_{gs} + \Pi_{qgs} - \alpha_s \rho_s \left[ \langle u''_{g,i} u''_{s,j} \rangle_s \frac{\partial U_{s,i}}{\partial x_j} + \langle u''_{g,j} u''_{s,i} \rangle_s \frac{\partial U_{g,i}}{\partial x_j} \right], \qquad (29) $$

where $\varepsilon_{gs}$ is the fluid-particle covariance dissipation rate due to viscous dissipation and crossing trajectory effects, which is modeled as:

$$ \varepsilon_{gs} = \frac{q_{gs}}{\tau_{gs}^t}. \qquad (30) $$

The interphase interaction term, $\Pi_{qgs}$, is written as:

$$ \Pi_{qgs} = -\alpha_s \rho_s \frac{1}{\tau_{gs}^F} \left[ (q_{gs} - 2k) + \frac{\alpha_s \rho_s}{\alpha_g \rho_g} (q_{gs} - 2\tilde{q}_s^2) \right]. \qquad (31) $$

where $\tilde{q}_s^2$ is the correlated part of the random part kinetic energy defined later. The first contribution on the right-hand side of the above equation, proportional to $(q_{gs} - 2k)$, represents the effect of particle entrainment by gas turbulence and is dominant, and generally positive, in dilute flows ($2k > q_{gs}$). The second term, proportional to $(q_{gs} - 2\tilde{q}_s^2)$, represents the effect of two-way coupling by particle agitation and is dominant and generally negative ($2\tilde{q}_s^2 < q_{gs}$), in high solid mass loaded flows.

The frictional tensor, $\phi_{s,ij}$, in Eq. (6), is defined according to the frictional model (Bennani et al. (2017)):

$$ \phi_{s,ij} = 2\mu_s^{fr} D_{s,ij} - P_s^{fr} \delta_{ij}. \qquad (32) $$

$D_{s,ij}$ is the particle shear tensor written as:

$$ D_{s,ij} = \frac{1}{2} \left[ \frac{\partial U_{s,i}}{\partial x_j} + \frac{U_{s,j}}{\partial x_i} - \frac{2}{3} \frac{\partial U_k}{\partial x_k} \delta_{ij} \right]. \qquad (33) $$

The frictional pressure, $P_s^{fr}$, is modeled according to Johnson and Jackson (1987) and Johnson et al. (1990):

$$ P_s^{fr} = \begin{cases} F r \frac{(\alpha_s - \alpha_s^{min})^r}{(\alpha_s^{max} - \alpha_s)^s}; & \alpha_s > \alpha_s^{min} \\ 0; & \text{else} \end{cases} \qquad (34) $$
$$ \varphi_c = \frac{3}{5} (1 + e_c)^2 (2e_c - 1) \qquad (25) $$

and

$$ \xi_c = \frac{(1 + e_c)(49 - 33e_c)}{100} \qquad (26) $$

$\varepsilon_s$, in Eq. 22, is the particle kinetic energy dissipation rate due to the inelastic collisions (Simonin et al. (2002)):

$$ \varepsilon_s = \frac{1}{3}(1 - e_c^2) \frac{\delta q_s^2}{\tau_s^c}. \qquad (27) $$

where $\delta q_s^2$ represents the uncorrelated part of the random particle kinetic energy, also named granular temperature ($\Theta_s = 2/3 \delta q_s^2$) (Fox (2014)). $\Pi_{qs}$, in Eq. (22), is the interphase turbulent kinetic energy transfer rate and is written as

$$ \Pi_{qs} = -\alpha_s \rho_s \frac{1}{\tau_{gs}^F} (2q_s^2 - q_{gs}), \qquad (28) $$

where $q_{gs} = \langle u''_{g,i} u''_{s,i} \rangle_s$ is the fluid-particle velocity covariance, which is solved by the following transport equation (Simonin (2000))

$$ \alpha_s \rho_s \left( \frac{\partial q_{gs}}{\partial t} + \partial U_{s,j} \frac{\partial q_{gs}}{\partial x_j} \right) = \frac{\partial}{\partial x_j} \left( \alpha_s \rho_s \frac{\nu_{gs}^t}{\sigma_k} \frac{\partial q_{gs}}{\partial x_j} \right) - \alpha_s \rho_s \varepsilon_{gs} + \Pi_{qgs} - \alpha_s \rho_s \left[ \langle u''_{g,i} u''_{s,j} \rangle_s \frac{\partial U_{s,i}}{\partial x_j} + \langle u''_{g,j} u''_{s,i} \rangle_s \frac{\partial U_{g,i}}{\partial x_j} \right], \qquad (29) $$

where $\varepsilon_{gs}$ is the fluid-particle covariance dissipation rate due to viscous dissipation and crossing trajectory effects, which is modeled as:

$$ \varepsilon_{gs} = \frac{q_{gs}}{\tau_{gs}^t}. \qquad (30) $$

The interphase interaction term, $\Pi_{qgs}$, is written as:

$$ \Pi_{qgs} = -\alpha_s \rho_s \frac{1}{\tau_{gs}^F} \left[ (q_{gs} - 2k) + \frac{\alpha_s \rho_s}{\alpha_g \rho_g} (q_{gs} - 2\tilde{q}_s^2) \right]. \qquad (31) $$

where $\tilde{q}_s^2$ is the correlated part of the random part kinetic energy defined later. The first contribution on the right-hand side of the above equation, proportional to $(q_{gs} - 2k)$, represents the effect of particle entrainment by gas turbulence and is dominant, and generally positive, in dilute flows ($2k > q_{gs}$). The second term, proportional to $(q_{gs} - 2\tilde{q}_s^2)$, represents the effect of two-way coupling by particle agitation and is dominant and generally negative ($2\tilde{q}_s^2 < q_{gs}$), in high solid mass loaded flows.

The frictional tensor, $\phi_{s,ij}$, in Eq. (6), is defined according to the frictional model (Bennani et al. (2017)):

$$ \phi_{s,ij} = 2\mu_s^{fr} D_{s,ij} - P_s^{fr} \delta_{ij}. \qquad (32) $$

$D_{s,ij}$ is the particle shear tensor written as:

$$ D_{s,ij} = \frac{1}{2} \left[ \frac{\partial U_{s,i}}{\partial x_j} + \frac{U_{s,j}}{\partial x_i} - \frac{2}{3} \frac{\partial U_k}{\partial x_k} \delta_{ij} \right]. \qquad (33) $$

The frictional pressure, $P_s^{fr}$, is modeled according to Johnson and Jackson (1987) and Johnson et al. (1990):

$$ P_s^{fr} = \begin{cases} Fr \frac{(\alpha_s - \alpha_s^{min})^r}{(\alpha_s^{max} - \alpha_s)^s}; & \alpha_s > \alpha_s^{min} \\ 0; & \text{else} \end{cases} \qquad (34) $$
where $Fr = 0.05$, $r = 2$ and $s = 5$ are model parameters, which may be varied depending on the types of particles. In this work, $\alpha_s^{min}$ is set equal to 0.55, which is an appropriate value for spherical particles. The frictional viscosity is modeled as follows (Srivastava and Sundaresan (2003); Bennani et al. (2017)):

$$ \mu_s^{fr} = \begin{cases} Fr \frac{\sqrt{2} P_p^{fr} \sin(\eta)}{2\sqrt{D_{s,ij} D_{s,ij} + \psi}}; & \alpha_s > \alpha_s^{min} \\ 0; & \text{else} \end{cases} \qquad (35) $$

where $\eta$ is the internal friction angle (25°) and $\psi = 2/3(q_s^2/d_s^2)$. The interphase momentum transfer between gas and solid after subtraction of the gas pressure gradient effect is written as:

$$ I_{s \to g, i} = -I_{g \to s, i} = \alpha_s \rho_s \frac{1}{\tau_{gs}^F} V_{r,i} \qquad (36) $$

on the basis of the mean relative velocity and the mean particle relaxation time, $\tau_{gs}^F$, which accounts for the drag effect on the particles:

$$ \frac{1}{\tau_{gs}^F} = \frac{3\rho_g \langle |v_r| \rangle_s}{4\rho_s d_s} C_D. \qquad (37) $$

This time is expressed using two different experimental laws, Wen and Yu and Ergun's law, according to the modeling proposed by Gobin et al. (2003):

$$ C_D = \begin{cases} C_{D,WY}; & \alpha_g \ge 0.7 \\ \min[C_{D,WY}; C_{D,Erg}]; & \alpha_g < 0.7 \end{cases} \qquad (38) $$

with

$$ C_{D,Erg} = 200 \frac{1 - \alpha_g}{Re_p} + \frac{7}{3} \qquad (39) $$

$$ C_{D,WY} = \begin{cases} \frac{24}{Re_p} [1 + 0.15 Re_p^{0.687}] \alpha_g^{-1.7}; & Re_p < 1000 \\ 0.44 \alpha_g^{-1.7}; & Re_p \ge 1000 \end{cases} \qquad (40) $$

Here, $Re_p$ is the particle Reynolds number defined as

$$ Re_p = \frac{\alpha_g \rho_g \langle |v_r| \rangle_s d_s}{\mu_g}. \qquad (41) $$

Finally, $\tau_s^c$ is the inter-particle collision time:

$$ \tau_s^c = \left( 6 \frac{\alpha_s g_0}{d_s} \sqrt{\frac{16}{\pi}} \frac{2}{3} \delta q_s^2 \right)^{-1} . \qquad (42) $$

In this study, two different models are used for the uncorrelated contribution of the random kinetic energy in the inter-particle collision time (Equation (42)) and in the kinetic energy dissipation by inelastic collision (Equation (27)). Indeed, according Février et al. (2005) and Fox (2014), we may assume that the random particle kinetic energy $q_s^2$ may be separated in two parts:

$$ q_s^2 = \tilde{q}_s^2 + \delta q_s^2 \qquad (43) $$

where $\tilde{q}_s^2$ is the correlated contribution, representing the collective fluctuating motion of the particles, and $\delta q_s^2$ is the uncorrelated contribution, representing the particle-particle relative
where $Fr = 0.05$, $r = 2$ and $s = 5$ are model parameters, which may be varied depending on the types of particles. In this work, $\alpha_s^{min}$ is set equal to 0.55, which is an appropriate value for spherical particles. The frictional viscosity is modeled as follows (Srivastava and Sundaresan (2003); Bennani et al. (2017)):

$$ \mu_s^{fr} = \begin{cases} Fr \frac{\sqrt{2} P_p^{fr} \sin(\eta)}{2\sqrt{D_{s,ij} D_{s,ij} + \psi}}; & \alpha_s > \alpha_s^{min} \\ 0; & \text{else} \end{cases} \qquad (35) $$

where $\eta$ is the internal friction angle (25°) and $\psi = 2/3(q_s^2/d_s^2)$. The interphase momentum transfer between gas and solid after subtraction of the gas pressure gradient effect is written as:

$$ I_{s \to g, i} = -I_{g \to s, i} = \alpha_s \rho_s \frac{1}{\tau_{gs}^F} V_{r,i} \qquad (36) $$

on the basis of the mean relative velocity and the mean particle relaxation time, $\tau_{gs}^F$, which accounts for the drag effect on the particles:

$$ \frac{1}{\tau_{gs}^F} = \frac{3\rho_g \langle |v_r| \rangle_s}{4\rho_s d_s} C_D. \qquad (37) $$

This time is expressed using two different experimental laws, Wen and Yu and Ergun's law, according to the modeling proposed by Gobin et al. (2003):

$$ C_D = \begin{cases} C_{D,WY}; & \alpha_g \ge 0.7 \\ \min[C_{D,WY}; C_{D,Erg}]; & \alpha_g < 0.7 \end{cases} \qquad (38) $$

with

$$ C_{D,Erg} = 200 \frac{1 - \alpha_g}{Re_p} + \frac{7}{3} \qquad (39) $$

$$ C_{D,WY} = \begin{cases} \frac{24}{Re_p} [1 + 0.15 Re_p^{0.687}] \alpha_g^{-1.7}; & Re_p < 1000 \\ 0.44 \alpha_g^{-1.7}; & Re_p \ge 1000 \end{cases} \qquad (40) $$

Here, $Re_p$ is the particle Reynolds number defined as

$$ Re_p = \frac{\alpha_g \rho_g \langle |v_r| \rangle_s d_s}{\mu_g}. \qquad (41) $$

Finally, $\tau_s^c$ is the inter-particle collision time:

$$ \tau_s^c = \left( 6 \frac{\alpha_s g_0}{d_s} \sqrt{\frac{16}{\pi}} \frac{2}{3} \delta q_s^2 \right)^{-1} . \qquad (42) $$

In this study, two different models are used for the uncorrelated contribution of the random kinetic energy in the inter-particle collision time (Equation (42)) and in the kinetic energy dissipation by inelastic collision (Equation (27)). Indeed, according Février et al. (2005) and Fox (2014), we may assume that the random particle kinetic energy $q_s^2$ may be separated in two parts:

$$ q_s^2 = \tilde{q}_s^2 + \delta q_s^2 \qquad (43) $$

where $\tilde{q}_s^2$ is the correlated contribution, representing the collective fluctuating motion of the particles, and $\delta q_s^2$ is the uncorrelated contribution, representing the particle-particle relative
fluctuating motion (Simonin et al. (2002)). For the uncorrelated model, the correlation effect of the neighboring particles is not taken into consideration and we get the following assumption:

$$ \begin{cases} \tilde{q}_s^2 = 0 \\ \delta q_s^2 = q_s^2 \end{cases} \qquad (44) $$

Such an assumption is corresponding to classic kinetic theory of granular flow (Gidaspow (1994)). For the correlated model, according to the works of Laviéville et al. (1995) and Simonin et al. (2002), it can be expressed as :

$$ \begin{cases} \tilde{q}_s^2 = \zeta_{gs}^2 q_s^2 \\ \delta q_s^2 = [1 - \zeta_{gs}^2] q_s^2 \end{cases} \qquad (45) $$

where $\zeta_{gs}^2$ represents a correlation coefficient and is written as

$$ \zeta_{gs}^2 = \frac{[q_{gs}]^2}{4kq_s^2}, \quad 0 < \zeta_{gs}^2 < 1. \qquad (46) $$

The correlated model will account for correlation between colliding particles due to the inter-action with the fluid turbulence (Février et al. (2005)). When an uncorrelated assumption is used, $\tilde{q}_s^2$ in Equation (31) turns to zero, and the contribution of two-way coupling to the inter-phase term is a destruction term directly proportional to the fluid-particle velocity covariance. When a correlated model is used, this contribution is obtained from the total particle kinetic energy, using the above correlation coefficient (Simonin et al. (2002)). We can notice that according to Equations (45) and (46), we may write, $q_{gs} - 2\tilde{q}_s^2 = q_{gs}(1 - q_{gs}/2k)$ showing that the sign of the two-way contribution is directly depending on the ratio between the fluid-particle velocity covariance and the fluid turbulent kinetic energy. In addition, when the value of $\zeta_{gs}^2$ tends towards zero, that is for very large Stokes numbers in particle-laden turbulent flows, the correlated model reverts to the uncorrelated one.

### 3. Experimental system and simulation setup

In this work, a double-loop CFB reactor system, corresponding to the experimental facility at SINTEF Energy Research (Trondheim, Norway), is adopted to investigate the hydrodynamics of the CLC unit. Two reactors, two cyclones, two loop seals and one lifter are designed and built for this facility, which has also been used to study CLC of gaseous fuels (Langørgen et al. (2017)). The operating schematic diagram is displayed in Figure 1. The dimensions of the CLC unit are reported in Table 1. In the experiments, ilmenite from Titania A/S in Norway (of bulk density 2600 kg/m³ and mean diameter ($D_{50}$) 90 µm) is used as oxygen carrier.

In the numerical simulation, the system is meshed by an O-grid method with approximately 0.7 million cells (the reference case), which is a suitable compromise between fine and coarse mesh considering both accuracy and calculation costs. No-slip or free-slip wall boundary conditions for the mean particle velocity and zero-flux boundary conditions for the particle kinetic energy are imposed (Fede et al. (2016)). Friction conditions are used for the gas phase, according to the $k-\epsilon$ modeling considered in this work. The operating temperature is set to 1273 K, according to the experiments. An overview of the CLC mesh is given in Figure 2.

Mass inventories are calculated from the experimental pressure-drop measurements, and are summarized in the Table 2. In the experiments, each loop seal was designed with three chambers: central, external and internal. The particles separated by the cyclone enter the
239 fluctuating motion (Simonin et al. (2002)). For the uncorrelated model, the correlation effect of the neighboring particles is not taken into consideration and we get the following assumption:

241

$$ \begin{cases} \tilde{q}_s^2 = 0 \\ \delta q_s^2 = q_s^2 \end{cases} \qquad (44) $$

242 Such an assumption is corresponding to classic kinetic theory of granular flow (Gidaspow (1994)). For the correlated model, according to the works of Laviéville et al. (1995) and 243 Simonin et al. (2002), it can be expressed as :

244

$$ \begin{cases} \tilde{q}_s^2 = \zeta_{gs}^2 q_s^2 \\ \delta q_s^2 = [1 - \zeta_{gs}^2] q_s^2 \end{cases} \qquad (45) $$

245

where $\zeta_{gs}^2$ represents a correlation coefficient and is written as

$$ \zeta_{gs}^2 = \frac{[q_{gs}]^2}{4kq_s^2}, \quad 0 < \zeta_{gs}^2 < 1. \qquad (46) $$

246 The correlated model will account for correlation between colliding particles due to the inter-action with the fluid turbulence (Février et al. (2005)). When an uncorrelated assumption is used, $\tilde{q}_s^2$ in Equation (31) turns to zero, and the contribution of two-way coupling to the inter-phase term is a destruction term directly proportional to the fluid-particle velocity covariance. When a correlated model is used, this contribution is obtained from the total particle kinetic energy, using the above correlation coefficient (Simonin et al. (2002)). We can notice that according to Equations (45) and (46), we may write, $q_{gs} - 2\tilde{q}_s^2 = q_{gs}(1 - q_{gs}/2k)$ showing that the sign of the two-way contribution is directly depending on the ratio between the fluid-particle velocity covariance and the fluid turbulent kinetic energy. In addition, when the value of $\zeta_{gs}^2$ tends towards zero, that is for very large Stokes numbers in particle-laden turbulent flows, the correlated model reverts to the uncorrelated one.

247

248

249

250

251

252

253

254

255

256

257

### 3. Experimental system and simulation setup

In this work, a double-loop CFB reactor system, corresponding to the experimental facility at SINTEF Energy Research (Trondheim, Norway), is adopted to investigate the hydrodynamics of the CLC unit. Two reactors, two cyclones, two loop seals and one lifter are designed and built for this facility, which has also been used to study CLC of gaseous fuels (Langørgen et al. (2017)). The operating schematic diagram is displayed in Figure 1. The dimensions of the CLC unit are reported in Table 1. In the experiments, ilmenite from Titania A/S in Norway (of bulk density 2600 kg/m³ and mean diameter ($D_{50}$) 90 µm) is used as oxygen carrier.

258

259

260

261

262

263

264

In the numerical simulation, the system is meshed by an O-grid method with approximately 0.7 million cells (the reference case), which is a suitable compromise between fine and coarse mesh considering both accuracy and calculation costs. No-slip or free-slip wall boundary conditions for the mean particle velocity and zero-flux boundary conditions for the particle kinetic energy are imposed (Fede et al. (2016)). Friction conditions are used for the gas phase, according to the $k-\epsilon$ modeling considered in this work. The operating temperature is set to 1273 K, according to the experiments. An overview of the CLC mesh is given in Figure 2.

265

266

267

268

269

270

271

Mass inventories are calculated from the experimental pressure-drop measurements, and are summarized in the Table 2. In the experiments, each loop seal was designed with three chambers: central, external and internal. The particles separated by the cyclone enter the

272

273

274
<figure><img src="image_2.png" /><figcaption>Figure 1: Scheme of the 150 $kW_{th}$ chemical looping combustion pilot at SINTEF, Norway.</figcaption></figure>

<figure><img src="image_3.png" /><figcaption>Figure 2: Structure of the 150 $kW_{th}$ chemical looping combustion pilot at SINTEF and mesh plan.</figcaption></figure>

central part of the loop seal. Then, the particles are transported to the other reactor through the external chamber, or re-circulated back into the original reactor through the internal chamber. For the current CLC experiments, the particle outlet leg connected with the internal chamber was shut down. For this reason, in the numerical simulation only the central and external chambers were considered (as shown in Figure A.22 for the FR loop seal). Therefore, only two third of the mass inventory of each loop seal was taken into account, in addition to the mass
<figure><img src="image_2.png" /><figcaption>Figure 1: Scheme of the 150 $kW_{th}$ chemical looping combustion pilot at SINTEF, Norway.</figcaption></figure>

<figure><img src="image_3.png" /><figcaption>Figure 2: Structure of the 150 $kW_{th}$ chemical looping combustion pilot at SINTEF and mesh plan.</figcaption></figure>

central part of the loop seal. Then, the particles are transported to the other reactor through the external chamber, or re-circulated back into the original reactor through the internal chamber. For the current CLC experiments, the particle outlet leg connected with the internal chamber was shut down. For this reason, in the numerical simulation only the central and external chambers were considered (as shown in Figure A.22 for the FR loop seal). Therefore, only two third of the mass inventory of each loop seal was taken into account, in addition to the mass
<table><caption>Table 1: Dimension of the CLC unit.</caption><thead><tr><th>Item</th><th>Value</th><th>Units</th></tr></thead><tbody><tr><td>Height of AR</td><td>6.0</td><td>m</td></tr><tr><td>Inner diameter of AR (cylindrical part)</td><td>23.0</td><td>cm</td></tr><tr><td>Height of the AR conical part</td><td>1.0</td><td>m</td></tr><tr><td>Bottom diameter of the AR conical part</td><td>15.0</td><td>cm</td></tr><tr><td>Height of FR (including lifter)</td><td>6.7</td><td>m</td></tr><tr><td>Inner diameter of FR (cylindrical part)</td><td>15.4</td><td>cm</td></tr><tr><td>Height of the FR conical part</td><td>1.0</td><td>m</td></tr><tr><td>Bottom diameter of the FR conical part</td><td>10.0</td><td>cm</td></tr></tbody></table>

<table><caption>Table 2: Mass inventories.</caption><thead><tr><th>Item</th><th>Value</th><th>Units</th></tr></thead><tbody><tr><td>Air reactor</td><td>18.1</td><td>kg</td></tr><tr><td>Fuel reactor</td><td>28.2</td><td>kg</td></tr><tr><td>Lifter</td><td>11.5</td><td>kg</td></tr><tr><td>AR loop seal</td><td>52.1 (37.7)</td><td>kg</td></tr><tr><td>FR loop seal</td><td>42.4 (30.4)</td><td>kg</td></tr><tr><td>Total</td><td>152.3(125.9)</td><td>kg</td></tr></tbody></table>

of the particles contained in each connecting pipe. These values are shown in parentheses in Table 2.

At the initial time, the solid phase is initialized by a solid volume fraction of 0.55, the particle diameter is set to 90 µm and the particle density to 4727 kg/m³ (reference case), corresponding to a mean voidage of 0.45 (cf. Section 4.4). Particle diameter and particle density are kept constant during the numerical simulation. The initial mass distribution in the CLC is set according to the experiments (cf. previous discussion), which is beneficial for shortening the computational time to reach a steady state.

In this study, only one solid phase is considered, the oxygen carrier, while the biomass is taken into account through its end products. Moreover, to estimate as best as possible the amount of gases in the system, the mass transfer of oxygen, from the solid in the fuel reactor, and to the solid in the air reactor, is considered as well. Three different ways of injection are tested. The numerical strategy to account for the change in gas flow rate due to reactions is detailed below.

The biomass is composed of volatiles, char and ash. Volatiles are released into the fuel reactor and take part in OC reduction reactions. Gas products and a part of volatiles are involved in char gasification reactions. Assuming complete reactions in the fuel reactor, and a statistically stationary state, one may estimate the mass flow rate of the whole mixture from the mass flow rate of the final products, CO₂ and H₂O, defined as

$$ \dot{Q}_{CO_2} = \dot{Q}_{bio} Y_{vol} X_{CO_2} + \dot{Q}_{bio} Y_{vol} X_{CO} \frac{W_{CO_2}}{W_{CO}} + \dot{Q}_{bio} Y_{vol} X_{CH_4} \frac{W_{CO_2}}{W_{CH_4}} + \dot{Q}_{bio} Y_{char} \frac{W_{CO_2}}{W_{char}}, \quad (47) $$

$$ \dot{Q}_{H_2O} = \dot{Q}_{bio} Y_{vol} X_{H_2O} + \dot{Q}_{bio} Y_{vol} X_{H_2} \frac{W_{H_2O}}{W_{H_2}} + \dot{Q}_{bio} Y_{vol} X_{CH_4} \frac{2W_{H_2O}}{W_{CH_4}}, \quad (48) $$
<table><caption>Table 1: Dimension of the CLC unit.</caption><thead><tr><th>Item</th><th>Value</th><th>Units</th></tr></thead><tbody><tr><td>Height of AR</td><td>6.0</td><td>m</td></tr><tr><td>Inner diameter of AR (cylindrical part)</td><td>23.0</td><td>cm</td></tr><tr><td>Height of the AR conical part</td><td>1.0</td><td>m</td></tr><tr><td>Bottom diameter of the AR conical part</td><td>15.0</td><td>cm</td></tr><tr><td>Height of FR (including lifter)</td><td>6.7</td><td>m</td></tr><tr><td>Inner diameter of FR (cylindrical part)</td><td>15.4</td><td>cm</td></tr><tr><td>Height of the FR conical part</td><td>1.0</td><td>m</td></tr><tr><td>Bottom diameter of the FR conical part</td><td>10.0</td><td>cm</td></tr></tbody></table>

<table><caption>Table 2: Mass inventories.</caption><thead><tr><th>Item</th><th>Value</th><th>Units</th></tr></thead><tbody><tr><td>Air reactor</td><td>18.1</td><td>kg</td></tr><tr><td>Fuel reactor</td><td>28.2</td><td>kg</td></tr><tr><td>Lifter</td><td>11.5</td><td>kg</td></tr><tr><td>AR loop seal</td><td>52.1 (37.7)</td><td>kg</td></tr><tr><td>FR loop seal</td><td>42.4 (30.4)</td><td>kg</td></tr><tr><td>Total</td><td>152.3(125.9)</td><td>kg</td></tr></tbody></table>

of the particles contained in each connecting pipe. These values are shown in parentheses in Table 2.

At the initial time, the solid phase is initialized by a solid volume fraction of 0.55, the particle diameter is set to 90 µm and the particle density to 4727 kg/m³ (reference case), corresponding to a mean voidage of 0.45 (cf. Section 4.4). Particle diameter and particle density are kept constant during the numerical simulation. The initial mass distribution in the CLC is set according to the experiments (cf. previous discussion), which is beneficial for shortening the computational time to reach a steady state.

In this study, only one solid phase is considered, the oxygen carrier, while the biomass is taken into account through its end products. Moreover, to estimate as best as possible the amount of gases in the system, the mass transfer of oxygen, from the solid in the fuel reactor, and to the solid in the air reactor, is considered as well. Three different ways of injection are tested. The numerical strategy to account for the change in gas flow rate due to reactions is detailed below.

The biomass is composed of volatiles, char and ash. Volatiles are released into the fuel reactor and take part in OC reduction reactions. Gas products and a part of volatiles are involved in char gasification reactions. Assuming complete reactions in the fuel reactor, and a statistically stationary state, one may estimate the mass flow rate of the whole mixture from the mass flow rate of the final products, CO₂ and H₂O, defined as

$$ \dot{Q}_{CO_2} = \dot{Q}_{bio} Y_{vol} X_{CO_2} + \dot{Q}_{bio} Y_{vol} X_{CO} \frac{W_{CO_2}}{W_{CO}} + \dot{Q}_{bio} Y_{vol} X_{CH_4} \frac{W_{CO_2}}{W_{CH_4}} + \dot{Q}_{bio} Y_{char} \frac{W_{CO_2}}{W_{char}}, \quad (47) $$

$$ \dot{Q}_{H_2O} = \dot{Q}_{bio} Y_{vol} X_{H_2O} + \dot{Q}_{bio} Y_{vol} X_{H_2} \frac{W_{H_2O}}{W_{H_2}} + \dot{Q}_{bio} Y_{vol} X_{CH_4} \frac{2W_{H_2O}}{W_{CH_4}}, \quad (48) $$
where $\dot{Q}_{bio}$ is the injection mass flow rate (kg/s) of biomass, $Y_{vol}$ is the mass fraction of volatiles, and $Y_{char}$ the mass fraction of char ($Y_{vol} + Y_{char} = 1 - Y_{ash}$). $X_{\beta}$ is the mass fraction of the species $\beta$ in the volatiles ($\sum_{\beta} X_{\beta} = 1$). The values of $Y_{vol}$ (0.845) and $Y_{char}$ (0.15), as well as the volatile composition ($X_{CO} = 0.5581$, $X_{CO_2} = 0.1594$, $X_{CH_4} = 0.1926$, $X_{H_2} = 0.0183$, and $X_{H_2O} = 0.0716$) are obtained from proximate analysis and heat and mass balance (Thunman et al. (2001)), assuming these five species as the primary volatiles. The mass flow rate of oxygen required to full conversion in the fuel reactor is therefore

$$ \dot{Q}_{O_2} = \dot{Q}_{bio} Y_{vol} X_{CO} \frac{0.5 W_{O_2}}{W_{CO}} + \dot{Q}_{bio} Y_{vol} X_{H_2} \frac{0.5 W_{O_2}}{W_{H_2}} + \dot{Q}_{bio} Y_{vol} X_{CH_4} \frac{2 W_{O_2}}{W_{CH_4}} + \dot{Q}_{bio} Y_{char} \frac{W_{O_2}}{W_{char}}. \quad (49) $$

This mass flow rate represents the amount of oxygen per unit time that is required in the air reactor to return the oxygen carrier to its original oxidization state. The mass flow rate to be added at the fuel reactor injection to reproduce biomass conversion and reduction reactions is $\dot{Q}_{CO_2} + \dot{Q}_{H_2O}$.

<table><caption>Table 3: Inlet mass flow rates: experiments and simulation with lateral injection of products.</caption><thead><tr><th>Item</th><th>Simulation</th><th colspan="2">Experiments</th></tr></thead><tbody><tr><td>AR primary gas inlet</td><td>124.49 kg/h</td><td>N<sub>2</sub>,O<sub>2</sub></td><td>146.67 kg/h Air</td></tr><tr><td>AR secondary gas inlet G1</td><td>14.46 kg/h</td><td>N<sub>2</sub>,O<sub>2</sub></td><td>17.04 kg/h Air</td></tr><tr><td>AR secondary gas inlet G2</td><td>24.75 kg/h</td><td>N<sub>2</sub>,O<sub>2</sub></td><td>29.16 kg/h Air</td></tr><tr><td>FR bottom inlet</td><td>11.66 kg/h</td><td>N<sub>2</sub></td><td>11.66 kg/h N<sub>2</sub></td></tr><tr><td>FR lateral inlet</td><td>49.27 kg/h</td><td>CO<sub>2</sub>,H<sub>2</sub>O</td><td>20.2 kg/h Bio</td></tr><tr><td>Lifter inlet</td><td>2.27 kg/h</td><td>N<sub>2</sub></td><td>2.27 kg/h N<sub>2</sub></td></tr><tr><td>AR loop seal, particle inlet leg</td><td>2.23 kg/h</td><td>N<sub>2</sub></td><td>2.23 kg/h N<sub>2</sub></td></tr><tr><td>AR loop seal, particle outlet leg</td><td>3.21 kg/h</td><td>N<sub>2</sub></td><td>3.21 kg/h N<sub>2</sub></td></tr><tr><td>FR loop seal, particle inlet leg</td><td>2.15 kg/h</td><td>N<sub>2</sub></td><td>2.15 kg/h N<sub>2</sub></td></tr><tr><td>FR loop seal, particle outlet leg</td><td>1.73 kg/h</td><td>N<sub>2</sub></td><td>1.73 kg/h N<sub>2</sub></td></tr></tbody></table>

Three ways of injecting such additional gases to mimic reactions are tested:

* Lateral injection of products: the CO₂-H₂O mixture, corresponding to the whole products from the full biomass conversion and reduction reactions (therefore accounting for the oxygen from ilmenite), is injected in the fuel reactor from the lateral inlet, which is at the same location as the fuel particle inlet, according to the experimental configuration; the oxygen consumed by the oxidation is directly removed from the AR inlet. In this case, the primary and secondary inlets to AR are reduced with 15% compared to the experimental values.

* Bottom injection of products: the CO₂-H₂O mixture is injected in the fuel reactor from the bottom inlet, together with the fluidizing gas (N₂); the inlet conditions in the AR are the same as above.

* Lateral injection of a part of the products while using source terms for the oxygen transfer: inlet conditions for the air reactor are kept the same as in the experiments. The change of flow rate due to the mass transfer between the oxygen carrier and the gas phase is taken into consideration by source terms for both the air and fuel reactors. This method of injection allows part of the gases to be distributed inside the reactors in proportion to
where $\dot{Q}_{bio}$ is the injection mass flow rate (kg/s) of biomass, $Y_{vol}$ is the mass fraction of volatiles, and $Y_{char}$ the mass fraction of char ($Y_{vol} + Y_{char} = 1 - Y_{ash}$). $X_{\beta}$ is the mass fraction of the species $\beta$ in the volatiles ($\sum_{\beta} X_{\beta} = 1$). The values of $Y_{vol}$ (0.845) and $Y_{char}$ (0.15), as well as the volatile composition ($X_{CO} = 0.5581$, $X_{CO_2} = 0.1594$, $X_{CH_4} = 0.1926$, $X_{H_2} = 0.0183$, and $X_{H_2O} = 0.0716$) are obtained from proximate analysis and heat and mass balance (Thunman et al. (2001)), assuming these five species as the primary volatiles. The mass flow rate of oxygen required to full conversion in the fuel reactor is therefore

$$ \dot{Q}_{O_2} = \dot{Q}_{bio} Y_{vol} X_{CO} \frac{0.5 W_{O_2}}{W_{CO}} + \dot{Q}_{bio} Y_{vol} X_{H_2} \frac{0.5 W_{O_2}}{W_{H_2}} + \dot{Q}_{bio} Y_{vol} X_{CH_4} \frac{2 W_{O_2}}{W_{CH_4}} + \dot{Q}_{bio} Y_{char} \frac{W_{O_2}}{W_{char}}. \quad (49) $$

This mass flow rate represents the amount of oxygen per unit time that is required in the air reactor to return the oxygen carrier to its original oxidization state. The mass flow rate to be added at the fuel reactor injection to reproduce biomass conversion and reduction reactions is $\dot{Q}_{CO_2} + \dot{Q}_{H_2O}$.

<table><caption>Table 3: Inlet mass flow rates: experiments and simulation with lateral injection of products.</caption><thead><tr><th>Item</th><th>Simulation</th><th>Experiments</th></tr></thead><tbody><tr><td>AR primary gas inlet</td><td>124.49 kg/h N<sub>2</sub>O<sub>2</sub></td><td>146.67 kg/h Air</td></tr><tr><td>AR secondary gas inlet G1</td><td>14.46 kg/h N<sub>2</sub>O<sub>2</sub></td><td>17.04 kg/h Air</td></tr><tr><td>AR secondary gas inlet G2</td><td>24.75 kg/h N<sub>2</sub>O<sub>2</sub></td><td>29.16 kg/h Air</td></tr><tr><td>FR bottom inlet</td><td>11.66 kg/h N<sub>2</sub></td><td>11.66 kg/h N<sub>2</sub></td></tr><tr><td>FR lateral inlet</td><td>49.27 kg/h CO<sub>2</sub>H<sub>2</sub>O</td><td>20.2 kg/h Bio</td></tr><tr><td>Lifter inlet</td><td>2.27 kg/h N<sub>2</sub></td><td>2.27 kg/h N<sub>2</sub></td></tr><tr><td>AR loop seal, particle inlet leg</td><td>2.23 kg/h N<sub>2</sub></td><td>2.23 kg/h N<sub>2</sub></td></tr><tr><td>AR loop seal, particle outlet leg</td><td>3.21 kg/h N<sub>2</sub></td><td>3.21 kg/h N<sub>2</sub></td></tr><tr><td>FR loop seal, particle inlet leg</td><td>2.15 kg/h N<sub>2</sub></td><td>2.15 kg/h N<sub>2</sub></td></tr><tr><td>FR loop seal, particle outlet leg</td><td>1.73 kg/h N<sub>2</sub></td><td>1.73 kg/h N<sub>2</sub></td></tr></tbody></table>

Three ways of injecting such additional gases to mimic reactions are tested:

* Lateral injection of products: the CO₂-H₂O mixture, corresponding to the whole products from the full biomass conversion and reduction reactions (therefore accounting for the oxygen from ilmenite), is injected in the fuel reactor from the lateral inlet, which is at the same location as the fuel particle inlet, according to the experimental configuration; the oxygen consumed by the oxidation is directly removed from the AR inlet. In this case, the primary and secondary inlets to AR are reduced with 15% compared to the experimental values.

* Bottom injection of products: the CO₂-H₂O mixture is injected in the fuel reactor from the bottom inlet, together with the fluidizing gas (N₂); the inlet conditions in the AR are the same as above.

* Lateral injection of a part of the products while using source terms for the oxygen transfer: inlet conditions for the air reactor are kept the same as in the experiments. The change of flow rate due to the mass transfer between the oxygen carrier and the gas phase is taken into consideration by source terms for both the air and fuel reactors. This method of injection allows part of the gases to be distributed inside the reactors in proportion to
the local amount of solid. In the fuel reactor, the source terms for CO₂ and H₂O in each computational cell are computed as follows:

$$ \Gamma_{CO_2}^{FR}(\mathbf{x}, t) = \dot{Q}_{O_2} Y_{CO_2}^* \frac{\alpha_{oc}(\mathbf{x}, t) \rho_{oc}(\mathbf{x}, t)}{m_{oc,FR}(t)}, \qquad \Gamma_{H_2O}^{FR} = \dot{Q}_{O_2} Y_{H_2O}^* \frac{\alpha_{oc}(\mathbf{x}, t) \rho_{oc}(\mathbf{x}, t)}{m_{oc,FR}(t)}, \quad (50) $$

where $Y_{CO_2}^*$ and $Y_{H_2O}^*$ are the mass fractions of CO₂ and H₂O in the FR products, according to a full conversion assumption (i.e. $Y_{CO_2}^* + Y_{H_2O}^* = 1.0$). $\alpha_{oc}(\mathbf{x}, t)$ and $\rho_{oc}(\mathbf{x}, t)$ are, respectively, the particle volume fraction and density in the corresponding computational cell. $m_{oc,FR}(t)$ is the instantaneous total solid mass in the fuel reactor.

For the air reactor, the source term is written as

$$ \Gamma_{O_2}^{AR} = -\dot{Q}_{O_2} \frac{\alpha_{oc}(\mathbf{x}, t) \rho_{oc}(\mathbf{x}, t)}{m_{oc,AR}(t)}, $$

where $m_{oc,AR}(t)$ is the instantaneous total mass of oxygen carrier in the air reactor. This term is negative because of the oxidation reaction, which leads to a mass transfer of oxygen from the gas to the solid phase.

The evaluation of the injection methods will be shown in section 4.1.

The gas and solid flow rates from experiments are given in Table 3, as well as the gas flow rates considered in the numerical simulation according to the lateral injection method (the first presented above). Loop-seal and lifter injection rates correspond to the experimental operating conditions, except for the injection temperature (here set to 1273 K).

# 4. Results and discussion

Hereafter, if not otherwise mentioned, results refer to the uncorrelated model case (cf. Equation (44)), with a free-slip mean particle velocity wall boundary condition. Results indicate that the system reaches a hydrodynamic steady state from about 15 seconds. Accordingly, to get statistics from the numerical simulations (pressure, velocity, solid mass flow rate, etc.), time-averaging of the results starts after 15 seconds of physical time.

## 4.1. Injection method

First, the three different methods of injection are tested to assess their effect on the numerical predictions. Time-averaged pressure profiles obtained by the numerical simulations are shown in Figure 3. In the figure, experimental pressure measurements are also shown. The given pressures are gauge pressures, with the atmosphere as the zero reference. The numerical results are almost identical to each other in the air reactor, as well as in the lifter. In the fuel reactor, however, the bed expansion is slightly higher between 0.5 and 3 meters when injecting from the bottom or when using local source terms to account for oxygen mass transfer. The results are consistent with each other, though. Indeed, injecting from the bottom, or directly introducing the gas in the reactor by source terms, should lead to a more homogeneous bed for which a larger expansion is expected. This is more pronounced when introducing the gas in each computational cell.

Globally, results indicate that the air reactor operates in a circulating regime, with a denser part at the bottom, as expected. In contrast, the pressure profiles in the fuel reactor exhibit two smoothly connected linear trends, typical of bubbling fluidized beds. The slope of the pressure in the upper part of the fuel reactor indicates however that a part of the solid leaves
the local amount of solid. In the fuel reactor, the source terms for CO₂ and H₂O in each computational cell are computed as follows:

$$ \Gamma_{CO_2}^{FR}(\mathbf{x}, t) = \dot{Q}_{O_2} Y_{CO_2}^* \frac{\alpha_{oc}(\mathbf{x}, t) \rho_{oc}(\mathbf{x}, t)}{m_{oc,FR}(t)}, \qquad \Gamma_{H_2O}^{FR} = \dot{Q}_{O_2} Y_{H_2O}^* \frac{\alpha_{oc}(\mathbf{x}, t) \rho_{oc}(\mathbf{x}, t)}{m_{oc,FR}(t)}, \quad (50) $$

where $Y_{CO_2}^*$ and $Y_{H_2O}^*$ are the mass fractions of CO₂ and H₂O in the FR products, according to a full conversion assumption (i.e. $Y_{CO_2}^* + Y_{H_2O}^* = 1.0$). $\alpha_{oc}(\mathbf{x}, t)$ and $\rho_{oc}(\mathbf{x}, t)$ are, respectively, the particle volume fraction and density in the corresponding computational cell. $m_{oc,FR}(t)$ is the instantaneous total solid mass in the fuel reactor.

For the air reactor, the source term is written as

$$ \Gamma_{O_2}^{AR} = -\dot{Q}_{O_2} \frac{\alpha_{oc}(\mathbf{x}, t) \rho_{oc}(\mathbf{x}, t)}{m_{oc,AR}(t)}, $$

where $m_{oc,AR}(t)$ is the instantaneous total mass of oxygen carrier in the air reactor. This term is negative because of the oxidation reaction, which leads to a mass transfer of oxygen from the gas to the solid phase.

The evaluation of the injection methods will be shown in section 4.1.

The gas and solid flow rates from experiments are given in Table 3, as well as the gas flow rates considered in the numerical simulation according to the lateral injection method (the first presented above). Loop-seal and lifter injection rates correspond to the experimental operating conditions, except for the injection temperature (here set to 1273 K).

# 4. Results and discussion

Hereafter, if not otherwise mentioned, results refer to the uncorrelated model case (cf. Equation (44)), with a free-slip mean particle velocity wall boundary condition. Results indicate that the system reaches a hydrodynamic steady state from about 15 seconds. Accordingly, to get statistics from the numerical simulations (pressure, velocity, solid mass flow rate, etc.), time-averaging of the results starts after 15 seconds of physical time.

## 4.1. Injection method

First, the three different methods of injection are tested to assess their effect on the numerical predictions. Time-averaged pressure profiles obtained by the numerical simulations are shown in Figure 3. In the figure, experimental pressure measurements are also shown. The given pressures are gauge pressures, with the atmosphere as the zero reference. The numerical results are almost identical to each other in the air reactor, as well as in the lifter. In the fuel reactor, however, the bed expansion is slightly higher between 0.5 and 3 meters when injecting from the bottom or when using local source terms to account for oxygen mass transfer. The results are consistent with each other, though. Indeed, injecting from the bottom, or directly introducing the gas in the reactor by source terms, should lead to a more homogeneous bed for which a larger expansion is expected. This is more pronounced when introducing the gas in each computational cell.

Globally, results indicate that the air reactor operates in a circulating regime, with a denser part at the bottom, as expected. In contrast, the pressure profiles in the fuel reactor exhibit two smoothly connected linear trends, typical of bubbling fluidized beds. The slope of the pressure in the upper part of the fuel reactor indicates however that a part of the solid leaves
<figure><img src="image_4.png" /><figcaption>Figure 3: Time-averaged pressure using three different injection methods.</figcaption></figure>

the reactor at the top. The fuel reactor operates therefore in a mixed regime, as observed in the experiments.

From now on, the lateral injection will be used for the numerical simulations. This injection method along with the uncorrelated model and the free-slip boundary condition represent our reference case. In the reference case, the particle density is set to 4727 kg/m³ (as already mentioned). The question of the particle density is discussed in Section 4.4.

## 4.2. Flow pattern

Instantaneous concentrations of the particulate phase (oxygen carrier) are computed and displayed in Figures 4 and 5, at different times, to provide information of the flow evolution in the system. In the air reactor, particles are fluidized and carried by the fluidization gas, then they are separated by the cyclone and fall down into the loop seal. After passing the loop seal, particles then enter the fuel reactor. The lifter is located between the air reactor and the fuel reactor and works as an additional connector for transporting particles from the fuel reactor to the air reactor. As can be seen, there are bubbles between the fuel reactor and the lifter. In a real reactive case, attention should be paid to the flow rate of fuel particles entering the lifter from the fuel reactor, which could dramatically affect the carbon capture efficiency, since such particles will be transported to the air reactor where they will burn, producing CO₂. Results show that the solid volume fraction is higher in the fuel reactor than in the air reactor because of the different gas velocity. The evolution in time of the solid volume fraction within the fuel reactor is shown in Figure A.23, on a plane located in the middle of the reactor. To get a better
<figure><img src="image_4.png" /><figcaption>Figure 3: Time-averaged pressure using three different injection methods.</figcaption></figure>

the reactor at the top. The fuel reactor operates therefore in a mixed regime, as observed in the experiments.

From now on, the lateral injection will be used for the numerical simulations. This injection method along with the uncorrelated model and the free-slip boundary condition represent our reference case. In the reference case, the particle density is set to 4727 kg/m³ (as already mentioned). The question of the particle density is discussed in Section 4.4.

## 4.2. Flow pattern

Instantaneous concentrations of the particulate phase (oxygen carrier) are computed and displayed in Figures 4 and 5, at different times, to provide information of the flow evolution in the system. In the air reactor, particles are fluidized and carried by the fluidization gas, then they are separated by the cyclone and fall down into the loop seal. After passing the loop seal, particles then enter the fuel reactor. The lifter is located between the air reactor and the fuel reactor and works as an additional connector for transporting particles from the fuel reactor to the air reactor. As can be seen, there are bubbles between the fuel reactor and the lifter. In a real reactive case, attention should be paid to the flow rate of fuel particles entering the lifter from the fuel reactor, which could dramatically affect the carbon capture efficiency, since such particles will be transported to the air reactor where they will burn, producing CO₂. Results show that the solid volume fraction is higher in the fuel reactor than in the air reactor because of the different gas velocity. The evolution in time of the solid volume fraction within the fuel reactor is shown in Figure A.23, on a plane located in the middle of the reactor. To get a better
<figure><img src="image_5.png" /><figcaption>Figure 4: Instantaneous solid volume fraction: air reactor view.</figcaption></figure>

<figure><img src="image_6.png" /><figcaption>Figure 5: Instantaneous solid volume fraction: fuel reactor and lifter view.</figcaption></figure>

386 look, the view is zoomed from 0 to 3 meters in height. The bubble formation, breakage and

387 the flow state can be observed.
<figure><img src="image_5.png" /><figcaption>Figure 4: Instantaneous solid volume fraction: air reactor view.</figcaption></figure>

<figure><img src="image_6.png" /><figcaption>Figure 5: Instantaneous solid volume fraction: fuel reactor and lifter view.</figcaption></figure>

386 look, the view is zoomed from 0 to 3 meters in height. The bubble formation, breakage and

387 the flow state can be observed.
<figure><img src="image_7.png" /><figcaption>Figure 6: Time-averaged pressure in the air reactor. Comparison between collision models (left) and between mean particle velocity boundary conditions while using the same uncorrelated model (right).</figcaption></figure>

### 4.3. Collision model and mean particle velocity wall boundary condition

In this section, the effects of the collision model and the mean particle wall velocity boundary condition on the numerical predictions are analyzed. To properly read the results that will follow, it is important to keep in mind that, in the numerical simulations, the pressure at the outlet was estimated from the experimental pressure in each reactor at the same height, subtracting the pressure between the reactor and its corresponding cyclone. Such an estimate was computed using the free-slip boundary condition, and the corresponding value set as an outlet condition for all the numerical simulations, including those that use the no-slip boundary condition. For the latter, a shift in the pressure profile with respect to the experimental measurement is therefore expected. This point does not deserve further analysis since it is due solely to the pressure outlet conditions.

Figure 6 compares numerical results and experimental measurements of the time-averaged pressure in the air reactor. In Figure 6 (left) both the uncorrelated and correlated model results are shown, using the same free-slip boundary condition. The results reveal only little difference between the two numerical model predictions in the air reactor. The pressure profiles in the fuel reactor and lifter are shown in Figure 7 (top). No appreciable difference is found between the uncorrelated and correlated model predictions in these zones.

The solid mass corresponding to each element of the CLC system is calculated by a volume integral using the time-averaged solid volume fraction, together with the constant particle density. Results are listed in Table 4. Comparison between correlated and uncorrelated models confirms very few differences in the solid mass distribution as well. They are slightly more pronounced in the air reactor, where the particulate phase is more dilute and the effects of a correlated contribution to the particle velocity fluctuation, due to the interactions with the fluid, should be more important than in the fuel reactor and lifter. This point will be discussed further in Section 4.6.

The effect of the mean particle velocity boundary condition on the numerical predictions is then analyzed. No-slip and free-slip wall boundary conditions correspond to the limit cases of maximum particle wall friction effect and pure elastic frictionless particle bouncing, respectively. The results indicate that the wall boundary conditions significantly affect the solid flow behavior. As shown by Figure 6 (right), the pressure obtained using the no-slip wall boundary condition corresponds to an increase of the solid entrainment in the air reactor. To exclude any substantial dependency on the coupling of the different parts of the system, an additional
<figure><img src="image_7.png" /><figcaption>Figure 6: Time-averaged pressure in the air reactor. Comparison between collision models (left) and between mean particle velocity boundary conditions while using the same uncorrelated model (right).</figcaption></figure>

### 4.3. Collision model and mean particle velocity wall boundary condition

In this section, the effects of the collision model and the mean particle wall velocity boundary condition on the numerical predictions are analyzed. To properly read the results that will follow, it is important to keep in mind that, in the numerical simulations, the pressure at the outlet was estimated from the experimental pressure in each reactor at the same height, subtracting the pressure between the reactor and its corresponding cyclone. Such an estimate was computed using the free-slip boundary condition, and the corresponding value set as an outlet condition for all the numerical simulations, including those that use the no-slip boundary condition. For the latter, a shift in the pressure profile with respect to the experimental measurement is therefore expected. This point does not deserve further analysis since it is due solely to the pressure outlet conditions.

Figure 6 compares numerical results and experimental measurements of the time-averaged pressure in the air reactor. In Figure 6 (left) both the uncorrelated and correlated model results are shown, using the same free-slip boundary condition. The results reveal only little difference between the two numerical model predictions in the air reactor. The pressure profiles in the fuel reactor and lifter are shown in Figure 7 (top). No appreciable difference is found between the uncorrelated and correlated model predictions in these zones.

The solid mass corresponding to each element of the CLC system is calculated by a volume integral using the time-averaged solid volume fraction, together with the constant particle density. Results are listed in Table 4. Comparison between correlated and uncorrelated models confirms very few differences in the solid mass distribution as well. They are slightly more pronounced in the air reactor, where the particulate phase is more dilute and the effects of a correlated contribution to the particle velocity fluctuation, due to the interactions with the fluid, should be more important than in the fuel reactor and lifter. This point will be discussed further in Section 4.6.

The effect of the mean particle velocity boundary condition on the numerical predictions is then analyzed. No-slip and free-slip wall boundary conditions correspond to the limit cases of maximum particle wall friction effect and pure elastic frictionless particle bouncing, respectively. The results indicate that the wall boundary conditions significantly affect the solid flow behavior. As shown by Figure 6 (right), the pressure obtained using the no-slip wall boundary condition corresponds to an increase of the solid entrainment in the air reactor. To exclude any substantial dependency on the coupling of the different parts of the system, an additional
<table><caption>Table 4: Solid mass distribution using different collision models and mean particle velocity wall boundary conditions (units: kg).</caption><thead><tr><th></th><th>CASE 1<br>uncorrelated<br>free-slip</th><th>CASE 2<br>correlated<br>free-slip</th><th>CASE 3<br>uncorrelated<br>no-slip</th></tr></thead><tbody><tr><td>Air reactor</td><td>15.656</td><td>15.194</td><td>12.662</td></tr><tr><td>Fuel reactor</td><td>24.597</td><td>24.889</td><td>26.768</td></tr><tr><td>Lifter</td><td>19.698</td><td>19.832</td><td>18.895</td></tr><tr><td>AR loop seal and cyclone</td><td>34.277</td><td>34.308</td><td>35.908</td></tr><tr><td>FR loop seal and cyclone</td><td>30.771</td><td>30.777</td><td>30.707</td></tr><tr><td>Total mass</td><td>125.00</td><td>125.00</td><td>125.00</td></tr></tbody></table>

<figure><img src="image_8.png" /><figcaption>Figure 7: Time-averaged pressure in the fuel reactor and lifter. Comparison between collision models (top), and between mean particle velocity boundary conditions while using the same uncorrelated model (bottom).</figcaption></figure>

420 simulation was carried out. Starting from the end (i.e. 40 seconds of physical time) of the simulation that uses uncorrelated model and no-slip boundary condition everywhere, the free-slip boundary condition was applied to the air reactor only. This simulation was run for an additional 45 seconds, and averages in time were computed from 65 to 85 seconds. The results (Figure 6 (right)) show that the pressure profile mostly returns to the one obtained when the free-slip boundary condition was used everywhere in the system. This demonstrates that the
<table><caption>Table 4: Solid mass distribution using different collision models and mean particle velocity wall boundary conditions (units: kg).</caption><thead><tr><th></th><th>CASE 1<br>uncorrelated<br>free-slip</th><th>CASE 2<br>correlated<br>free-slip</th><th>CASE 3<br>uncorrelated<br>no-slip</th></tr></thead><tbody><tr><td>Air reactor</td><td>15.656</td><td>15.194</td><td>12.662</td></tr><tr><td>Fuel reactor</td><td>24.597</td><td>24.889</td><td>26.768</td></tr><tr><td>Lifter</td><td>19.698</td><td>19.832</td><td>18.895</td></tr><tr><td>AR loop seal and cyclone</td><td>34.277</td><td>34.308</td><td>35.908</td></tr><tr><td>FR loop seal and cyclone</td><td>30.771</td><td>30.777</td><td>30.707</td></tr><tr><td>Total mass</td><td>125.00</td><td>125.00</td><td>125.00</td></tr></tbody></table>

<figure><img src="image_8.png" /><figcaption>Figure 7: Time-averaged pressure in the fuel reactor and lifter. Comparison between collision models (top), and between mean particle velocity boundary conditions while using the same uncorrelated model (bottom).</figcaption></figure>

420 simulation was carried out. Starting from the end (i.e. 40 seconds of physical time) of the simulation that uses uncorrelated model and no-slip boundary condition everywhere, the free-slip boundary condition was applied to the air reactor only. This simulation was run for an additional 45 seconds, and averages in time were computed from 65 to 85 seconds. The results (Figure 6 (right)) show that the pressure profile mostly returns to the one obtained when the free-slip boundary condition was used everywhere in the system. This demonstrates that the
changes in pressure predictions are mainly due to the flow behavior within the air reactor rather than on the coupling of the entire system. We can conclude that the no-slip condition in the AR is the reason for a larger extension of the linear pressure gradient, corresponding to a reduced acceleration region of the solid phases, which leads to a more efficient entrainment by the gas flow.

<figure><img src="image_9.png" /><figcaption>Figure 8: Effect of the mean particle velocity boundary condition on the entrainment estimated as the ratio between the solid mass flow rate and the solid mass.</figcaption></figure>

Only smaller differences are observed when comparing free-slip and no-slip boundary condition results in the fuel reactor and lifter (Figure 7 (bottom)). The no-slip boundary condition leads to slightly higher expansion of the denser part of the bed in the fuel reactor (roughly estimated as the point of intersection of the tangents corresponding to the two linear pressure distributions in the pressure profiles). This observation is consistent with the results of the work of Fede et al. (2016) who showed that a no-slip boundary condition acts at reducing the downward solid mass flux at the walls in a dense fluidized bed, leading to a more expanded bed. The no-slip boundary condition also affects the slope of the pressure profile in the lifter. However, Figure 7 (bottom right) shows that, for the lifter, the effect of the coupling is more important than the boundary condition itself. Globally, in the fuel reactor, the agreement between the numerical results and experimental data is good, except in the penultimate measurement point, while an overestimation of the pressure is observed at the bottom of the lifter.

In contrast to what is observed when comparing the two collision models, the mass distribution changes considerably with the boundary condition (see Table 4). In particular, the mass increases in the fuel reactor and decreases in the air reactor when the no-slip boundary condition is used. In order to estimate the effect of the boundary conditions on the entrainment, which is strictly related to the mass distribution, the ratio between the solid mass flow rate and the solid mass in each relevant part of the system is computed. Molodtsov (2003) (and references cited in) showed indeed that a fully developed gas-particle flow in dilute regime (typical of circulating fluidized beds) exhibits a linear dependency of the solid flux on the solid concentration at a given superficial gas velocity, according to the regime. Since the gas flow rate in the air reactor is almost the same, regardless of the boundary condition, the solid mass flow rate and solid mass are expected to be linearly related, at the operating conditions considered here. Therefore, also their ratio should be almost the same if the boundary condition had no effect on the numerical predictions. Figure 8 shows instead that this ratio increases with the no-slip boundary condition in the air reactor, confirming the conclusions drawn from the pressure profiles.
changes in pressure predictions are mainly due to the flow behavior within the air reactor rather than on the coupling of the entire system. We can conclude that the no-slip condition in the AR is the reason for a larger extension of the linear pressure gradient, corresponding to a reduced acceleration region of the solid phases, which leads to a more efficient entrainment by the gas flow.

<figure><img src="image_9.png" /><figcaption>Figure 8: Effect of the mean particle velocity boundary condition on the entrainment estimated as the ratio between the solid mass flow rate and the solid mass.</figcaption></figure>

Only smaller differences are observed when comparing free-slip and no-slip boundary condition results in the fuel reactor and lifter (Figure 7 (bottom)). The no-slip boundary condition leads to slightly higher expansion of the denser part of the bed in the fuel reactor (roughly estimated as the point of intersection of the tangents corresponding to the two linear pressure distributions in the pressure profiles). This observation is consistent with the results of the work of Fede et al. (2016) who showed that a no-slip boundary condition acts at reducing the downward solid mass flux at the walls in a dense fluidized bed, leading to a more expanded bed. The no-slip boundary condition also affects the slope of the pressure profile in the lifter. However, Figure 7 (bottom right) shows that, for the lifter, the effect of the coupling is more important than the boundary condition itself. Globally, in the fuel reactor, the agreement between the numerical results and experimental data is good, except in the penultimate measurement point, while an overestimation of the pressure is observed at the bottom of the lifter.

In contrast to what is observed when comparing the two collision models, the mass distribution changes considerably with the boundary condition (see Table 4). In particular, the mass increases in the fuel reactor and decreases in the air reactor when the no-slip boundary condition is used. In order to estimate the effect of the boundary conditions on the entrainment, which is strictly related to the mass distribution, the ratio between the solid mass flow rate and the solid mass in each relevant part of the system is computed. Molodtsov (2003) (and references cited in) showed indeed that a fully developed gas-particle flow in dilute regime (typical of circulating fluidized beds) exhibits a linear dependency of the solid flux on the solid concentration at a given superficial gas velocity, according to the regime. Since the gas flow rate in the air reactor is almost the same, regardless of the boundary condition, the solid mass flow rate and solid mass are expected to be linearly related, at the operating conditions considered here. Therefore, also their ratio should be almost the same if the boundary condition had no effect on the numerical predictions. Figure 8 shows instead that this ratio increases with the no-slip boundary condition in the air reactor, confirming the conclusions drawn from the pressure profiles.
<figure><img src="image_10.png" /><figcaption>Figure 9: Radial profiles of the time-averaged solid volume fraction in air reactor, fuel reactor and lifter at different heights, depending on the mean particle velocity boundary condition.</figcaption></figure>

The radial profiles of the time-averaged solid volume fraction in the system are shown in Figure 9. No significant effect of the boundary condition is observed in the fuel reactor and lifter, away from the injections. In contrast, results indicate that in the air reactor the no-slip boundary condition leads to lower values of the solid volume fraction at the wall than the free-slip condition. A possible explanation is the additional production of the particle fluctuant kinetic energy, $q_s^2$, due to the larger values of the particle velocity gradient imposed at the wall by the no-slip condition. As a matter of fact, such larger values of $q_s^2$ close to the wall lead to a “turbophoresis” effect that pushes the particles back towards the core of the flow. The time-averaged particle fluctuant kinetic energy is displayed in Figure 10 on a selected plane in the center of the system.

The radial profiles of the time-averaged solid vertical velocity are shown in Figure 11. In the fuel reactor, as expected, particles flow up in the center and down near the wall. Negative solid velocities are found in the air reactor as well, also at the top. For the lifter, the trend is similar, while the values are lower. The effect of the mean particle velocity wall boundary condition on the particle velocity is not really conclusive. Radial profiles are not symmetric in the air reactor, and it is unclear whether this asymmetry depends on the convergence of the numerical simulation or on the influence of the injections, even at these heights. We can however observe a decrease of the solid axial velocity in the air reactor when a no-slip condition is used at these two heights.

In conclusion, results suggest that while a no-slip condition can be considered satisfactory
<figure><img src="image_10.png" /><figcaption>Figure 9: Radial profiles of the time-averaged solid volume fraction in air reactor, fuel reactor and lifter at different heights, depending on the mean particle velocity boundary condition.</figcaption></figure>

The radial profiles of the time-averaged solid volume fraction in the system are shown in Figure 9. No significant effect of the boundary condition is observed in the fuel reactor and lifter, away from the injections. In contrast, results indicate that in the air reactor the no-slip boundary condition leads to lower values of the solid volume fraction at the wall than the free-slip condition. A possible explanation is the additional production of the particle fluctuant kinetic energy, $q_s^2$, due to the larger values of the particle velocity gradient imposed at the wall by the no-slip condition. As a matter of fact, such larger values of $q_s^2$ close to the wall lead to a “turbophoresis” effect that pushes the particles back towards the core of the flow. The time-averaged particle fluctuant kinetic energy is displayed in Figure 10 on a selected plane in the center of the system.

The radial profiles of the time-averaged solid vertical velocity are shown in Figure 11. In the fuel reactor, as expected, particles flow up in the center and down near the wall. Negative solid velocities are found in the air reactor as well, also at the top. For the lifter, the trend is similar, while the values are lower. The effect of the mean particle velocity wall boundary condition on the particle velocity is not really conclusive. Radial profiles are not symmetric in the air reactor, and it is unclear whether this asymmetry depends on the convergence of the numerical simulation or on the influence of the injections, even at these heights. We can however observe a decrease of the solid axial velocity in the air reactor when a no-slip condition is used at these two heights.

In conclusion, results suggest that while a no-slip condition can be considered satisfactory
<figure><img src="image_11.png" /><figcaption>Figure 10: Time-averaged particle fluctuant kinetic energy, $q_s^2$, on a plane in the middle of the system using free-slip (left) and no-slip (right) mean particle velocity wall boundary conditions.</figcaption></figure>

in a dense regime (fuel reactor, lifter), its use is very questionable in dilute zones (air reactor).

4.4. Particle density and solid mass inventory

The bulk density of the oxygen carrier is given by the experiments and has a value of 2600 kg/m³. To recover the value of the particle density (needed for the numerical simulations) one must know the mean packed-bed voidage, which depends on several parameters including the particle shape. Results presented so far were obtained with a particle density of 4727 kg/m³, corresponding to a mean voidage of 0.45. This value is consistent with the specific gravity of ilmenite found in the literature. However, the oxygen carrier is not pure ilmenite and its composition also changes with redox cycles. Abad et al. (2011) reported a lower true density for the ilmenite oxygen carrier, also depending on the particle state (pre-oxidized or activated), and an increasing porosity of the most oxidized state with particle activation. Since the exact value to be attributed to a spherical particle modeling the real material (including pores) is a priori unknown, additional numerical simulations were carried out using a lower density to investigate the influence of the particle density on the numerical predictions. A particle density corresponding to the maximum packing (0.64) was considered, i.e. 4062 kg/m³. The numerical predictions of the pressure using the two different particle densities are shown in Figure 12.
<figure><img src="image_11.png" /><figcaption>Figure 10: Time-averaged particle fluctuant kinetic energy, $q_s^2$, on a plane in the middle of the system using free-slip (left) and no-slip (right) mean particle velocity wall boundary conditions.</figcaption></figure>

in a dense regime (fuel reactor, lifter), its use is very questionable in dilute zones (air reactor).

4.4. *Particle density and solid mass inventory*

The bulk density of the oxygen carrier is given by the experiments and has a value of 2600 kg/m³. To recover the value of the particle density (needed for the numerical simulations) one must know the mean packed-bed voidage, which depends on several parameters including the particle shape. Results presented so far were obtained with a particle density of 4727 kg/m³, corresponding to a mean voidage of 0.45. This value is consistent with the specific gravity of ilmenite found in the literature. However, the oxygen carrier is not pure ilmenite and its composition also changes with redox cycles. Abad et al. (2011) reported a lower true density for the ilmenite oxygen carrier, also depending on the particle state (pre-oxidized or activated), and an increasing porosity of the most oxidized state with particle activation. Since the exact value to be attributed to a spherical particle modeling the real material (including pores) is a priori unknown, additional numerical simulations were carried out using a lower density to investigate the influence of the particle density on the numerical predictions. A particle density corresponding to the maximum packing (0.64) was considered, i.e. 4062 kg/m³. The numerical predictions of the pressure using the two different particle densities are shown in Figure 12.
<figure><img src="image_12.png" /><figcaption>Figure 11: Radial profiles of the time-averaged solid vertical velocity in air reactor, fuel reactor and lifter at different heights, depending on the mean particle velocity boundary condition.</figcaption></figure>

The corresponding integrated mass distribution is given in Table 5. The results show that the pressure is overestimated by the lowest particle density in the air reactor. Globally, the mass in the air reactor is higher, and comes mainly from the loop seals and lifter (see Table 5). In the fuel reactor, the pressure is also overestimated by the lowest particle density with respect to the highest density. A greater difference is observed in the lifter where the density makes the slope of the pressure to change. The pressure is better predicted at the bottom than at the top in this case (because of the coupling effect).

In order to investigate the effect of the mass inventory, an additional numerical simulation was carried out. In this simulation the mass was decreased in order to improve the numerical predictions obtained with the lowest density. Comparing the two cases with the same particle density and different mass inventory, it turns out that the smaller the mass, the smaller the pressure at the bottom of the air reactor, as well as in the fuel reactor, as expected. Further, the results show that decreasing the mass inventory gives a better agreement with the experimental data in the lifter, but not in the other reactors. Noteworthy is that the total solid inventory does not change the slope of the pressure profile in the lifter, since this connection operates almost filled with particles, in a very dense regime, and therefore its total mass is primarily determined by the particle density.

### 4.5. Mesh refinement (sub-grid scale effects)
<figure><img src="image_12.png" /><figcaption>Figure 11: Radial profiles of the time-averaged solid vertical velocity in air reactor, fuel reactor and lifter at different heights, depending on the mean particle velocity boundary condition.</figcaption></figure>

The corresponding integrated mass distribution is given in Table 5. The results show that the pressure is overestimated by the lowest particle density in the air reactor. Globally, the mass in the air reactor is higher, and comes mainly from the loop seals and lifter (see Table 5). In the fuel reactor, the pressure is also overestimated by the lowest particle density with respect to the highest density. A greater difference is observed in the lifter where the density makes the slope of the pressure to change. The pressure is better predicted at the bottom than at the top in this case (because of the coupling effect).

In order to investigate the effect of the mass inventory, an additional numerical simulation was carried out. In this simulation the mass was decreased in order to improve the numerical predictions obtained with the lowest density. Comparing the two cases with the same particle density and different mass inventory, it turns out that the smaller the mass, the smaller the pressure at the bottom of the air reactor, as well as in the fuel reactor, as expected. Further, the results show that decreasing the mass inventory gives a better agreement with the experimental data in the lifter, but not in the other reactors. Noteworthy is that the total solid inventory does not change the slope of the pressure profile in the lifter, since this connection operates almost filled with particles, in a very dense regime, and therefore its total mass is primarily determined by the particle density.

4.5. *Mesh refinement (sub-grid scale effects)*
<figure><img src="image_13.png" /><figcaption>Figure 12: Time-averaged pressure depending on the particle density and solid inventory.</figcaption></figure>

Drag modeling is a crucial aspect of the closure assumptions for the accurate prediction of fluidized beds using the Euler-Euler simulation approach. The main issue is the ability of the numerical simulation to take into account solids segregation effects, such as the formation of clusters in circulating fluidized beds, which can occur at very small length scales but with a very strong effect on the macroscopic hydrodynamics, in particular on the entrainment of the solid by the gas flow. Thus, the question of the drag closure model is directly related to the refinement of the mesh towards the characteristic length scale typical of the clustering effect. As pointed out by Igci and Sundaresan (2011b) and Ozel et al. (2013), if the mesh is not sufficiently refined, the Euler-Euler equations of momentum and random kinetic energy have to be supplemented by additional terms accounting for the clustering of particles at the sub-grid scale. The dominant effect is the overestimation of the drag term, which can be corrected using different approaches such as the energy minimization multi-scale (EMMS) approach (Li and Kwauk (1994)), or the sub-grid scale drift velocity modeling (Igci and Sundaresan (2011a); Parmentier et al. (2012); Ozel et al. (2013)). While previous studies have shown that the effect of the sub-grid drag modeling depends on the mesh size and particle characteristics (see, e.g., Wang et al. (2009)), unfortunately there are still no universal dimensionless parameters that allow a priori assessment of sub-grid scale effects and the need to consider sub-grid drag closures. However, it is found that the sub-grid scale effect decreases with the particle inertia and is much less effective for Geldart-B particles, such as those considered in the present study.

In order to estimate the unresolved clustering effect, which is expected to lead to an overesti-
<figure><img src="image_13.png" /><figcaption>Figure 12: Time-averaged pressure depending on the particle density and solid inventory.</figcaption></figure>

Drag modeling is a crucial aspect of the closure assumptions for the accurate prediction of fluidized beds using the Euler-Euler simulation approach. The main issue is the ability of the numerical simulation to take into account solids segregation effects, such as the formation of clusters in circulating fluidized beds, which can occur at very small length scales but with a very strong effect on the macroscopic hydrodynamics, in particular on the entrainment of the solid by the gas flow. Thus, the question of the drag closure model is directly related to the refinement of the mesh towards the characteristic length scale typical of the clustering effect. As pointed out by Igci and Sundaresan (2011b) and Ozel et al. (2013), if the mesh is not sufficiently refined, the Euler-Euler equations of momentum and random kinetic energy have to be supplemented by additional terms accounting for the clustering of particles at the sub-grid scale. The dominant effect is the overestimation of the drag term, which can be corrected using different approaches such as the energy minimization multi-scale (EMMS) approach (Li and Kwauk (1994)), or the sub-grid scale drift velocity modeling (Igci and Sundaresan (2011a); Parmentier et al. (2012); Ozel et al. (2013)). While previous studies have shown that the effect of the sub-grid drag modeling depends on the mesh size and particle characteristics (see, e.g., Wang et al. (2009)), unfortunately there are still no universal dimensionless parameters that allow a priori assessment of sub-grid scale effects and the need to consider sub-grid drag closures. However, it is found that the sub-grid scale effect decreases with the particle inertia and is much less effective for Geldart-B particles, such as those considered in the present study.

In order to estimate the unresolved clustering effect, which is expected to lead to an overesti-
<table><caption>Table 5: Solid mass in the system with different density and/or total inventories (units: kg).</caption><thead><tr><th rowspan="2"></th><th>CASE 1</th><th>CASE 4</th><th>CASE 5</th></tr><tr><th>ρ<sub>p</sub> = 4727 kg/m<sup>3</sup><br>m<sub>total</sub>=125 kg</th><th>ρ<sub>p</sub> = 4062 kg/m<sup>3</sup><br>m<sub>total</sub>=125 kg</th><th>ρ<sub>p</sub> = 4062 kg/m<sup>3</sup><br>m<sub>total</sub>=108 kg</th></tr></thead><tbody><tr><td>Air reactor</td><td>15.656</td><td>23.787</td><td>15.257</td></tr><tr><td>Fuel reactor</td><td>24.597</td><td>25.217</td><td>20.739</td></tr><tr><td>Lifter</td><td>19.698</td><td>16.596</td><td>16.072</td></tr><tr><td>AR loop seal and cyclone</td><td>34.277</td><td>31.329</td><td>28.885</td></tr><tr><td>FR loop seal and cyclone</td><td>30.771</td><td>28.071</td><td>27.047</td></tr><tr><td>Total mass</td><td>125.00</td><td>125.00</td><td>108.00</td></tr></tbody></table>

<table><caption>Table 6: Solid mass in the system depending on the mesh refinement (units: kg).</caption><thead><tr><th rowspan="2"></th><th>CASE 1</th><th>CASE 6</th></tr><tr><th>Reference mesh</th><th>Refined mesh</th></tr></thead><tbody><tr><td>Air reactor</td><td>15.656</td><td>16.508</td></tr><tr><td>Fuel reactor</td><td>24.597</td><td>25.564</td></tr><tr><td>Lifter</td><td>19.698</td><td>19.798</td></tr><tr><td>AR loop seal and cyclone</td><td>34.277</td><td>33.323</td></tr><tr><td>FR loop seal and cyclone</td><td>30.771</td><td>29.807</td></tr><tr><td>Total mass</td><td>125.00</td><td>125.00</td></tr></tbody></table>

mation of the solid entrainment, and to assess the need for a sub-grid drag model, an additional numerical study with a refined mesh was carried out. The numerical simulation was performed by decreasing the mesh size in each direction by a factor of two, in both the air reactor and fuel reactor.

Figure 13 compares the time-averaged pressure predictions obtained using the reference mesh (716 312 cells) and the refined mesh (2 723 176 cells). The corresponding mass distribution in each part of the CLC unit is provided in Table 6.

The results show that the effect of the refinement is low in the fuel reactor and unperceivable in the lifter. Concerning the air reactor, a better agreement between numerical simulations and experiments is observed when a finer mesh is used. However, the difference in pressure predictions is rather small. Additionally, no substantial difference in the entrainment of the solid by the gas flow is found neither in the air reactor nor in the fuel reactor (see, Figure 14 (left)). The mesh refinement leads to a slightly lower bed expansion in the fuel reactor, as expected, but globally the solid volume fraction predictions are very close, as shown in Figure 14 (right). In this figure, the profiles are computed by spatially averaging the local time-averaged solid volume fraction. Finally, comparing the two numerical predictions leads to the conclusions that the sub-grid clustering effect is negligible for such highly inertial particles and that we can overcome the use of a sub-grid drag model in our study.

## 4.6. Analysis of the CLC behavior

Back to the reference case (CASE 1), Figure 15 (left) shows the time evolution of the pressure in the two reactors, at a selected height for each one. At these locations, the instantaneous pressure is averaged in space (over a plane normal to the vertical direction). The dashed/dot-dashed lines represent the mean experimental values. From about 15 seconds, numerical results
<table><caption>Table 5: Solid mass in the system with different density and/or total inventories (units: kg).</caption><thead><tr><th rowspan="2"></th><th>CASE 1</th><th>CASE 4</th><th>CASE 5</th></tr><tr><th>ρ<sub>p</sub> = 4727 kg/m<sup>3</sup><br>m<sub>total</sub>=125 kg</th><th>ρ<sub>p</sub> = 4062 kg/m<sup>3</sup><br>m<sub>total</sub>=125 kg</th><th>ρ<sub>p</sub> = 4062 kg/m<sup>3</sup><br>m<sub>total</sub>=108 kg</th></tr></thead><tbody><tr><th>Air reactor</th><td>15.656</td><td>23.787</td><td>15.257</td></tr><tr><th>Fuel reactor</th><td>24.597</td><td>25.217</td><td>20.739</td></tr><tr><th>Lifter</th><td>19.698</td><td>16.596</td><td>16.072</td></tr><tr><th>AR loop seal and cyclone</th><td>34.277</td><td>31.329</td><td>28.885</td></tr><tr><th>FR loop seal and cyclone</th><td>30.771</td><td>28.071</td><td>27.047</td></tr><tr><th>Total mass</th><td>125.00</td><td>125.00</td><td>108.00</td></tr></tbody></table>

<table><caption>Table 6: Solid mass in the system depending on the mesh refinement (units: kg).</caption><thead><tr><th rowspan="2"></th><th>CASE 1</th><th>CASE 6</th></tr><tr><th>Reference mesh</th><th>Refined mesh</th></tr></thead><tbody><tr><th>Air reactor</th><td>15.656</td><td>16.508</td></tr><tr><th>Fuel reactor</th><td>24.597</td><td>25.564</td></tr><tr><th>Lifter</th><td>19.698</td><td>19.798</td></tr><tr><th>AR loop seal and cyclone</th><td>34.277</td><td>33.323</td></tr><tr><th>FR loop seal and cyclone</th><td>30.771</td><td>29.807</td></tr><tr><th>Total mass</th><td>125.00</td><td>125.00</td></tr></tbody></table>

mation of the solid entrainment, and to assess the need for a sub-grid drag model, an additional numerical study with a refined mesh was carried out. The numerical simulation was performed by decreasing the mesh size in each direction by a factor of two, in both the air reactor and fuel reactor.

Figure 13 compares the time-averaged pressure predictions obtained using the reference mesh (716 312 cells) and the refined mesh (2 723 176 cells). The corresponding mass distribution in each part of the CLC unit is provided in Table 6.

The results show that the effect of the refinement is low in the fuel reactor and unperceivable in the lifter. Concerning the air reactor, a better agreement between numerical simulations and experiments is observed when a finer mesh is used. However, the difference in pressure predictions is rather small. Additionally, no substantial difference in the entrainment of the solid by the gas flow is found neither in the air reactor nor in the fuel reactor (see, Figure 14 (left)). The mesh refinement leads to a slightly lower bed expansion in the fuel reactor, as expected, but globally the solid volume fraction predictions are very close, as shown in Figure 14 (right). In this figure, the profiles are computed by spatially averaging the local time-averaged solid volume fraction. Finally, comparing the two numerical predictions leads to the conclusions that the sub-grid clustering effect is negligible for such highly inertial particles and that we can overcome the use of a sub-grid drag model in our study.

## 4.6. Analysis of the CLC behavior

Back to the reference case (CASE 1), Figure 15 (left) shows the time evolution of the pressure in the two reactors, at a selected height for each one. At these locations, the instantaneous pressure is averaged in space (over a plane normal to the vertical direction). The dashed/dot-dashed lines represent the mean experimental values. From about 15 seconds, numerical results
<figure><img src="image_14.png" /><figcaption>Figure 13: Time-averaged pressure depending on the mesh refinement.</figcaption></figure>

<figure><img src="image_15.png" /><figcaption>Figure 14: Effect of the mesh refinement on the entrainment (left). Time-averaged vertical profile of the solid volume fraction depending on the mesh refinement (right).</figcaption></figure>

stabilize around a constant mean value in both the fuel and air reactors. At the selected location, the predicted pressure is very close to the experimental data in the air reactor, while it is slightly lower in the fuel reactor, which is consistent with the time-averaged profiles previously analyzed. Further, in the fuel reactor the pressure fluctuates violently around a mean value. This confirms that in this dense region, at the bottom, the fuel reactor operates rather like a bubbling fluidized
<figure><img src="image_14.png" /><figcaption>Figure 13: Time-averaged pressure depending on the mesh refinement.</figcaption></figure>

<figure><img src="image_15.png" /><figcaption>Figure 14: Effect of the mesh refinement on the entrainment (left). Time-averaged vertical profile of the solid volume fraction depending on the mesh refinement (right).</figcaption></figure>

stabilize around a constant mean value in both the fuel and air reactors. At the selected location, the predicted pressure is very close to the experimental data in the air reactor, while it is slightly lower in the fuel reactor, which is consistent with the time-averaged profiles previously analyzed. Further, in the fuel reactor the pressure fluctuates violently around a mean value. This confirms that in this dense region, at the bottom, the fuel reactor operates rather like a bubbling fluidized
<figure><img src="image_16.png" /><figcaption>Figure 15: Time evolutions of pressure in the fuel and air reactors (left). Profiles of the time-averaged pressure in the different parts of the CLC pilot (right).</figcaption></figure>

bed.

Figure 15 (right) shows the time-averaged pressure distribution along the height, according to the location in the CLC system. There is an obvious decrease of pressure with height in the lifter, and the same is true also for the fuel and air reactors. All pressure measurements depend on the amount of particle loading. The pressure balance of the current interconnected reactor system reveals that pressure is largest in the lifter and smallest in the FR cyclone. The coupling of the different unit components is clearly identified by the figure. One can recognize the connections between the bottom of the fuel reactor and the lifter, and between the top of the lifter and the air reactor, as shown in Figure 2. The pressure is the same at each connecting location. The pressure distribution therefore depends on the mass inventory in each part of the CLC unit, but also on the coupling effect of the entire system. A change in pressure that occurs in one part will lead to a pressure modification throughout the whole system.

The mass flow rates obtained from the numerical simulations are displayed in Figure 16 (left). The mass flow rate of the solid leaving the air reactor fluctuates around a mean value, which is close to the value expected from the experiments. Results confirm that the air reactor operates as a circulating fluidized bed and that a substantial quantity of oxygen carrier leaves the air reactor from the top, according to the CLC concept and design. The fuel reactor was expected to operate in a mixed regime with most of the oxygen carriers entering the air reactor through the lifter. This regime is consistent with the profile of the mean pressure in the fuel reactor. Results show however that about half of the solid is transported from the fuel reactor to the air reactor through the lifter, while about half leaves the fuel reactor from the top and enters the air reactor through the corresponding cyclone and loop seal. This amount is higher than expected based on the design values of the experiments (~30%). This point will be investigated in the future, under reactive conditions.

In the numerical simulation, the mass of solid in each part of the system was initialized using values estimated from the pressure drop measurements in the experiments. In order to check the accuracy of such an estimation method, the time evolutions of the mass obtained from the pressure in the numerical simulation, in the two reactors, are plotted and displayed in Figure 16 (right). The time averaged results should be compared with the mean values given in Table 6 (CASE 1), computed by a volume integral using the time-averaged solid volume fraction, together with the constant particle density. The results show that the solid mass is
<figure><img src="image_16.png" /><figcaption>Figure 15: Time evolutions of pressure in the fuel and air reactors (left). Profiles of the time-averaged pressure in the different parts of the CLC pilot (right).</figcaption></figure>

bed.

Figure 15 (right) shows the time-averaged pressure distribution along the height, according to the location in the CLC system. There is an obvious decrease of pressure with height in the lifter, and the same is true also for the fuel and air reactors. All pressure measurements depend on the amount of particle loading. The pressure balance of the current interconnected reactor system reveals that pressure is largest in the lifter and smallest in the FR cyclone. The coupling of the different unit components is clearly identified by the figure. One can recognize the connections between the bottom of the fuel reactor and the lifter, and between the top of the lifter and the air reactor, as shown in Figure 2. The pressure is the same at each connecting location. The pressure distribution therefore depends on the mass inventory in each part of the CLC unit, but also on the coupling effect of the entire system. A change in pressure that occurs in one part will lead to a pressure modification throughout the whole system.

The mass flow rates obtained from the numerical simulations are displayed in Figure 16 (left). The mass flow rate of the solid leaving the air reactor fluctuates around a mean value, which is close to the value expected from the experiments. Results confirm that the air reactor operates as a circulating fluidized bed and that a substantial quantity of oxygen carrier leaves the air reactor from the top, according to the CLC concept and design. The fuel reactor was expected to operate in a mixed regime with most of the oxygen carriers entering the air reactor through the lifter. This regime is consistent with the profile of the mean pressure in the fuel reactor. Results show however that about half of the solid is transported from the fuel reactor to the air reactor through the lifter, while about half leaves the fuel reactor from the top and enters the air reactor through the corresponding cyclone and loop seal. This amount is higher than expected based on the design values of the experiments (~30%). This point will be investigated in the future, under reactive conditions.

In the numerical simulation, the mass of solid in each part of the system was initialized using values estimated from the pressure drop measurements in the experiments. In order to check the accuracy of such an estimation method, the time evolutions of the mass obtained from the pressure in the numerical simulation, in the two reactors, are plotted and displayed in Figure 16 (right). The time averaged results should be compared with the mean values given in Table 6 (CASE 1), computed by a volume integral using the time-averaged solid volume fraction, together with the constant particle density. The results show that the solid mass is
<figure><img src="image_17.png" /><figcaption>Figure 16: Time evolution of the mass flow rates of solid leaving the reactors from their top outlets (left). Time evolution of the solid mass in the air and fuel reactors (right).</figcaption></figure>

overestimated, especially in the fuel reactor, when computed from the pressure measurements. Figure 16 (right) also shows that the solid mass in the fuel reactor stabilizes quite soon around a mean value close to the initial one. In contrast, after an initial increase, the solid mass in the air reactor decreases with time, and reaches a steady state only after 15 seconds of simulation. Further, results also show larger fluctuations in the fuel reactor compared to the air reactor, while the frequency is quite similar.

Figure 17 shows the time evolution of both the gas and solid velocities at a given height in the two reactors. Under the current conditions, the gas velocity is approximately 35-40% greater than the solid velocity in both the fuel and air reactors. Further, according to the numerical predictions, velocities in the air reactor are greater than in the fuel reactor, as expected. Both the gas and solid velocities fluctuate wildly due to the intense interaction between the two phases. In the experiments, the gas velocity was measured at the fuel reactor exit. A mean value is therefore available from the experiments for comparison. Figure 17 (right) shows that the numerical prediction matches well the experimental result.

Radial profiles of the time-averaged solid volume fraction in the two reactors and in the lifter are shown in Figure 18 at different heights. Values are plotted on a line through the center of the reactor in a radial direction. Profiles extend differently in the reactors due to the conical structure at the bottom of each. Radial coordinates are normalized by $R$, which is the maximum radius of the corresponding reactor. For the air reactor it corresponds to $R = 0.115$ m (the radius ranges from 0.077 m (at $H = 0$ m) to 0.115 m (at $H > 1$ m)). For the fuel reactor, $R = 0.077$ m (the radius spanning from 0.05 m (at $H = 0$ m) to 0.077 m (at $H > 1$ m)). For the lifter, two cylindrical zones are gradually connected by a conformal mesh. The radius of each is 0.051 m (at $-0.61 < H \le 0$ m) and 0.039 m (at $H > 0$ m).

The air reactor exhibits the well-known core-annulus flow structure almost at all locations, corresponding to accumulation of particles near the wall and a more dilute regime in the center. In the air reactor, the difference between the solid volume fractions at different heights is quite small, except at the wall. Profiles are not symmetrical close to the injections. In the air reactor, the secondary gas injections are located at different heights (around $H = 0.5$ m and $H = 0.95$ m). Also the connecting parts between the air reactor and lifter ($H = 0.9$ m) affects the velocity distribution. The fuel reactor behavior is closer to that of a dense fluidized bed, but with lower solid volume fractions. The fuel reactor is indeed working in a mixed regime,
<figure><img src="image_17.png" /><figcaption>Figure 16: Time evolution of the mass flow rates of solid leaving the reactors from their top outlets (left). Time evolution of the solid mass in the air and fuel reactors (right).</figcaption></figure>

overestimated, especially in the fuel reactor, when computed from the pressure measurements. Figure 16 (right) also shows that the solid mass in the fuel reactor stabilizes quite soon around a mean value close to the initial one. In contrast, after an initial increase, the solid mass in the air reactor decreases with time, and reaches a steady state only after 15 seconds of simulation. Further, results also show larger fluctuations in the fuel reactor compared to the air reactor, while the frequency is quite similar.

Figure 17 shows the time evolution of both the gas and solid velocities at a given height in the two reactors. Under the current conditions, the gas velocity is approximately 35-40% greater than the solid velocity in both the fuel and air reactors. Further, according to the numerical predictions, velocities in the air reactor are greater than in the fuel reactor, as expected. Both the gas and solid velocities fluctuate wildly due to the intense interaction between the two phases. In the experiments, the gas velocity was measured at the fuel reactor exit. A mean value is therefore available from the experiments for comparison. Figure 17 (right) shows that the numerical prediction matches well the experimental result.

Radial profiles of the time-averaged solid volume fraction in the two reactors and in the lifter are shown in Figure 18 at different heights. Values are plotted on a line through the center of the reactor in a radial direction. Profiles extend differently in the reactors due to the conical structure at the bottom of each. Radial coordinates are normalized by $R$, which is the maximum radius of the corresponding reactor. For the air reactor it corresponds to $R = 0.115$ m (the radius ranges from 0.077 m (at $H = 0$ m) to 0.115 m (at $H > 1$ m)). For the fuel reactor, $R = 0.077$ m (the radius spanning from 0.05 m (at $H = 0$ m) to 0.077 m (at $H > 1$ m)). For the lifter, two cylindrical zones are gradually connected by a conformal mesh. The radius of each is 0.051 m (at $-0.61 < H \le 0$ m) and 0.039 m (at $H > 0$ m).

The air reactor exhibits the well-known core-annulus flow structure almost at all locations, corresponding to accumulation of particles near the wall and a more dilute regime in the center. In the air reactor, the difference between the solid volume fractions at different heights is quite small, except at the wall. Profiles are not symmetrical close to the injections. In the air reactor, the secondary gas injections are located at different heights (around $H = 0.5$ m and $H = 0.95$ m). Also the connecting parts between the air reactor and lifter ($H = 0.9$ m) affects the velocity distribution. The fuel reactor behavior is closer to that of a dense fluidized bed, but with lower solid volume fractions. The fuel reactor is indeed working in a mixed regime,
<figure><img src="image_18.png" /><figcaption>Figure 17: Time evolution of the solid and gas velocities in the air and fuel reactors.</figcaption></figure>

bubbling and circulating, at these operating conditions. The lifter transports the particles from the fuel reactor to the air reactor by an overall upward movement, operating with high solid concentrations, especially close to the wall.

<figure><img src="image_19.png" /><figcaption>Figure 18: Radial profiles of the time-averaged solid volume fraction in the air reactor, fuel reactor and lifter, at different heights.</figcaption></figure>
<figure><img src="image_18.png" /><figcaption>Figure 17: Time evolution of the solid and gas velocities in the air and fuel reactors.</figcaption></figure>

bubbling and circulating, at these operating conditions. The lifter transports the particles from the fuel reactor to the air reactor by an overall upward movement, operating with high solid concentrations, especially close to the wall.

<figure><img src="image_19.png" /><figcaption>Figure 18: Radial profiles of the time-averaged solid volume fraction in the air reactor, fuel reactor and lifter, at different heights.</figcaption></figure>
Figure 19 (top left) shows a scatterplot of the correlation coefficient used for modeling correlated and uncorrelated contributions of the particle kinetic energy in the frame of the correlated collision model. As shown in Eq. (46), such a coefficient is related to the ratio between $q_{gs}$ and $\sqrt{4kq_s^2}$. From the figure we can observe that most of the instantaneous values are located in the range from $10^{-5}$ to $10^{-1}$. The value of the gas-particle correlation coefficient, $\zeta_{gs}$, is far less than 1, which leads to a very low correlation effect between neighboring particles due to their inertia with respect to the gas turbulent flow. As a result, the inter-particle collision time for correlated and uncorrelated model is nearly the same. Therefore, for the current case, correlated and uncorrelated models produce similar predictions.

<figure><img src="image_20.png" /><figcaption>Figure 19: Instantaneous fluid-particle correlation coefficient, particle fluctuant kinetic energy and inter-particle collision time versus the solid volume fraction.</figcaption></figure>

A scatterplot of the particle random kinetic energy versus the solid volume fraction is shown in Figure 19 (top right). Results show that in the air reactor particles are more fluctuating than in the fuel reactor, and much more than in the lifter where most of the movement is represented by a collective transport. The inter-particle collision time is shown in Figure 19 (bottom). Its dependency on the solid volume fraction is inherent to the model. Results additionally show that for a given value of the solid volume fraction, the collision time is lower in the air reactor where agitation is larger, i.e. collision frequency is higher in the part of the CLC corresponding to stronger particle fluctuations. The inter-particle collision time takes large values at the maximum compaction because the particle random kinetic energy tends to zero in such zones (see Figure 19 (top right)).
Figure 19 (top left) shows a scatterplot of the correlation coefficient used for modeling correlated and uncorrelated contributions of the particle kinetic energy in the frame of the correlated collision model. As shown in Eq. (46), such a coefficient is related to the ratio between $q_{gs}$ and $\sqrt{4kq_s^2}$. From the figure we can observe that most of the instantaneous values are located in the range from $10^{-5}$ to $10^{-1}$. The value of the gas-particle correlation coefficient, $\zeta_{gs}$, is far less than 1, which leads to a very low correlation effect between neighboring particles due to their inertia with respect to the gas turbulent flow. As a result, the inter-particle collision time for correlated and uncorrelated model is nearly the same. Therefore, for the current case, correlated and uncorrelated models produce similar predictions.

<figure><img src="image_20.png" /><figcaption>Figure 19: Instantaneous fluid-particle correlation coefficient, particle fluctuant kinetic energy and inter-particle collision time versus the solid volume fraction.</figcaption></figure>

A scatterplot of the particle random kinetic energy versus the solid volume fraction is shown in Figure 19 (top right). Results show that in the air reactor particles are more fluctuating than in the fuel reactor, and much more than in the lifter where most of the movement is represented by a collective transport. The inter-particle collision time is shown in Figure 19 (bottom). Its dependency on the solid volume fraction is inherent to the model. Results additionally show that for a given value of the solid volume fraction, the collision time is lower in the air reactor where agitation is larger, i.e. collision frequency is higher in the part of the CLC corresponding to stronger particle fluctuations. The inter-particle collision time takes large values at the maximum compaction because the particle random kinetic energy tends to zero in such zones (see Figure 19 (top right)).
Two additional quantities are examined, which are the ratio between the gas-particle velocity covariance, $q_{gs}$, and twice the gas ($k$) or particle ($q_s^2$) kinetic energies. These ratios are relevant in the interpretation of the flow behavior. Results are given in Figure 20. A first information is obtained looking at the term in Eq. (28), which represents the effect of the interphase kinetic energy exchange on the evolution of the particle random kinetic energy (Eq. (22)). Figure 20 (right) shows that in all the relevant parts of the CLC the ratio $q_{gs}/2q_s^2$ is generally smaller than unity, and even smaller in the air reactor, except in very dense zones where the solid volume fraction tends to the maximum compaction. In this case $q_s^2$ is very small because of the larger dissipation in such zones. A ratio $q_{gs}/2q_s^2$ smaller than unity means that the term in Eq. (28) is a negative quantity, i.e. the particle agitation is not due to the entrainment by the turbulence, which on the contrary acts at dissipating the particle fluctuations. The same destruction effect is found in the balance equation of the gas turbulent kinetic energy (Eq. (11)), based on the interphase coupling term (Eq. (13)) and the results of Figure 20 (left). In this case, the effect of the scalar product of the drift with the relative velocity is found to be lower than the other contributions in the coupling term. More complicated is instead the interpretation of the results on the evolution of the covariance itself (Eq. (29)). Looking at the source term in Eq. (31), it comes out that the first contribution is positive. The second one depends instead on the intensity of the correlated part of the particle kinetic energy, which may be related to the ratio between the gas-particle velocity covariance and the gas kinetic energy as follow: $(q_{sg} - 2\tilde{q}_s^2) = q_{sg}(1 - q_{sg}/2k)$ (Section 2). If $q_{sg} < 2k$, a negative sign can be anticipated, which means that the second term in Eq. (31) acts at dissipating the gas-particle velocity covariance. Globally, the source term corresponding to Eq. (31) is a destruction term when $\alpha_s\rho_s/\alpha_g\rho_g > 1$.

Finally, Figure 21 shows snapshots of some relevant quantities of the gas-particle flow. Some of the most important information is that the gas turbulent kinetic energy is largely dissipated by the particle two-way coupling effect and is much smaller than the particle fluctuant kinetic energy. The gas turbulent viscosity is found much lower than its laminar counterpart ($1.7 \times 10^{-4} \text{ m}^2/\text{s}$) revealing that, in the CLC at the current operating conditions, the gas turbulence predicted by the $k-\epsilon$ model has no effect on the gas flow prediction. In addition, the correlation coefficient based on the fluid-particle velocity covariance (Eq. (46)) is very small showing that both gas and particle fluctuating velocities are uncorrelated. As a consequence, the proposed correlated model predicts that, in such a flow configuration, the random velocities of neighboring discrete particles are largely uncorrelated and the total predicted particle fluctuant kinetic energy may be recognized as the granular temperature.

# 5. Conclusion

A model based on the Euler-Euler approach is adopted in this study to predict the hydrodynamic behavior of a chemical looping combustion system. Three-dimensional unsteady numerical simulations of a 150 kW$_{th}$ pilot (operating at SINTEF, Trondheim, Norway) were carried out using NEPTUNE_CFD, with the main goal to gain insight in the local and instantaneous flow behavior and operating characteristics. In the original experiments, the CLC pilot operated with ilmenite as oxygen carrier and biomass (wood pellet) as fuel. In this numerical study, biomass was not considered as an additional solid phase and gases from biomass conversion and redox reactions were accounted for by adjusting the injection conditions. Moreover, an isothermal flow was assumed, since the 150 kW$_{th}$ CLC system operates in almost uniform temperature conditions, according to the experiments. The numerical geometry was built according to the experimental facility and was discretized by using a numerical mesh corresponding to a
Two additional quantities are examined, which are the ratio between the gas-particle velocity covariance, $q_{gs}$, and twice the gas ($k$) or particle ($q_s^2$) kinetic energies. These ratios are relevant in the interpretation of the flow behavior. Results are given in Figure 20. A first information is obtained looking at the term in Eq. (28), which represents the effect of the interphase kinetic energy exchange on the evolution of the particle random kinetic energy (Eq. (22)). Figure 20 (right) shows that in all the relevant parts of the CLC the ratio $q_{gs}/2q_s^2$ is generally smaller than unity, and even smaller in the air reactor, except in very dense zones where the solid volume fraction tends to the maximum compaction. In this case $q_s^2$ is very small because of the larger dissipation in such zones. A ratio $q_{gs}/2q_s^2$ smaller than unity means that the term in Eq. (28) is a negative quantity, i.e. the particle agitation is not due to the entrainment by the turbulence, which on the contrary acts at dissipating the particle fluctuations. The same destruction effect is found in the balance equation of the gas turbulent kinetic energy (Eq. (11)), based on the interphase coupling term (Eq. (13)) and the results of Figure 20 (left). In this case, the effect of the scalar product of the drift with the relative velocity is found to be lower than the other contributions in the coupling term. More complicated is instead the interpretation of the results on the evolution of the covariance itself (Eq. (29)). Looking at the source term in Eq. (31), it comes out that the first contribution is positive. The second one depends instead on the intensity of the correlated part of the particle kinetic energy, which may be related to the ratio between the gas-particle velocity covariance and the gas kinetic energy as follow: $(q_{sg} - 2\tilde{q}_s^2) = q_{sg}(1 - q_{sg}/2k)$ (Section 2). If $q_{sg} < 2k$, a negative sign can be anticipated, which means that the second term in Eq. (31) acts at dissipating the gas-particle velocity covariance. Globally, the source term corresponding to Eq. (31) is a destruction term when $\alpha_s \rho_s / \alpha_g \rho_g > 1$.

Finally, Figure 21 shows snapshots of some relevant quantities of the gas-particle flow. Some of the most important information is that the gas turbulent kinetic energy is largely dissipated by the particle two-way coupling effect and is much smaller than the particle fluctuant kinetic energy. The gas turbulent viscosity is found much lower than its laminar counterpart ($1.7 \times 10^{-4} \text{ m}^2/\text{s}$) revealing that, in the CLC at the current operating conditions, the gas turbulence predicted by the $k-\epsilon$ model has no effect on the gas flow prediction. In addition, the correlation coefficient based on the fluid-particle velocity covariance (Eq. (46)) is very small showing that both gas and particle fluctuating velocities are uncorrelated. As a consequence, the proposed correlated model predicts that, in such a flow configuration, the random velocities of neighboring discrete particles are largely uncorrelated and the total predicted particle fluctuant kinetic energy may be recognized as the granular temperature.

# 5. Conclusion

A model based on the Euler-Euler approach is adopted in this study to predict the hydrodynamic behavior of a chemical looping combustion system. Three-dimensional unsteady numerical simulations of a 150 kW$_{th}$ pilot (operating at SINTEF, Trondheim, Norway) were carried out using NEPTUNE_CFD, with the main goal to gain insight in the local and instantaneous flow behavior and operating characteristics. In the original experiments, the CLC pilot operated with ilmenite as oxygen carrier and biomass (wood pellet) as fuel. In this numerical study, biomass was not considered as an additional solid phase and gases from biomass conversion and redox reactions were accounted for by adjusting the injection conditions. Moreover, an isothermal flow was assumed, since the 150 kW$_{th}$ CLC system operates in almost uniform temperature conditions, according to the experiments. The numerical geometry was built according to the experimental facility and was discretized by using a numerical mesh corresponding to a
<figure><img src="image_21.png" /></figure>

<figure><img src="image_22.png" /></figure>

<figure><img src="image_23.png" /><figcaption>Figure 20: Instantaneous ratio between the gas-particle velocity covariance and twice the turbulent kinetic energy (on the left) or particle kinetic energy (on the right), versus the solid volume fraction.</figcaption></figure>
<figure><img src="image_21.png" /></figure>

<figure><img src="image_22.png" /></figure>

<figure><img src="image_23.png" /><figcaption>Figure 20: Instantaneous ratio between the gas-particle velocity covariance and twice the turbulent kinetic energy (on the left) or particle kinetic energy (on the right), versus the solid volume fraction.</figcaption></figure>
<figure><img src="image_24.png" /><figcaption>Figure 21: Instantaneous visualization of gas kinetic energy, gas turbulent viscosity, particle fluctuant kinetic energy and gas-particle velocity covariance.</figcaption></figure>

suitable compromise between fine and coarse meshes, considering both the accuracy and computational costs. Results about the pressure in the different parts of the pilot showed a general agreement between numerical predictions and experimental data, proving that the simplifying assumptions considered in this study allow to reproduce satisfactorily the flow regime. The hydrodynamics of the process was therefore investigated in detail, in particular studying the solid mass flow rates, the gas and solid velocities and the particle distribution in the relevant parts of the CLC system. Numerical simulations showed that the air reactor operates in a circulating bed regime, while the fuel reactor works in a mixed regime, in between a dense and a circulating fluidized bed. Numerical simulations also showed that the gas turbulence is negligible at this operating condition and weakly correlated with the particle fluctuating motion. So, according to the modeling approach, agitation between neighboring particles was found rather uncorrelated and for these reasons, both the uncorrelated and correlated collision models led to almost the same results. The effects of the two limit-case wall boundary conditions (free-slip and no-slip) for the mean particle velocity were also analyzed. According to the solid
<figure><img src="image_24.png" /><figcaption>Figure 21: Instantaneous visualization of gas kinetic energy, gas turbulent viscosity, particle fluctuant kinetic energy and gas-particle velocity covariance.</figcaption></figure>

suitable compromise between fine and coarse meshes, considering both the accuracy and computational costs. Results about the pressure in the different parts of the pilot showed a general agreement between numerical predictions and experimental data, proving that the simplifying assumptions considered in this study allow to reproduce satisfactorily the flow regime. The hydrodynamics of the process was therefore investigated in detail, in particular studying the solid mass flow rates, the gas and solid velocities and the particle distribution in the relevant parts of the CLC system. Numerical simulations showed that the air reactor operates in a circulating bed regime, while the fuel reactor works in a mixed regime, in between a dense and a circulating fluidized bed. Numerical simulations also showed that the gas turbulence is negligible at this operating condition and weakly correlated with the particle fluctuating motion. So, according to the modeling approach, agitation between neighboring particles was found rather uncorrelated and for these reasons, both the uncorrelated and correlated collision models led to almost the same results. The effects of the two limit-case wall boundary conditions (free-slip and no-slip) for the mean particle velocity were also analyzed. According to the solid
circulation, it was found that a no-slip condition in the air reactor leads to an increase of the global circulation rate. The reason is not completely understood. The asymmetry of the solid velocity radial profiles in the air reactor makes the back-mixing analysis inconclusive at this stage. Further studies are needed to clarify this point. The results however suggested that a no-slip boundary condition can be considered satisfactory in a dense regime, but its use should be avoided in dilute zones, such as in the air reactor. The question of the wall boundary conditions for the solid phase is an important point that deserves to be investigated further, and it is left as a future work. Indeed, more appropriate boundary conditions should be used to represent the behavior of the different particle-wall interactions in the presence of both dense and dilute regimes. Globally, the current study assessed satisfactorily the isothermal, non-reactive modeling approach regarding the hydrodynamic predictions of a reactive unit. This allows the design phase to deal with the reactive aspects at a later time.

## Acknowledgment

This work is supported by Chinese-European Emission-Reducing Solutions (CHEERS) under the European Union's Horizon 2020 research and innovation program (No 764697). For the numerical simulations at IMFT, it was granted access to the HPC resources of CALMIP supercomputing center under the allocation P19017 and CINES supercomputing center under the allocation A0082B10864. CALMIP and CINES are gratefully acknowledged. The authors wish to thank Ing. Hervé Neau (CoSiNus) for his help with the NEPTUNE_CFD code.

## Nomenclature

### Latin Symbols

<table><tr><td>C<sub>D</sub></td><td>drag coefficient</td></tr><tr><td>d<sub>s</sub></td><td>particle diameter</td></tr><tr><td>D<sub>s,ij</sub></td><td>particle shear tensor</td></tr><tr><td>e<sub>c</sub></td><td>normal restitution coefficient</td></tr><tr><td>g</td><td>gravity</td></tr><tr><td>g<sub>0</sub></td><td>radial distribution function</td></tr><tr><td>I</td><td>interphase momentum transfer</td></tr><tr><td>k</td><td>gas turbulent kinetic energy</td></tr><tr><td>P</td><td>pressure</td></tr><tr><td>P<sup>fr</sup><sub>s</sub></td><td>frictional pressure</td></tr><tr><td>q<sub>gs</sub></td><td>fluid-particle velocity covariance</td></tr><tr><td>q<sup>2</sup><sub>s</sub></td><td>particle fluctuant kinetic energy</td></tr><tr><td>q̃<sup>2</sup><sub>s</sub></td><td>correlated particle kinetic energy</td></tr><tr><td>R<sub>g,ij</sub></td><td>turbulent-Reynolds stress tensor</td></tr><tr><td>Re<sub>p</sub></td><td>particle Reynolds number</td></tr><tr><td>R<sub>s,ij</sub></td><td>particle kinetic stress tensor</td></tr><tr><td>u''</td><td>velocity fluctuation</td></tr><tr><td>U</td><td>mean velocity</td></tr><tr><td>v<sub>r</sub></td><td>instantaneous relative velocity</td></tr><tr><td>V<sub>r</sub></td><td>relative velocity</td></tr><tr><td>V<sub>d</sub></td><td>drift velocity</td></tr></table>

### Greek letters
circulation, it was found that a no-slip condition in the air reactor leads to an increase of the global circulation rate. The reason is not completely understood. The asymmetry of the solid velocity radial profiles in the air reactor makes the back-mixing analysis inconclusive at this stage. Further studies are needed to clarify this point. The results however suggested that a no-slip boundary condition can be considered satisfactory in a dense regime, but its use should be avoided in dilute zones, such as in the air reactor. The question of the wall boundary conditions for the solid phase is an important point that deserves to be investigated further, and it is left as a future work. Indeed, more appropriate boundary conditions should be used to represent the behavior of the different particle-wall interactions in the presence of both dense and dilute regimes. Globally, the current study assessed satisfactorily the isothermal, non-reactive modeling approach regarding the hydrodynamic predictions of a reactive unit. This allows the design phase to deal with the reactive aspects at a later time.

## Acknowledgment

This work is supported by Chinese-European Emission-Reducing Solutions (CHEERS) under the European Union's Horizon 2020 research and innovation program (No 764697). For the numerical simulations at IMFT, it was granted access to the HPC resources of CALMIP supercomputing center under the allocation P19017 and CINES supercomputing center under the allocation A0082B10864. CALMIP and CINES are gratefully acknowledged. The authors wish to thank Ing. Hervé Neau (CoSiNus) for his help with the NEPTUNE_CFD code.

## Nomenclature

### Latin Symbols

<table><tr><td>C<sub>D</sub></td><td>drag coefficient</td></tr><tr><td>d<sub>s</sub></td><td>particle diameter</td></tr><tr><td>D<sub>s,ij</sub></td><td>particle shear tensor</td></tr><tr><td>e<sub>c</sub></td><td>normal restitution coefficient</td></tr><tr><td>g</td><td>gravity</td></tr><tr><td>g<sub>0</sub></td><td>radial distribution function</td></tr><tr><td>I</td><td>interphase momentum transfer</td></tr><tr><td>k</td><td>gas turbulent kinetic energy</td></tr><tr><td>P</td><td>pressure</td></tr><tr><td>P<sup>fr</sup><sub>s</sub></td><td>frictional pressure</td></tr><tr><td>q<sub>gs</sub></td><td>fluid-particle velocity covariance</td></tr><tr><td>q<sup>2</sup><sub>s</sub></td><td>particle fluctuant kinetic energy</td></tr><tr><td>q̃<sup>2</sup><sub>s</sub></td><td>correlated particle kinetic energy</td></tr><tr><td>R<sub>g,ij</sub></td><td>turbulent-Reynolds stress tensor</td></tr><tr><td>Re<sub>p</sub></td><td>particle Reynolds number</td></tr><tr><td>R<sub>s,ij</sub></td><td>particle kinetic stress tensor</td></tr><tr><td>u''</td><td>velocity fluctuation</td></tr><tr><td>U</td><td>mean velocity</td></tr><tr><td>v<sub>r</sub></td><td>instantaneous relative velocity</td></tr><tr><td>V<sub>r</sub></td><td>relative velocity</td></tr><tr><td>V<sub>d</sub></td><td>drift velocity</td></tr></table>

### Greek letters
<table><tr><td>α</td><td>volume fraction</td></tr><tr><td>δ<sub>ij</sub></td><td>Kronecker symbol</td></tr><tr><td>δq<sup>2</sup><sub>s</sub></td><td>uncorrelated contribution of the particle kinetic energy</td></tr><tr><td>ε</td><td>gas turbulent dissipation rate</td></tr><tr><td>ε<sub>gs</sub></td><td>fluid-particle covariance dissipation rate</td></tr><tr><td>ζ<sup>2</sup><sub>gs</sub></td><td>correlation coefficient</td></tr><tr><td>η</td><td>internal friction angle</td></tr><tr><td>Θ<sub>s</sub></td><td>granular temperature</td></tr><tr><td>Θ<sub>g,ij</sub></td><td>viscous stress tensor</td></tr><tr><td>Θ<sub>s,ij</sub></td><td>collisional stress tensor</td></tr><tr><td>κ<sup>eff</sup><sub>s</sub></td><td>particle effective diffusivity</td></tr><tr><td>κ<sup>kin</sup><sub>s</sub></td><td>particle kinetic diffusivity</td></tr><tr><td>κ<sup>col</sup><sub>s</sub></td><td>particle collisional diffusivity</td></tr><tr><td>μ<sub>g</sub></td><td>laminar dynamic viscosity</td></tr><tr><td>μ<sup>fr</sup><sub>s</sub></td><td>frictional viscosity</td></tr><tr><td>ν<sup>t</sup><sub>g</sub></td><td>turbulent kinematic viscosity</td></tr><tr><td>ν<sup>kin</sup><sub>s</sub></td><td>particle kinetic viscosity</td></tr><tr><td>ν<sup>col</sup><sub>s</sub></td><td>particle collisional viscosity</td></tr><tr><td>ν<sup>t</sup><sub>gs</sub></td><td>turbulent gas-particle viscosity</td></tr><tr><td>Π<sub>qgs</sub></td><td>interphase gas-particle covariance interaction term</td></tr><tr><td>Π<sub>qs</sub></td><td>interphase turbulent kinetic energy transfer rate</td></tr><tr><td>Π<sup>k</sup><sub>s→g</sub></td><td>interphase turbulent kinetic energy interaction term</td></tr><tr><td>Π<sup>ε</sup><sub>s→g</sub></td><td>interphase turbulent dissipation rate interaction term</td></tr><tr><td>ρ</td><td>density</td></tr><tr><td>Σ<sub>ij</sub></td><td>stress tensor</td></tr><tr><td>τ<sup>t</sup><sub>g</sub></td><td>fluid turbulent time scale</td></tr><tr><td>τ<sup>t</sup><sub>gs</sub></td><td>eddy-particle interaction time</td></tr><tr><td>τ<sup>F</sup><sub>gs</sub></td><td>mean particle relaxation time</td></tr><tr><td>τ<sup>c</sup><sub>s</sub></td><td>inter-particle collision time</td></tr><tr><td>φ<sub>s,ij</sub></td><td>frictional tensor</td></tr></table>

Abbreviation

<table><tr><td>AR</td><td>air reactor</td></tr><tr><td>CFB</td><td>circulating fluidized bed</td></tr><tr><td>CFD</td><td>computational fluid dynamics</td></tr><tr><td>CLC</td><td>chemical looping combustion</td></tr><tr><td>CLOU</td><td>chemical looping with oxygen uncoupling</td></tr><tr><td>DEM</td><td>Discrete Element Method</td></tr><tr><td>EMMS</td><td>energy minimization multi-scale</td></tr><tr><td>FR</td><td>fuel reactor</td></tr><tr><td>HPC</td><td>high performance computing</td></tr><tr><td>OC</td><td>oxygen carrier</td></tr></table>

724 References

725 Abad, A., Adánez, J., Cuadrat, A., García-Labiano, F., Gayán, P., de Diego, L. F., 2011.

726 Kinetics of redox reactions of ilmenite for chemical-looping combustion. Chemical Engineering

727 Science 66 (4), 689–702.
<table><tr><td>α</td><td>volume fraction</td></tr><tr><td>δ<sub>ij</sub></td><td>Kronecker symbol</td></tr><tr><td>δq<sup>2</sup><sub>s</sub></td><td>uncorrelated contribution of the particle kinetic energy</td></tr><tr><td>ε</td><td>gas turbulent dissipation rate</td></tr><tr><td>ε<sub>gs</sub></td><td>fluid-particle covariance dissipation rate</td></tr><tr><td>ζ<sup>2</sup><sub>gs</sub></td><td>correlation coefficient</td></tr><tr><td>η</td><td>internal friction angle</td></tr><tr><td>Θ<sub>s</sub></td><td>granular temperature</td></tr><tr><td>Θ<sub>g,ij</sub></td><td>viscous stress tensor</td></tr><tr><td>Θ<sub>s,ij</sub></td><td>collisional stress tensor</td></tr><tr><td>κ<sup>eff</sup><sub>s</sub></td><td>particle effective diffusivity</td></tr><tr><td>κ<sup>kin</sup><sub>s</sub></td><td>particle kinetic diffusivity</td></tr><tr><td>κ<sup>col</sup><sub>s</sub></td><td>particle collisional diffusivity</td></tr><tr><td>μ<sub>g</sub></td><td>laminar dynamic viscosity</td></tr><tr><td>μ<sup>fr</sup><sub>s</sub></td><td>frictional viscosity</td></tr><tr><td>ν<sup>t</sup><sub>g</sub></td><td>turbulent kinematic viscosity</td></tr><tr><td>ν<sup>kin</sup><sub>s</sub></td><td>particle kinetic viscosity</td></tr><tr><td>ν<sup>col</sup><sub>s</sub></td><td>particle collisional viscosity</td></tr><tr><td>ν<sup>t</sup><sub>gs</sub></td><td>turbulent gas-particle viscosity</td></tr><tr><td>Π<sub>qgs</sub></td><td>interphase gas-particle covariance interaction term</td></tr><tr><td>Π<sub>qs</sub></td><td>interphase turbulent kinetic energy transfer rate</td></tr><tr><td>Π<sup>k</sup><sub>s→g</sub></td><td>interphase turbulent kinetic energy interaction term</td></tr><tr><td>Π<sup>ε</sup><sub>s→g</sub></td><td>interphase turbulent dissipation rate interaction term</td></tr><tr><td>ρ</td><td>density</td></tr><tr><td>Σ<sub>ij</sub></td><td>stress tensor</td></tr><tr><td>τ<sup>t</sup><sub>g</sub></td><td>fluid turbulent time scale</td></tr><tr><td>τ<sup>t</sup><sub>gs</sub></td><td>eddy-particle interaction time</td></tr><tr><td>τ<sup>F</sup><sub>gs</sub></td><td>mean particle relaxation time</td></tr><tr><td>τ<sup>c</sup><sub>s</sub></td><td>inter-particle collision time</td></tr><tr><td>φ<sub>s,ij</sub></td><td>frictional tensor</td></tr></table>

Abbreviation

<table><tr><td>AR</td><td>air reactor</td></tr><tr><td>CFB</td><td>circulating fluidized bed</td></tr><tr><td>CFD</td><td>computational fluid dynamics</td></tr><tr><td>CLC</td><td>chemical looping combustion</td></tr><tr><td>CLOU</td><td>chemical looping with oxygen uncoupling</td></tr><tr><td>DEM</td><td>Discrete Element Method</td></tr><tr><td>EMMS</td><td>energy minimization multi-scale</td></tr><tr><td>FR</td><td>fuel reactor</td></tr><tr><td>HPC</td><td>high performance computing</td></tr><tr><td>OC</td><td>oxygen carrier</td></tr></table>

724 References

725 Abad, A., Adánez, J., Cuadrat, A., García-Labiano, F., Gayán, P., de Diego, L. F., 2011.

726 Kinetics of redox reactions of ilmenite for chemical-looping combustion. Chemical Engineering

727 Science 66 (4), 689–702.
728

Abad, A., Gayán, P., Pérez-Vega, R., García-Labiano, F., de Diego, L., Mendiara, T., Izquierdo, M., Adánez, J., 2020. Evaluation of different strategies to improve the efficiency of coal conversion in a 50 kWth chemical looping combustion unit. Fuel 271, 117514.

729

730

731

Abad, A., Pérez-Vega, R., de Diego, L. F., García-Labiano, F., Gayán, P., Adánez, J., 2015. Design and operation of a 50 kWth Chemical Looping Combustion (CLC) unit for solid fuels. Applied Energy 157, 295 – 303.

732

733

734

Bennani, L., Neau, H., Baudry, C., Laviéville, J., Fede, P., Simonin, O., 2017. Numerical simulation of unsteady dense granular flows with rotating geometries. Chemical Engineering Research and Design 120, 333 – 347.

735

736

737

Berguerand, N., Lyngfelt, A., 2008. Design and operation of a 10 kWth chemical-looping combustor for solid fuels – Testing with South African coal. Fuel 87 (12), 2713 – 2726.

738

739

Boelle, A., Balzer, G., Simonin, O., 1995. Second-order prediction of the particle-phase stress tensor of inelastic spheres in simple shear dense suspensions. American Society of Mechanical Engineers, Fluids Engineering Division (Publication) FED 228, 9–18.

740

741

742

Cao, Y., Casenas, B., Pan, W.-P., 2006. Investigation of chemical looping combustion by solid fuels. 2. Redox reaction kinetics and product characterization with coal, biomass, and solid waste as solid fuels and CuO as an oxygen carrier. Energy & Fuels 20 (5), 1845–1854.

743

744

745

Chen, X., Ma, J., Tian, X., Wan, J., Zhao, H., 2019. CPFD simulation and optimization of a 50 kWth dual circulating fluidized bed reactor for chemical looping combustion of coal. International Journal of Greenhouse Gas Control 90, 102800.

746

747

748

Cloete, S., Johansen, S. T., Amini, S., 2012. Performance evaluation of a complete lagrangian KTGF approach for dilute granular flow modelling. Powder Technology 226, 43–52.

749

750

Cundall, P. A., Strack, O. D. L., 1979. A discrete numerical model for granular assemblies. Géotechnique 29 (1), 47–65.

751

752

Di Renzo, A., Napolitano, E. S., Di Maio, F. P., 2021. Coarse-grain DEM modelling in fluidized bed simulation: A review. Processes 9 (2).

753

754

EDF R&D, 2017. NEPTUNE_CFD Version 4.0.1 User Guide. EDF R&D, Fluid Dynamics, Power Generation and Environment Department, Multi-Phase en Reactive Flow Group, 6 Quai Watier, 78401 Chatou CEDEX, France.

755

756

757

Fede, P., Simonin, O., Ingram, A., 2016. 3D numerical simulation of a lab-scale pressurized dense fluidized bed focussing on the effect of the particle-particle restitution coefficient and particle-wall boundary conditions. Chemical Engineering Science 142, 215–235.

758

759

760

Février, P., Simonin, O., Squires, K. D., 2005. Partitioning of particle velocities in gas-solid turbulent flows into a continuous field and a spatially uncorrelated random distribution: theoretical formalism and numerical study. Journal of Fluid Mechanics 533, 1–46.

761

762

763

Fox, R. O., 2014. On multiphase turbulence models for collisional fluid–particle flows. Journal of Fluid Mechanics 742, 368–424.

764

765

Gidaspow, D., 1994. In: Multiphase flow and fluidization: continuum and kinetic theory descriptions. Academic press, San Diego.

766
728

729

730

Abad, A., Gayán, P., Pérez-Vega, R., García-Labiano, F., de Diego, L., Mendiara, T., Izquierdo, M., Adánez, J., 2020. Evaluation of different strategies to improve the efficiency of coal conversion in a 50 kWth chemical looping combustion unit. Fuel 271, 117514.

731

732

733

Abad, A., Pérez-Vega, R., de Diego, L. F., García-Labiano, F., Gayán, P., Adánez, J., 2015. Design and operation of a 50 kWth Chemical Looping Combustion (CLC) unit for solid fuels. Applied Energy 157, 295 – 303.

734

735

736

Bennani, L., Neau, H., Baudry, C., Laviéville, J., Fede, P., Simonin, O., 2017. Numerical simulation of unsteady dense granular flows with rotating geometries. Chemical Engineering Research and Design 120, 333 – 347.

737

738

Berguerand, N., Lyngfelt, A., 2008. Design and operation of a 10 kWth chemical-looping combustor for solid fuels – Testing with South African coal. Fuel 87 (12), 2713 – 2726.

739

740

741

Boelle, A., Balzer, G., Simonin, O., 1995. Second-order prediction of the particle-phase stress tensor of inelastic spheres in simple shear dense suspensions. American Society of Mechanical Engineers, Fluids Engineering Division (Publication) FED 228, 9–18.

742

743

744

Cao, Y., Casenas, B., Pan, W.-P., 2006. Investigation of chemical looping combustion by solid fuels. 2. Redox reaction kinetics and product characterization with coal, biomass, and solid waste as solid fuels and CuO as an oxygen carrier. Energy & Fuels 20 (5), 1845–1854.

745

746

747

Chen, X., Ma, J., Tian, X., Wan, J., Zhao, H., 2019. CPFD simulation and optimization of a 50 kWth dual circulating fluidized bed reactor for chemical looping combustion of coal. International Journal of Greenhouse Gas Control 90, 102800.

748

749

Cloete, S., Johansen, S. T., Amini, S., 2012. Performance evaluation of a complete lagrangian KTGF approach for dilute granular flow modelling. Powder Technology 226, 43–52.

750

751

Cundall, P. A., Strack, O. D. L., 1979. A discrete numerical model for granular assemblies. Géotechnique 29 (1), 47–65.

752

753

Di Renzo, A., Napolitano, E. S., Di Maio, F. P., 2021. Coarse-grain DEM modelling in fluidized bed simulation: A review. Processes 9 (2).

754

755

756

EDF R&D, 2017. NEPTUNE_CFD Version 4.0.1 User Guide. EDF R&D, Fluid Dynamics, Power Generation and Environment Department, Multi-Phase en Reactive Flow Group, 6 Quai Watier, 78401 Chatou CEDEX, France.

757

758

759

Fede, P., Simonin, O., Ingram, A., 2016. 3D numerical simulation of a lab-scale pressurized dense fluidized bed focussing on the effect of the particle-particle restitution coefficient and particle-wall boundary conditions. Chemical Engineering Science 142, 215–235.

760

761

762

Février, P., Simonin, O., Squires, K. D., 2005. Partitioning of particle velocities in gas-solid turbulent flows into a continuous field and a spatially uncorrelated random distribution: theoretical formalism and numerical study. Journal of Fluid Mechanics 533, 1–46.

763

764

Fox, R. O., 2014. On multiphase turbulence models for collisional fluid-particle flows. Journal of Fluid Mechanics 742, 368–424.

765

766

Gidaspow, D., 1994. In: Multiphase flow and fluidization: continuum and kinetic theory descriptions. Academic press, San Diego.
767

Gobin, A., Neau, H., Simonin, O., Llinas, J., Reiling, V., Sélo, J., 2003. Fluid dynamic numerical simulation of a gas phase polymerization reactor. *International Journal for Numerical Methods in Fluids* 43 (10-11), 1199-1220.

768

769

770

771

Goniva, C., Kloss, C., Deen, N. G., Kuipers, J. A., Pirker, S., 2012. Influence of rolling friction on single spout fluidized bed simulation. *Particuology* 10 (5), 582-591.

772

773

774

775

Hamidouche, Z., Masi, E., Fede, P., Ansart, R., Neau, H., Hemati, M., Simonin, O., 2018. Chapter two - numerical simulation of multiphase reactive flows. In: *Bridging Scales in Modelling and Simulation of Non-Reacting and Reacting Flows. Part I. Vol. 52 of Advances in Chemical Engineering*. Academic Press, pp. 51-124.

776

777

778

Hamidouche, Z., Masi, E., Fede, P., Simonin, O., Mayer, K., Penthor, S., 2019. Unsteady three-dimensional theoretical model and numerical simulation of a 120-kW chemical looping combustion pilot plant. *Chemical Engineering Science* 193, 102-119.

779

780

781

782

Haus, J., Lindmüller, L., Dymala, T., Jarolin, K., Feng, Y., Hartge, E.-U., Heinrich, S., Werther, J., 2020. Increasing the efficiency of chemical looping combustion of biomass by a dual-stage fuel reactor design to reduce carbon capture costs. *Mitigation and Adaptation Strategies for Global Change* 25 (6), 969-986.

783

784

Igci, Y., Sundaresan, S., 2011a. Constitutive models for filtered two-fluid models of fluidized gas-particle flows. *Industrial & Engineering Chemistry Research* 50 (23), 13190-13201.

785

786

Igci, Y., Sundaresan, S., 2011b. Verification of filtered two-fluid models for gas-particle flows in risers. *AIChE Journal* 57 (10), 2691-2707.

787

788

789

Jenkins, J. T., Richman, M. W., 1986. Grad's 13-moment system for a dense gas of inelastic spheres. In: *The Breadth and Depth of Continuum Mechanics*. Springer Berlin Heidelberg, Berlin, Heidelberg, pp. 647-669.

790

791

Johnson, P. C., Jackson, R., 1987. Frictional-collisional constitutive relations for granular materials, with application to plane shearing. *Journal of Fluid Mechanics* 176, 67-93.

792

793

Johnson, P. C., Nott, P., Jackson, R., 1990. Frictional-collisional equations of motion for particulate flows and their application to chutes. *Journal of Fluid Mechanics* 210, 501-535.

794

795

796

Kim, H. R., Wang, D., Zeng, L., Bayham, S., Tong, A., Chung, E., Kathe, M. V., Luo, S., McGiveron, O., Wang, A., Sun, Z., Chen, D., Fan, L., 2013. Coal direct chemical looping combustion process: Design and operation of a 25-kWth sub-pilot unit. *Fuel* 108, 370-384.

797

798

Kolehmainen, J., Ozel, A., Sundaresan, S., 2018. Eulerian modelling of gas-solid flows with triboelectric charging. *Journal of Fluid Mechanics* 848, 340-369.

799

800

801

Langørgen, Ø., Saanum, I., Haugen, N. E. L., 2017. Chemical looping combustion of methane using a copper-based oxygen carrier in a 150 kW reactor system. *Energy Procedia* 114, 352-360.

802

803

804

Laviéville, J., Deutsch, E., Simonin, O., 1995. Large eddy simulation of interactions between colliding particles and a homogeneous isotropic turbulence field. *American Society of Mechanical Engineers, Fluids Engineering Division (Publication) FED 228*, 347-357.
767

Gobin, A., Neau, H., Simonin, O., Llinas, J., Reiling, V., Sélo, J., 2003. Fluid dynamic numerical simulation of a gas phase polymerization reactor. *International Journal for Numerical Methods in Fluids* 43 (10-11), 1199-1220.

768

769

770

Goniva, C., Kloss, C., Deen, N. G., Kuipers, J. A., Pirker, S., 2012. Influence of rolling friction on single spout fluidized bed simulation. *Particuology* 10 (5), 582-591.

771

772

Hamidouche, Z., Masi, E., Fede, P., Ansart, R., Neau, H., Hemati, M., Simonin, O., 2018. Chapter two - numerical simulation of multiphase reactive flows. In: *Bridging Scales in Modelling and Simulation of Non-Reacting and Reacting Flows. Part I. Vol. 52 of Advances in Chemical Engineering*. Academic Press, pp. 51-124.

773

774

775

776

Hamidouche, Z., Masi, E., Fede, P., Simonin, O., Mayer, K., Penthor, S., 2019. Unsteady three-dimensional theoretical model and numerical simulation of a 120-kW chemical looping combustion pilot plant. *Chemical Engineering Science* 193, 102-119.

777

778

779

Haus, J., Lindmüller, L., Dymala, T., Jarolin, K., Feng, Y., Hartge, E.-U., Heinrich, S., Werther, J., 2020. Increasing the efficiency of chemical looping combustion of biomass by a dual-stage fuel reactor design to reduce carbon capture costs. *Mitigation and Adaptation Strategies for Global Change* 25 (6), 969-986.

780

781

782

783

Igci, Y., Sundaresan, S., 2011a. Constitutive models for filtered two-fluid models of fluidized gas-particle flows. *Industrial & Engineering Chemistry Research* 50 (23), 13190-13201.

784

785

Igci, Y., Sundaresan, S., 2011b. Verification of filtered two-fluid models for gas-particle flows in risers. *AIChE Journal* 57 (10), 2691-2707.

786

787

Jenkins, J. T., Richman, M. W., 1986. Grad's 13-moment system for a dense gas of inelastic spheres. In: *The Breadth and Depth of Continuum Mechanics*. Springer Berlin Heidelberg, Berlin, Heidelberg, pp. 647-669.

788

789

790

Johnson, P. C., Jackson, R., 1987. Frictional-collisional constitutive relations for granular materials, with application to plane shearing. *Journal of Fluid Mechanics* 176, 67-93.

791

792

Johnson, P. C., Nott, P., Jackson, R., 1990. Frictional-collisional equations of motion for particulate flows and their application to chutes. *Journal of Fluid Mechanics* 210, 501-535.

793

794

Kim, H. R., Wang, D., Zeng, L., Bayham, S., Tong, A., Chung, E., Kathe, M. V., Luo, S., McGiveron, O., Wang, A., Sun, Z., Chen, D., Fan, L., 2013. Coal direct chemical looping combustion process: Design and operation of a 25-kWth sub-pilot unit. *Fuel* 108, 370-384.

795

796

797

Kolehmainen, J., Ozel, A., Sundaresan, S., 2018. Eulerian modelling of gas-solid flows with triboelectric charging. *Journal of Fluid Mechanics* 848, 340-369.

798

799

Langørgen, Ø., Saanum, I., Haugen, N. E. L., 2017. Chemical looping combustion of methane using a copper-based oxygen carrier in a 150 kW reactor system. *Energy Procedia* 114, 352-360.

800

801

802

Laviéville, J., Deutsch, E., Simonin, O., 1995. Large eddy simulation of interactions between colliding particles and a homogeneous isotropic turbulence field. *American Society of Mechanical Engineers, Fluids Engineering Division (Publication) FED 228*, 347-357.

803

804
Leion, H., Mattisson, T., Lyngfelt, A., 2008. Solid fuels in chemical-looping combustion. International Journal of Greenhouse Gas Control 2 (2), 180–193.

Li, J., Kwauk, M., 1994. Particle-fluid two-phase flow : The Energy-Minimization Multi-Scale Method. Metallurgical Industry Press, Beijing, China.

Li, J., Zhang, H., Gao, Z., Fu, J., Ao, W., Dai, J., 2017. CO2 capture with chemical looping combustion of gaseous fuels: An overview. Energy & Fuels 31 (4), 3475–3524.

Linderholm, C., Schmitz, M., Knutsson, P., Lyngfelt, A., 2016. Chemical-looping combustion in a 100-kW unit using a mixture of ilmenite and manganese ore as oxygen carrier. Fuel 166, 533–542.

Lyngfelt, A., 2014. Chemical-looping combustion of solid fuels – status of development. Applied Energy 113, 1869 – 1873.

Lyngfelt, A., Leckner, B., Mattisson, T., 2001. A fluidized-bed combustion process with inherent CO2 separation; application of chemical-looping combustion. Chemical Engineering Science 56 (10), 3101–3113.

Lyngfelt, A., Linderholm, C., 2017. Chemical-looping combustion of solid fuels – status and recent progress. Energy Procedia 114, 371 – 386, 13th International Conference on Greenhouse Gas Control Technologies, GHGT-13, 14-18 November 2016, Lausanne, Switzerland.

Ma, J., Tian, X., Wang, C., Chen, X., Zhao, H., 2018. Performance of a 50 kWth coal-fuelled chemical looping combustor. International Journal of Greenhouse Gas Control 75, 98–106.

Mahalatkar, K., Kuhlman, J., Huckaby, E. D., O’Brien, T., 2011. CFD simulation of a chemical-looping fuel reactor utilizing solid fuel. Chemical Engineering Science 66 (16), 3617–3627.

Markström, P., Linderholm, C., Lyngfelt, A., 2013. Chemical-looping combustion of solid fuels – design and operation of a 100 kW unit with bituminous coal. International Journal of Greenhouse Gas Control 15, 150 – 162.

Mattisson, T., Keller, M., Linderholm, C., Moldenhauer, P., Rydén, M., Leion, H., Lyngfelt, A., 2018. Chemical-looping technologies using circulating fluidized bed systems: Status of development. Fuel Processing Technology 172, 1–12.

May, J., Alobaid, F., Ohlemüller, P., Stroh, A., Ströhle, J., Epple, B., 2018. Reactive two-fluid model for chemical-looping combustion – simulation of fuel and air reactors. International Journal of Greenhouse Gas Control 76, 175 – 192.

Molodtsof, Y., 2003. Theoretical analysis of the flow regimes and their characteristics in vertically flowing gas–solids suspensions. Chemical Engineering Journal 96 (1), 133–143.

Montilla, C., Ansart, R., Simonin, O., 2020. Modelling of the mean electric charge transport equation in a mono-dispersed gas–particle flow. Journal of Fluid Mechanics 902, A12.

Neau, H., Pigou, M., Fede, P., Ansart, R., Baudry, C., Mérigoux, N., Laviéville, J., Fournier, Y., Renon, N., Simonin, O., 2020. Massively parallel numerical simulation using up to 36,000 cpu cores of an industrial-scale polydispersed reactive pressurized fluidized bed with a mesh of one billion cells. Powder Technology 366, 906 – 924.
Leion, H., Mattisson, T., Lyngfelt, A., 2008. Solid fuels in chemical-looping combustion. International Journal of Greenhouse Gas Control 2 (2), 180–193.

Li, J., Kwauk, M., 1994. Particle-fluid two-phase flow : The Energy-Minimization Multi-Scale Method. Metallurgical Industry Press, Beijing, China.

Li, J., Zhang, H., Gao, Z., Fu, J., Ao, W., Dai, J., 2017. CO2 capture with chemical looping combustion of gaseous fuels: An overview. Energy & Fuels 31 (4), 3475–3524.

Linderholm, C., Schmitz, M., Knutsson, P., Lyngfelt, A., 2016. Chemical-looping combustion in a 100-kW unit using a mixture of ilmenite and manganese ore as oxygen carrier. Fuel 166, 533–542.

Lyngfelt, A., 2014. Chemical-looping combustion of solid fuels – status of development. Applied Energy 113, 1869 – 1873.

Lyngfelt, A., Leckner, B., Mattisson, T., 2001. A fluidized-bed combustion process with inherent CO2 separation; application of chemical-looping combustion. Chemical Engineering Science 56 (10), 3101–3113.

Lyngfelt, A., Linderholm, C., 2017. Chemical-looping combustion of solid fuels – status and recent progress. Energy Procedia 114, 371 – 386, 13th International Conference on Greenhouse Gas Control Technologies, GHGT-13, 14-18 November 2016, Lausanne, Switzerland.

Ma, J., Tian, X., Wang, C., Chen, X., Zhao, H., 2018. Performance of a 50 kWth coal-fuelled chemical looping combustor. International Journal of Greenhouse Gas Control 75, 98–106.

Mahalatkar, K., Kuhlman, J., Huckaby, E. D., O'Brien, T., 2011. CFD simulation of a chemical-looping fuel reactor utilizing solid fuel. Chemical Engineering Science 66 (16), 3617–3627.

Markström, P., Linderholm, C., Lyngfelt, A., 2013. Chemical-looping combustion of solid fuels – design and operation of a 100 kW unit with bituminous coal. International Journal of Greenhouse Gas Control 15, 150 – 162.

Mattisson, T., Keller, M., Linderholm, C., Moldenhauer, P., Rydén, M., Leion, H., Lyngfelt, A., 2018. Chemical-looping technologies using circulating fluidized bed systems: Status of development. Fuel Processing Technology 172, 1–12.

May, J., Alobaid, F., Ohlemüller, P., Stroh, A., Ströhle, J., Epple, B., 2018. Reactive two-fluid model for chemical-looping combustion – simulation of fuel and air reactors. International Journal of Greenhouse Gas Control 76, 175 – 192.

Molodtsof, Y., 2003. Theoretical analysis of the flow regimes and their characteristics in vertically flowing gas-solids suspensions. Chemical Engineering Journal 96 (1), 133–143.

Montilla, C., Ansart, R., Simonin, O., 2020. Modelling of the mean electric charge transport equation in a mono-dispersed gas-particle flow. Journal of Fluid Mechanics 902, A12.

Neau, H., Pigou, M., Fede, P., Ansart, R., Baudry, C., Mérigoux, N., Laviéville, J., Fournier, Y., Renon, N., Simonin, O., 2020. Massively parallel numerical simulation using up to 36,000 cpu cores of an industrial-scale polydispersed reactive pressurized fluidized bed with a mesh of one billion cells. Powder Technology 366, 906 – 924.
843

844

845

Ozel, A., Fede, P., Simonin, O., 2013. Development of filtered Euler–Euler two-phase model for circulating fluidised bed: High resolution simulation, formulation and a priori analyses. International Journal of Multiphase Flow 55, 43–63.

846

847

Parker, J. M., 2014. CFD model for the simulation of chemical looping combustion. Powder Technology 265, 47–53.

848

849

Parmentier, J.-F., Simonin, O., Delsart, O., 2012. A functional subgrid drift velocity model for filtered drag prediction in dense fluidized bed. AIChE Journal 58 (4), 1084–1098.

850

851

852

Penthor, S., Stollhof, M., Pröll, T., Hofbauer, H., 2016. Detailed fluid dynamic investigations of a novel fuel reactor concept for chemical looping combustion of solid fuels. Powder Technology 287, 61–69.

853

854

855

Pérez-Astray, A., Mendiara, T., de Diego, L., Abad, A., García-Labiano, F., Izquierdo, M., Adánez, J., 2020. Improving the oxygen demand in biomass CLC using manganese ores. Fuel 274, 117803.

856

857

858

859

Pérez-Vega, R., Abad, A., Gayán, P., García-Labiano, F., Izquierdo, M. T., de Diego, L. F., Adánez, J., 2020. Coal combustion via Chemical Looping assisted by Oxygen Uncoupling with a manganese-iron mixed oxide doped with titanium. Fuel Processing Technology 197, 106184.

860

861

862

Pirker, S., Kahrimanovic, D., Kloss, C., Popoff, B., Braun, M., 2010. Simulating coarse particle conveying by a set of Eulerian, Lagrangian and hybrid particle models. Powder Technology 204 (2), 203–213.

863

864

Pröll, T., Kolbitsch, P., Bolhàr-Nordenkampf, J., Hofbauer, H., 2009. A novel dual circulating fluidized bed system for chemical looping processes. AIChE Journal 55 (12), 3255–3266.

865

866

867

Reinking, Z., Shim, H.-S., Whitty, K. J., Lighty, J. S., 2019. Computational simulation of a 100 kW dual circulating fluidized bed reactor processing coal by chemical looping with oxygen uncoupling. International Journal of Greenhouse Gas Control 90, 102795.

868

869

Schneiderbauer, S., Pirker, S., 2014. Filtered and heterogeneity-based subgrid modifications for gas-solid drag and solid stresses in bubbling fluidized beds. AIChE Journal 60 (3), 839–854.

870

871

872

Shao, Y., Agarwal, R. K., Wang, X., Jin, B., 2020. Numerical simulation of a 3D full loop iG-CLC system including a two-stage counter-flow moving bed air reactor. Chemical Engineering Science 217, 115502.

873

874

875

Shao, Y., Agarwal, R. K., Wang, X., Jin, B., 2021. Review of computational fluid dynamics studies on chemical looping combustion. Journal of Energy Resources Technology 143 (8), 080802.

876

877

Shen, L., Wu, J., Xiao, J., Song, Q., Xiao, R., 2009. Chemical-looping combustion of biomass in a 10 kWth reactor with iron oxide as an oxygen carrier. Energy & Fuels 23 (5), 2498–2505.

878

879

880

881

Simonin, O., 2000. Statistical and continuum modelling of turbulent reactive particulate flows. part 1: Theoretical derivation of dispersed Eulerian modelling from probability density function kinetic equation. In: Lecture Series. Vol. 6 of Theoretical and Experimental Modeling of Particulate Flows. von Karman Institute for Fluid Dynamics Rhode Saint Genèse, Belgium.