Numerical Simulation of Chemical Looping and Calcium Looping Combustion Processes for Carbon Capture

Abstract

Efficient carbon capture and storage (CCS) technologies are needed to address the rising carbon emissions from power generation using fossil fuels that have been linked to global warming and climate change. Chemical looping combustion (CLC) is one such technology that has shown great promise due to its potential for high-purity carbon capture at low cost. Another CCS technology that has garnered interest in recent years is calcium looping (CaL), which utilizes calcium oxide and the carbonation-calcination equilibrium reactions to capture CO2 from the flue stream of fossil fuel power plants. Computational fluid dynamics (CFD) simulations of two CLC reactors are presented in this chapter, along with system level simulations of CaL for postcombustion carbon capture. CFD simulation of a CLC reactor based on a dual fluidized bed reactor is developed using the Eulerian approach to characterize the chemical reactions in the system. The solid phase consists of a Fe-based oxygen carrier while the gaseous fuel used is syngas. Later, the detailed hydrodynamics in a CLC system designed for solid coal fuel is presented based on a cold flow experimental setup at National Energy Technology Laboratory using the Lagrangian particle-tracking method. The process simulation of CaL using Aspen Plus shows an increasing marginal energy penalty associated with an increase in the CO2 capture efficiency, which suggests a limit on the maximum carbon capture efficiency in practical applications of CaL before the energy penalty becomes too large.

1. Introduction

The relationship between the global surface temperature of the Earth and the concentration of CO₂ was identified by Arrhenius as early as 1896 [1]. Since then, the concentration of CO₂ in the atmosphere has risen from 280 ppm to around 400 ppm today, largely due to carbon emissions from fossil fueled power plants and other anthropogenic sources. As such, there is an imminent need for high-efficiency carbon capture and storage (CCS) technologies to avoid the “unequivocal warming of the global climate system” [2]. Chemical looping combustion (CLC) is one such CCS technology that is well-suited for high-efficiency, low-cost carbon capture. In the CLC process, fuel combustion takes place in the fuel reactor using oxygen supplied by a metal oxide oxygen carrier. The reduced oxygen carrier is pneumatically transported to the air reactor where it is reoxidized in air; it then circulates back into the fuel reactor to complete the loop. The typical CLC setup employing two dual fluidized bed reactors is shown in Figure 1(a) [3]. Alternatively, a single packed bed reactor can be swapped between fuel reactor and air reactor configurations through a high temperature gas switching system as shown in Figure 1(b) [3].

Figure 1.

Schematic representation of a chemical looping combustion system with (a) dual interconnected fluidized beds and (b) packed bed with alternating flow [3].

Since the fuel combustion in a CLC system takes place in the absence of air, the flue stream from the fuel reactor is not diluted or contaminated by other gases such as nitrogen. Hence, CLC can produce a high-purity stream of CO₂ available for capture at the fuel reactor without the need for the energy expensive gas separation process required by other CCS technologies such as oxy-fuel combustion. The only energy cost of separation associated with CLC is the cost of solid recirculation. Previous works based on energy and exergy analysis have demonstrated that CLC systems can achieve power efficiencies greater than 50% along with nearly complete CO₂ capture [4–8]. The low-cost carbon capture associated with CLC has a direct bearing on the cost of electricity, as confirmed by a techno-economic study published recently by the National Energy Technology Laboratory (NETL) that reported that the cost of electricity for a CLC plant using Fe₂O₃ oxygen carrier is $115.1 per MWh compared with $137.3 per MWh for a conventional pulverized coal power plant with a postcombustion amine-based CO₂ absorption system [9].

Unlike CLC where the combustion process is altered to facilitate carbon capture, CCS technologies can also be retrofitted into existing power plants postcombustion setting. The use of calcium oxide (CaO) as the sorbent to capture CO₂ from the flue stream of a power plant has been proposed in several studies [10–13]. The calcium looping (CaL) process consists of two interconnected reactors called the carbonator and the calciner in which the carbonation and calcination reactions take place respectively. The overall carbonation-calcination equilibrium reaction is given by

${\text{CaO}}_{(\text{s})}+{\text{ CO}}_{\text{2}(\text{g})}\leftrightharpoons {\text{CaCO}}_{\text{3}(\text{s})}$E1

The carbonation reaction entraps the CO₂ from the flue gas stream using the CaO sorbent to form calcium carbonate (CaCO₃). The flue stream exiting the carbonator is CO₂-lean and can be exhausted into the atmosphere. The solid CaCO₃ from the carbonator is transported to the calciner where it is heated to decompose into CaO and CO₂. The calciner produces a stream of pure CO₂ that is sent for pressurized storage for subsequent sequestration or use. The CaO is transported back into the carbonator to complete the loop. The typical setup for calcium looping (CaL) mirrors that of a chemical looping combustion plant, with the carbonator and calciner analogous to the air and fuel reactors in CLC and the CaO/CaCO₃ equivalent to the oxidized/reduced metal oxide oxygen carrier. A schematic representation of the CaL process is shown in Figure 2. The make-up CaO flow and solids disposal are required to maintain reactivity of the sorbent; if the CaO was not replenished, the reaction rates would degenerate over time as the loop ran its course. It is noted that CaL can also be employed to remove carbon from the fuel prior to combustion and generate a stream of H₂ for clean burning. However, the precombustion CaL scenario requires significant alteration to the combustion process and loses the advantage of being able to be retrofitted into existing power plants; it also has a higher energy penalty compared with the postcombustion CaL process considered in this chapter.

Figure 2.

Schematic representation of a calcium looping system with interconnected reactors.

Computational fluid dynamics (CFD) simulations of chemical looping combustion provide an efficient means to analyze the performance and characterize the fluidization behavior and chemical kinetics in the system. Past CFD simulations of CLC have primarily focused on the single reactor setup similar to Figure 1(b), such as the work of Mahalatkar et al. [14, 15] demonstrating the capability of CFD to match the reaction mechanics inside a CLC fuel reactor. However, the single reactor setup cannot be operated with solid fuels since there is no way to remove the coal from the system prior to switching to the air reactor configuration. Given that coal is projected to remain the dominant fossil fuel in the near future [16], it is necessary to establish a credible CFD simulation for CLC based on the dual fluidized bed setup shown in Figure 1(a).

A detailed three-dimensional CFD simulation of the experimental CLC reactor of Abad et al. [17] is presented in this chapter employing the Eulerian or continuum approach to model the solid phase. Later, the cold flow behavior of the fully-looped CLC system at NETL is modeled employing the Lagrangian discrete element method (DEM) as a follow up to the CFD simulation conducted for this case previously by Parker [18]. The CFD simulations discussed in this work are among the very few present in the current literature of the complete circulating dual fluidized bed setup for CLC. Finally, a system level simulation of postcombustion CaL in Aspen Plus is presented to evaluate the energy penalty associated with the carbon capture. The estimation of the energy penalty in CaL is of great interest in the field of CCS; the goal of an ideal carbon capture process is to consume the least amount of energy while achieving a high CO₂ capture efficiency.

2. Eulerian simulation of the experimental CLC reactor of Abad et al.

In this section, the laboratory scale experiment of Abad et al. [17] is used as a basis to perform a detailed CFD simulation of a CLC system using the Eulerian multifluid approach. It is one of the few CFD models in the literature of a complete circulating dual fluidized bed setup with chemical reactions. The fluidization behavior in both air and fuel reactor beds and the circulation of the oxygen carrier between the beds is investigated and compared with the experiment. The simulation results of the chemical reactions in the fuel reactor are evaluated against the outlet concentrations of the flue gases in the experiment.

The experiment uses the two-compartment fluidized bed design proposed by Chong et al. [19] and further investigated by Fang et al. [20]. The experimental reactor setup is illustrated in Figure 3. Dimensions and additional details can be found in the work of Abad et al. [17]. The experiment used a Fe-based oxygen carrier consisting of 60% Fe₂O₃ by mass and 40% Al₂O₃ designated as F6A1100. The gaseous fuels used in the experiment are natural gas, consisting of primarily CH₄ and syngas consisting of a mixture of 50% CO and 50% H₂. The fluidization velocity in the air reactor is greater than the terminal velocity of the oxygen carrier particles and carries the particles upward. The flow then undergoes a sudden expansion (i.e., a reduction in velocity) in the particle separator at the top of the reactor, which causes the particles to fall into the down-comer and enter the fuel reactor. The fuel reactor constitutes a bubbling fluidized bed given the fluidization velocity is smaller than the terminal velocity of the particles. The flue streams from both the reactors are connected to a gas analyzer to measure the concentrations of the outlet gases.

2.1. Numerical solution procedure

A numerical model of the experimental CLC reactor of Abad et al. [17] is developed using the commercial CFD software ANSYS Fluent, release version 14.5 [21, 22]. The complete Navier-Stokes equations of fluid dynamics are solved to account for the chemical active multiphase flow with heat transfer. The Eulerian two-fluid model is used to approximate the solid phase as a secondary fluid phase by averaging particle variables such as mass, velocity, temperature, etc. over a region that is large compared with the particle size. Interactions between the solid and gas phases are included in the model via constitutive equations for solid phase pressure and viscosity that are provided by the kinetic theory of granular flow, an extension of the classical kinetic gas theory that includes inelastic interparticle interactions [23, 24].

2.1.1. Eulerian two-fluid model equations

For multiphase simulations, the Navier-Stokes equations are modified by including the phase volume fraction α to account for the presence of other phases [22]. The continuity equation for phase q is given as

$\frac{\partial }{\partial t}\left({\alpha }_q{\rho }_q\right)+\nabla \cdot \left({\alpha }_q{\rho }_q{u}_q\right)={{\displaystyle \sum }}^{\text{}}\left({\dot{m}}_{pq}-{\dot{m}}_{qp}\right)$E2

where ${\dot{m}}_{pq}$ is the mass transfer rate from the pth phase to the qth phase. The current simulation contains one gas phase comprising six species (CO, CO₂, H₂, H₂O, O₂, and N₂) and one solid phase comprising two species (Fe₂O₃ and Fe₃O₄). The species transport equation is given as

$\frac{\partial }{\partial t}\left({\alpha }_q{\rho }_q{Y}_{iq}\right)+\nabla \cdot \left({\alpha }_q{\rho }_q{u}_q{Y}_{iq}\right)={{\displaystyle \sum }}^{\text{}}({\dot{m}}_{ij}^{qp}-{\dot{m}}_{ji}^{pq})$E3

where Y_iq is the mass fraction of the species i in the qth phase and ${\dot{m}}_{ij}^{qp}$ is the mass transfer rate from the jth species of the pth phase to the ith species in the qth phase.

The momentum equation for the gas phase is given as

$\frac{\partial }{\partial t}\left({\alpha }_g{\rho }_g{u}_g\right)+\nabla \cdot \left({\alpha }_g{\rho }_g{u}_g{u}_g\right)=-{\alpha }_g\nabla p+\nabla \cdot {\overline{\overline{\tau }}}_g+{\alpha }_g{\rho }_{g}g+{{\displaystyle \sum }}^{\text{}}\left(R_{sg}+{\dot{m}}_{sg}{u}_{sg}-{\dot{m}}_{gs}{u}_{gs}\right).$E4

where the terms in the summation are source terms added to the standard form of the Navier-Stokes momentum equations to account for the momentum transfer between the solid phase and the gas phase. Specifically, R_sg=β_sg(u_s−u_g) is the momentum transfer due to interphase drag and the other terms are due to the transfer of mass. The momentum equation for the solid phase follows from the momentum equation for the gas phase with the source term for interphase drag being equal but opposite.

$\frac{\partial }{\partial t}\left({\alpha }_s{\rho }_s{u}_s\right)+\nabla \cdot \left({\alpha }_s{\rho }_s{u}_s{u}_s\right)=-{\alpha }_s\nabla p+\nabla \cdot {\overline{\overline{\tau }}}_s+{\alpha }_s{\rho }_{s}g+{{\displaystyle \sum }}^{\text{}}\left(R_{rs}+{\dot{m}}_{rs}{u}_{rs}-{\dot{m}}_{sr}{u}_{sr}\right)$E5

The gas inside a CLC reactor can be considered as an incompressible fluid. As such, the fluid stress tensor ${\overline{\overline{\tau }}}_g$ is simply

For the flow conditions in a fuel reactor, the gas can be considered as an incompressible fluid. The fluid stress tensor is simply the Cauchy stress tensor with zero bulk viscosity.

${\overline{\overline{\text{τ}}}}_g={\alpha }_g{\mu }_g\left(\nabla u_g+\nabla u_g^{\text{T}}\right).$E6

On the other hand, the granular solid stress tensor considers all terms in the Cauchy stress tensor,

${\overline{\overline{\tau }}}_s=-p_s\overline{\overline{I}}+{\alpha }_s{\mu }_s\left(\nabla u_s+\nabla u_s^{\text{T}}\right)+{\alpha }_s{\lambda }_s\left(\nabla \cdot u_s\right)\overline{\overline{I}}$E7

where p_s is the solids pressure, μ_s is the granular viscosity, and λ_s is the granular bulk viscosity. The definition of these terms and the interphase exchange coefficient β_s provide the basis for the Eulerian approach for multiphase flow simulation. The solids pressure and granular bulk viscosity used in the present work are according to Lun et al. [25]; the granular viscosity is according to Gidaspow [26].

The experimental reactor of Abad et al. [17] includes a wide range of solid loadings, from the densely-packed loop-seal and down-comer, the bubbling fluidized bed in the fuel reactor, and the dilute regions in the particle separator. The Gidaspow drag model [26] is selected to model the solid-gas interaction because it accounts for the differences in the behavior in dilute and dense regions by switching between the drag prediction of the Ergun equation [27] and the drag model of Wen and Yu [28] based on the solids volume fraction α_s. For α_s > 0.8, the Gidaspow model for the exchange coefficient β_s gives

${\beta }_{sg}=\frac{3}{4}{C}_D\frac{{\alpha }_s{\alpha }_g{\rho }_g|u_s-u_g|}{d_s}{\alpha }_g^{-2.65};C_D=\frac{24}{{\alpha }_g{\mathrm{Re}}_s}\left[1+0.15{\left({\alpha }_g{\mathrm{Re}}_s\right)}^{0.687}\right]$E8

Conversely, for α_s ≤ 0.8,

${\beta }_{sg}=150\frac{{\alpha }_s\left(1-{\alpha }_g\right){\mu }_g}{{\alpha }_g{d}_s^2}+1.75\frac{{\rho }_g{\alpha }_s|u_s-u_g|}{d_s}$E9

where d_s is the particle diameter and Re_s is the Reynolds number based on d_s.

Finally, the energy equation for phase q is expressed in terms of the enthalpy as

$\frac{\partial }{\partial t}\left({\alpha }_q{\rho }_q{h}_q\right)+\nabla \cdot \left({\alpha }_q{u}_q{h}_q\right)={\alpha }_q\frac{\partial p}{\partial t}+\nabla \cdot \left({\overline{\overline{\tau }}}_q\cdot u_q\right)-\nabla \cdot q_q+S_q+{{\displaystyle \sum }}^{\text{}}{Q}_{pq}$E10

where h_q and q_q are the specific enthalpy and heat flux of phase q respectively. As with the continuity and momentum equations, the source terms S_q and Q_pq are included in Eq. (10) to account for the transfer of enthalpy between the solid and gas phases. S_q is the enthalpy source due to chemical reaction and Q_pq is the heat transfer as a function of the temperature difference between the pth phase and the qth phase, given in terms of the interphase heat transfer coefficient h_pq as

$Q_{pq}=h_{pq}\left(T_p-T_q\right)$E11

The interphase heat transfer coefficient is based on the Nusselt number, Nu, the ratio of convective to conductive heat transfer across a boundary.

$h_{pq}=\frac{6{\kappa }_q{\alpha }_p{\alpha }_q{\text{Nu}}_p}{d_p^2}$E12

where κ_q is the thermal conductivity of phase q. In the current simulation, the Nusselt number for the solid phase is obtained using the empirical relation of Gunn [29] given in terms of the Prandtl number Pr as

${\text{Nu}}_s=\left(7-10{\alpha }_q+5{\alpha }_q^2\right)\left(1+0.7{\mathrm{Re}}_s^{0.2}{\mathrm{Pr}}^{1/3}\right)+(1.33-2.4{\alpha }_q+1.2{\alpha }_q^2){\mathrm{Re}}_s^{0.7}{\mathrm{Pr}}^{1/3}$E13

2.1.2. Chemical reaction scheme and rates

The oxygen carrier used in the numerical simulation is F6A1100, consisting of 60% Fe₂O₃ and 40% inert Al₂O₃ by mass, following the experiment of Abad et al. [17]. Of the two gaseous fuels used in the experiment, only syngas has been considered in the current work because the chemical kinetics for the reaction of Fe₂O₃ with the nonmethane components of natural gas is not available. The metal oxide reduction reactions that are used in the simulation are

${\text{3Fe}}_{\text{2}}{\text{O}}_{\text{3}}+\text{ CO }\to {\text{2Fe}}_{\text{3}}{\text{O}}_{\text{4}}+{\text{ CO}}_{\text{2}}$E14

${\text{3Fe}}_{\text{2}}{\text{O}}_{\text{3}}+{\text{ H}}_{\text{2}}\to {\text{2Fe}}_{\text{3}}{\text{O}}_{\text{4}}+{\text{ H}}_{\text{2}}\text{O}$E15

Exact reaction rates for the reduction of F6A1100 with CO and H₂ are not available in the literature; the reaction rates for the reduction of hematite (Fe₂O₃) with CO and H₂ obtained from the work of Mattisson et al. [30] and further developed by Mahalatkar et al. [14] are used as an approximation. Based on these papers, the reaction rates $\dot{m}$ (in kg/s per cell volume or kg/(m³ s)) of the fuel gases with Fe₂O₃ are given by

${\dot{m}}_{{\text{H}}_2}=\frac{k_{{\text{H}}_2}{R}_o}{2M{W}_{{\text{O}}_2}}{\rho }_{avg}{\alpha }_s\left(Y_{{\text{Fe}}_2{\text{O}}_3}+Y_{{\text{Fe}}_3{\text{O}}_4}\frac{{\nu }_{{\text{Fe}}_2{\text{O}}_3}M{W}_{{\text{Fe}}_2{\text{O}}_3}}{{\nu }_{{\text{Fe}}_3{\text{O}}_4}M{W}_{{\text{Fe}}_3{\text{O}}_4}}\right){\left(1-X\right)}^{2/3}M{W}_{{\text{H}}_2}$E16

And

${\dot{m}}_{\text{CO}}=\frac{k_{\text{CO}}{R}_o}{2M{W}_{{\text{O}}_2}}{\rho }_{avg}{\alpha }_s\left(Y_{{\text{Fe}}_2{\text{O}}_3}+Y_{{\text{Fe}}_3{\text{O}}_4}\frac{{\nu }_{{\text{Fe}}_2{\text{O}}_3}M{W}_{{\text{Fe}}_2{\text{O}}_3}}{{\nu }_{{\text{Fe}}_3{\text{O}}_4}M{W}_{{\text{Fe}}_3{\text{O}}_4}}\right){\left(1-X\right)}^{2/3}M{W}_{\text{CO}}$E17

where k is the nominal reaction rate based on the Arrhenius rate, R_o is the oxygen carrying capacity, MW is the molecular weight (in kg/kmol), Y is the mass fraction, ν is the stoichiometric coefficient, and X is the conversion fraction based on the fully reduced state; in each case, the subscript identifies the species under consideration. More details of the reaction rate derivation can be found in the work of Mahalatkar et al. [14]. The reaction rates identified in Eqs. (16) (17) are implemented into the numerical simulation through separate user-defined functions.

2.2. Three-dimensional simulation of Abad et al. experiment

The results of a two-dimensional simulation of the experiment of Abad et al. [17] can be found in the authors’ previous work [31]. Although the two-dimensional model successfully captured the salient features of the fluidization behavior in the dual fluidized bed system, it was unable to produce the expected concentrations of CO₂ and H₂O in the fuel reactor because of the inadequacy of the two-dimensional simulation in modeling the gaseous diffusion, which is an inherently three-dimensional process. Therefore, a three-dimensional simulation of Abad et al.’s experiment [17] is performed to produce a more accurate match for the chemical reactions between the simulation and experiment. The computational domain is an exact representation of the experimental reactor shown in Figure 3. A structured mesh with around 45,000 elements is used to model the geometry; the grid is relatively fine in the lower part of the reactor where the solids loading is densely-packed (minimum cell volume of around 10⁻⁸ m) and coarser in the dilute upper part (cell volume of around 10⁻⁷ m). The mesh is shown in Figure 4.

Figure 3.

Sketch of experimental reactor [17]: (1) air reactor, (2) down-comer, (3) fuel reactor, (4) slot, (5) gas distributor plate, (6) wind box, (7) reactor part, (8) particle separator, and (9) leaning wall. The symbols (×) and (o) indicate fluidization in the down-comer and slot.

Figure 4.

Computational domain and grid for three-dimensional CFD simulation with detailed view of lower part.

The oxygen carrier used by Abad et al. [17] has a density of 2150 kg/m³ with a diameter of 90–212 μm; the average value of 150 μm is used in the simulation. The batch processing results of Abad et al. [17] are used as a basis for the simulation. The initial solids loading in the bed is about 300 g, of which 110 g is in the fuel reactor, in line with the experiment. The temperature in the system is set at 850°C or 1123 K based on the experimental conditions. It is expected that the densely-packed solids in the down-comer and slot regions will be enough to keep the leakage to a minimum without the need for an active pressure differential between the reactors. As such, both reactors are set at atmospheric pressure. The initial oxygen carrier mass in the fuel reactor is sufficient to react with all the injected fuel, so the fuel conversion is not affected by reoxidation in the air reactor. The CFD simulation is thus considerably simplified by setting the fluidization gas in the air reactor to an inert gas (in this case, nitrogen). The secondary phase mass fraction is set to zero at both fuel and air reactor inlets because no new oxygen carrier is added. The numerical parameters used in the CFD simulation are summarized in Table 1. The simulation was run on a Dell workstation with a quad-core Intel Xeon CPU for around four days to complete 30 s of simulation. The contours of the mass fraction of CO₂ are shown in Figure 5.

Primary phase	Fuel-gas mixture
Secondary phase	Oxygen carrier (F6AL1100)
Average particle diameter	150 μm
Average particle density	2150 kg/m^3
Initial bed mass	∼180 g
Fluidizing gas composition in fuel reactor	50% CO, 50% H2
Fluidizing gas composition in air reactor	100% N2
Inlet boundary condition in fuel reactor	Velocity inlet with velocity 0.1 m/s
Inlet boundary condition in air reactor	Velocity inlet with velocity 0.5 m/s
Outlet boundary condition in fuel reactor	Pressure outlet at atmospheric pressure
Outlet boundary condition in air reactor	Pressure outlet at atmospheric pressure
Operating temperature	1123 K
Solids pressure	Lun et al. [25]
Granular bulk viscosity	Lun et al. [25]
Granular viscosity	Gidaspow [26]
Drag law	Gidaspow [26]
Heat transfer coefficient	Gunn [29]
Numerical scheme	Phase-coupled SIMPLE
Time step size	0.0005 s
Iterations per time step	20

Table 1.

Modeling parameters for three-dimensional CFD simulation of the Abad et al. experiment.

Figure 5.

Contours of CO₂ mass fraction for the first 10 seconds of three-dimensional simulation showing the increased diffusion and absence of the vortex pattern compared with the two-dimensional case.

As expected, the three-dimensional simulation exhibits greater diffusion compared with the two-dimensional case presented in reference [31]. The local mass fraction of CO₂ at the base of the bed where the injected CO first comes into contact with the oxygen carrier and begins to react is around 15%. From there, the CO₂ diffuses through the fuel reactor more homogeneously as it travels toward the outlet; the vortex patterns observed in the two-dimensional case [31] are notably absent. The absence of the vortex pattern can be quantitatively confirmed by the plots of the mass fractions of CO₂ and H₂O at the fuel reactor outlet as shown in Figure 6, which also includes the two-dimensional results from reference [31]. The outlet mass fractions of both gases are initially lower because the gases have to diffuse through the existing N₂ in the fuel reactor instead of displacing it. The large fluctuations in the outlet mass fraction caused by pockets of reversed flow are also eliminated. The outlet mass fractions of both gases keep increasing as the fuel reaction produces more and more; by 20 s, both CO₂ and H₂O have exceeded their stagnation values from the two-dimensional simulation. The outlet mass fraction of H₂O reaches the expected value from Abad et al.’s experiment around 30 s [17].

Figure 6.

Mass fractions of CO₂ and H₂O at the fuel reactor outlet for the three-dimensional simulation of the CLC reactor of Abad et al. [17].

The three-dimensional simulation shows a significant increase in the mass fraction of the flue gases at the fuel reactor outlet. However, despite the improvement, the mass fraction of CO₂ still shows a significant discrepancy from the experimental value, which may be due to various external factors. It is known that significant apparent diffusion can occur in gases when they travel through pipes [32]. In the experiment [17], the gas streams from the reactors were pipelined to an electric cooler and then to the gas analyzer. Hence, it is reasonable to expect that the concentrations measured at the gas analyzer may be different from the concentrations present right at the fuel reactor outlet. It should also be noted that the reaction rate kinetics used in the simulation were based on the experimental study of Mattisson et al. [30] using hematite (Fe₂O₃), while the oxygen carrier used by Abad et al. [17] was F6A1100 comprising 60% Fe₂O₃ and 40% Al₂O₃. One of the reasons F6A1100 is preferred over hematite as the oxygen carrier for CLC operation is its improved reactivity owing to an increased apparent surface area due to the presence of the porous Al₂O₃[33]. Thus, it stands to reason that the experiment would show a higher concentration of the reaction products compared with the current simulation. Similar discrepancies between experiment and simulation have been previously noted in the work of Mahalatkar et al. [15] where the reaction rates obtained from experiment for the char gasification reaction had to be doubled in the simulation in order to match the result. Further research is required to determine more accurate empirical formulas for the reduction of F6A1100 specifically to improve the accuracy of the results of the CFD simulation.

3. Lagrangian simulation of complete CLC reactor at NETL

The Eulerian two-fluid model can accurately capture the bulk behavior of the solid phase in the dual fluidized bed reactor for CLC using gaseous fuels. However, given coal is projected to remain one of the dominant fossil fuels in the near future, the concept of coal-direct chemical looping combustion (CD-CLC) with in situ gasification has garnered significant interest in recent years. In the CD-CLC process, the oxygen carrier needs to be formed into particles with a relatively larger diameter compared with pulverized coal for easier separation. The work of Gryczka et al. [34] with larger particles has suggested that accurate numerical representation of particle dynamics is not likely to be achieved using the granular solid phase approximation due to “the inadequacies of the continuum model.” The inaccuracy arises from the nonphysical closure terms used in the Eulerian model such as the frictional solids viscosity or the solids pressure based on the kinetic theory of granular flow. Thus, for a more detailed understanding of the hydrodynamics in a multiphase system inside a CLC reactor, the Lagrangian particle-based approach is employed to study the CD-CLC reactor system at NETL, previously investigated by Parker [18].

The geometry of the CD-CLC system used at NETL comprises an air reactor, cyclone, loop-seal, and fuel reactor, as shown in Figure 7(a). In the cold flow experiment, the oxygen carrier particles start from the bottom of the air reactor and move up along the riser and into the cyclone. In the cyclone, the particles are separated from the air stream and drop into the loop-seal due to the gravity. After passing through the slightly fluidized loop-seal, the particles move into the fuel reactor. The oxygen carrier particles exit the fuel reactor and return to the air reactor through an L-valve. The reactor dimensions and additional descriptions of the various components of the CD-CLC setup at NETL can be found in the work of Parker [18]. The computational grid, shown in Figure 7(b), is an exact representation of the geometry. The oxygen carrier particles in the cold flow experiment are primarily ilmenite (FeTiO₃) with some uncombined TiO₂ and Fe₂O₃ as well. The ilmenite particles used in the experiment had a size distribution of 13-320 μm. The particle size used in the numerical simulation corresponds to the median particle size of 150 μm with the average density of ilmenite of 4450 kg/m³.

Figure 7.

(a) Geometry of the coal-direct chemical-looping combustion system at NETL [18], and (b) computational grid with static pressure measurement locations marked S1–S6.

3.1. Numerical solution procedure

The equations for mass and momentum conservation for the fluid phase are identical to those used in the Eulerian model given in Eqs. (2) and (4) with the exception that the source term in Eq. (4) for the solid-gas momentum exchange term, R_sg, is obtained from the average of the drag forces acting on all the discrete particles in a given computational cell. The shear stress term in the momentum equation is given in Eq. (6). Since the current simulation considers a cold flow with no species transport, the energy and species conservation equations are not applicable in this case.

3.1.1. Particle equations

In the Lagrangian approach, the motion of each solid particle is calculated by summing the forces acting on the particle and applying Newton’s second law of motion. The force balance equation is given by

$m_s\frac{\partial u_s}{\partial t}={{\displaystyle \sum }}^{\text{}}{F}_i=F_{gra}+F_{buo}+F_{drag}+F_{pre}+F_{Saf}+F_{Mag}+F_{con}$E18

The individual force terms in Eq. (18) are, in order, the bulk forces due to gravity and buoyancy, the hydrodynamic force due to drag, the pressure force due to the pressure gradient, the Saffman lift force due to fluid shear, and the Magnus force due to particle spin, and the contact force due to particle-particle or particle-wall collision. Given the large difference between the particle and fluid density, the pressure force can be dropped from Eq. (18) without loss of accuracy; the Magnus force can also be dropped because of negligible particle rotation. In this work, this contact force F_con is computed using the soft-sphere model, which decouples its normal and tangential components [22]. The normal force is given by

$F_{con}^n=\left(k\delta +\gamma \left(u_{12}e\right)\right)e$E19

where k is the spring constant of the particle, δ is the overlap between the particle pair involved in the collision as illustrated in Figure 8, γ is the damping coefficient, u₁₂ is the relative velocity vector of the colliding pair, and e is the unit vector. For large values of k, the results of the soft-sphere model are interchangeable with those obtained using a hard-sphere model [35]. The tangential contact force is calculated based on the normal force as ${\mathit{F}}_{con}^t=\mu {\mathit{F}}_{con}^n$ where the coefficient of friction μ is given as a function of the relative tangential velocity ν_r by

$\mu \left(v_r\right)=\{\begin{array}{cc}\begin{array}{l}{\mu }_{stick}+\left({\mu }_{stick}-{\mu }_{glide}\right)\left(v_r/v_{glide}-2\right)\left(v_r/v_{glide}\right)\hfill \\ {\mu }_{glide}\hfill \end{array}& \begin{array}{c}\text{if }{v}_r<v_{glide}\\ \text{if }{v}_r\ge v_{glide}\end{array}\end{array}$E20

Figure 8.

Schematic of particle collision model for DEM.

3.1.2. Solid-gas momentum exchange

The momentum exchange between the solid and fluid phase is a crucial element for modeling multiphase flow using the coupled CFD/DEM solver. The transfer of momentum from the fluid to a solid particle as it moves through each cell in the computational domain is attributed to the drag force given by

$F_{drag}=F_D\left(u_f-u_p\right)$E21

where u_f is the fluid velocity, u_p is the particle velocity, and F_D is the net drag coefficient,

$F_D=\frac{18\mu }{{\rho }_p{d}_p^2}\frac{C_{D}R{e}_p}{24}$E22

where μ, ρ, and d_p are the viscosity of the gas and the density and diameter of the solid particle, respectively. C_D and Re_p are the particle drag coefficient for a sphere and the relative Reynolds number based on the particle diameter respectively. Re_p is defined as

$R{e}_p=\frac{{\rho }_f{d}_p\left|u_f-u_p\right|}{\mu }$E23

The drag coefficient can be modeled using various empirical relations. The spherical or Stokes drag law is chosen in this work for its simplicity.

3.1.3. Parcel conceptcel concept

The computational cost of the DEM approach is driven by the number of collisions between particles; to track each individual particle in a CLC system using the DEM approach is extremely computationally demanding since the total number of particles increases drastically as the particle size becomes smaller. Hence, the parcel methodology first proposed by Patankar and Joseph [36] is employed in this work to overcome the high computational cost.

According to Patankar and Joseph [36], one parcel can represent a group of particles with the same properties such as density and size. The mass used in collisions is that of the whole parcel rather than a single particle. By summing the mass and volume of each individual particle in the parcel, the total mass m_p and volume V_p of the parcel can be obtained. The radius of the parcel is then determined by the mass of the entire parcel and the particle density. For a given point in the fluid flow, the driving force of a parcel due to fluid forces is assumed to be the same as the sum of the fluid force acting on the group of particles it represents.

$F_{f,p}={\displaystyle {\sum }_{i=1}^{N_p}{F}_{f,i}}$E24

where N_p is the number of particles in the parcel, and F_f is the net fluid force acting on a parcel p or particle i depending on the subscript. The acceleration due to inter-particle collision forces and particle-wall collisions forces is computed based on the mass properties of the parcel.

3.2. Cold flow simulation of NETL CD-CLC reactor system

The boundary conditions for the cold flow simulation are obtained from Parker [18] and are summarized in Table 2. Given the high gas velocity required in the air reactor to carry the particles to the top of the reactor and into the cyclone, the flow is turbulent. However, it is well-established that for gas-solid flows, the effect of turbulence is increasingly negligible compared with the effect of the solids for solid volume fractions above 0.001 [37]. For the present simulation of a fluidized bed with densely packed regions, the effect of turbulence can be ignored without loss of accuracy, in line with the work of Parker [18]. The particles are injected into the air reactor, loop-seal, and air reactor and are allowed to settle prior to the start of the simulation. A total of 717,879 particles in total are injected into the system—73,360 particles in the air reactor, 365,057 in the fuel reactor, and 279,462 in the loop-seal. After the particles are settled in each bed, the CFD/DEM model is run to simulate 360 ms of cold flow operation. The development of particle movement is shown in Figure 9 at 10 ms intervals with the particles colored by velocity magnitude.

Unit	Boundary	Gas	Flow rate (m/s)
Air reactor	Fluidizing air	Air	20.0
Fuel reactor	Fluidizing gas	N2	4.0
Loop-seal	Fluidizing gas	N2	2.0
L-valve	Stripper (upper)	N2	0.5
	Aeration (middle)	N2	1.0
	Eductor (lower)	N2	1.0

Table 2.

Boundary conditions for cold flow simulation.

As shown in Figure 9, the particles in the air reactor reach the top of the riser at around 190 ms, and then travel horizontally along the pipe toward the cyclone aided by two secondary gas injections on the side of the air reactor. After another 40 ms, the particles enter the cyclone and start to drop down to the loop-seal. Due to an erroneously high gas velocity in the loop-seal and fuel reactor during the initial startup, the particles in these chambers are also shot up to the top. Once the gas injections are reduced to their correct values of 2 m/s and 4 m/s, respectively at 210 ms, the particles settle down again. From this point onward, it is expected that the particles in the loop-seal will drop into the fuel reactor and the L-valve, and finally be pushed back into the air reactor by the gas injection in the L-valve.

Figure 9.

Particle tracks colored by velocity magnitude for the first 360 ms of cold flow simulation.

Three hundred sixty millisecond is not sufficient time to see the complete particle recirculation; however, the pressure contours in Figure 10 indicate the development of favorable pressure gradients for particle recirculation as the simulation time increases. The static pressure in the system is evaluated at surfaces S-1 to S-6 shown in Figure 7(b) to quantify the pressure gradients observed in Figure 10; the static pressures variation at 360 ms is presented in Figure 11. The arrows indicate the particle movement direction. It can be observed from Figure 11 that there is a consistent positive pressure differential between surface S-1 (air reactor bed) through to S-5 (fuel reactor bed), which confirms that particle continuous recirculation can occur between these surfaces. It is noted that the pressure gradient between S-5 and S-1 via S-6 (L-valve) is an adverse gradient, which explains why the L-valve gas injection feeds seen in Figure 7(a) are required to ensure particle circulation around the loop.

Figure 10.

Pressure contours for cold flow inside the CLC apparatus.

Figure 11.

Static pressure at surfaces S1–S6 in the CD-CLC system shown in Figure 7 at t = 360 ms.

4. Process simulation of calcium looping

The calcium looping process offers a solution for capturing CO₂ from existing power plants. In this section, CaL is modeled at the system level using the process simulation software Aspen Plus. In order to investigate the energy penalty associated with a CaL system, the overall heat production from a power plant without and with calcium looping is determined.

4.1. Simulation of calcium looping with postcombustion capture

In the postcombustion capture scenario, the carbonator and calciner are included downstream of the combustion process to capture CO₂ from the flue gases generated by the combustion of coal. The materials used in the process simulation include conventional and non-conventional components. Pure materials, including all possible chemical compounds comprising the elements C, N, O, H, S, and Cl that might be produced during the chemical reactions, are designated as conventional. Properties for conventional materials are obtained from the Aspen Plus data bank. Mixtures such as coal and ash are designated as non-conventional solids.

4.1.1. Process simulation setupcess simulation setup

The doctoral work of Sivalingam [38] is used as a basis for the process simulation of calcium looping presented in this chapter. Illinois #6 coal is used in the simulation in line with the work of Sivalingam [38]. The RYIELD reactor block is employed in Aspen Plus to decompose the nonconventional material coal into its constituent conventional materials—ash, H₂O, C, H₂, N₂, Cl₂, S, and O₂. Mass percentages for the component yields are set based on the proximate and ultimate analysis of the Illinois #6 coal given in Table 3.

Proximate analysis (wt.%)			Ultimate analysis (wt.%)
Moisture	Volatile matter	Fixed carbon	Ash	C	H	N	Cl	S	O
11.12	34.99	44.19	9.70	80.51	5.68	1.58	0.37	3.17	8.69

Table 3.

Physical and chemical properties of Illinois #6 coal [38].

The outlet stream from the RYIELD reactor goes into a burner, modeled as a RGIBBS reactor, along with air for combustion. The RGIBBS reactor automatically calculates the combustion products at equilibrium such that the Gibbs free energy is minimized. The air flow rate into the RGIBBS reactor is set at the minimum value where the carbon is completely combusted. The calculation for the proper amount of air is discussed in Section 4.1.2. The CO₂-rich flue gas after combustion then undergoes the calcium looping process. The temperature of the flue stream is maintained at 150°C in accordance with the lower limits on power plant flue gas temperatures provided by Feron [39].

The carbonator refers to the reactor where the carbonation reaction takes place. The RSTOIC reactor block is used in Aspen Plus to model the carbonator. The pressure is set at 1 bar and the temperature is set at 650°C. The RSTOIC is a reactor in which the user can define the specific reaction that occurs. The carbonation reaction is given by

${\text{CaO}}_{(\text{s})}+{\text{ CO}}_{\text{2}(\text{g})}\to {\text{CaCO}}_{\text{3}(\text{s})}$E25

In real situations, CaO and CO₂ do not react completely with each other. The amount of CaO that can react with CO₂ is constrained by the surface area of CaO particles and by the extent of the solid-gas mixing as a result of the fluidization behavior in the reactor. These effects are incorporated into Aspen Plus by defining the conversion fraction for one of the reactants, CaO. The dependence of the CO₂ capture efficiency of the carbonator on the sorbent flow ratios is based on the work of Abanades et al. [40]; the data were converted into a table by Sivalingam [38] and is presented in Table 4. $F_{{\text{CO}}_2}$ is the mole flow rate of CO₂, F_R is the mole flow rate of recycled (or looped) CaO, and F₀ refers to the make-up flow of CaO.

$F_0/F_{{\text{CO}}_2}$	$F_{\text{R}}/F_{{\text{CO}}_2}=3$	$F_{\text{R}}/F_{{\text{CO}}_2}=4$	$F_{\text{R}}/F_{{\text{CO}}_2}=5$
0.05	0.63	0.81	0.99
0.10	0.76	0.95	0.99

Table 4.

CO₂ capture efficiency for different flow ratios of CaO and CO₂ [38].

It is not possible to model the make-up flow of CaO as a variable in Aspen Plus. Instead, a fixed of $F_0/F_{{\text{CO}}_2}=0.1$ is chosen with three values of $F_{\text{R}}/F_{{\text{CO}}_2}=0.1$ to model three values of CO₂ capture efficiency in the range of 50–100%. Multiple cases are run in Aspen Plus using different values for the CaO conversion fraction until the correct CO₂ capture efficiency is achieved. As shown in Table 5, each specified CaO conversion fraction corresponds to a range of CO₂ capture efficiencies. The simulation results for the various CaO conversion fractions are plotted in Figure 12. The large symbols in Figure 12 refer to those cases whose results fit the data of Sivalingam [38] exactly; these are denoted as the results obtained from experimental data in the following discussion. The small symbols represent all the trial cases conducted in Aspen Plus and are denoted as extrapolated results.

$F_R/F_{{\text{CO}}_2}=3$	CaO conversion fraction	CO2 capture efficiency
3	0.33	0.66–0.86
4	0.25	0.86–0.97
5	0.20	0.97–0.99

Table 5.

Range of CO₂ capture efficiencies for each CaO conversion fraction.

Downstream of the carbonator, the solids mixture (primarily CaCO₃ with some CaO depending on the inlet flow rate of CaO) and the CO₂-lean flue gas is cooled back to the 150°C and sent to the calciner to regenerate the CaO. Similar to the carbonator, the RSTOIC reactor block is employed for calciner in ASPEN Plus. The calcination reaction that takes place in the calciner is given by

${\text{CaCO}}_{\text{3}(\text{s})}\to {\text{CaO}}_{(\text{s})}+{\text{ CO}}_{\text{2}(\text{g})}$E26

Figure 12.

Range of CO₂ capture efficiencies for various CaO conversion fractions.

The temperature in the calciner is 900°C and the pressure is 1 bar in accordance with Sivalingam [38]. The calcination reaction is a complete reaction, so the conversion fraction of CaCO₃ is set at 1. The reactor blocks used in Aspen Plus for calcium looping with postcombustion capture are summarized in Table 6 along with their functions and reactions; the final flow sheet is shown in Figure 13.

Name	Reactor	Function	Reaction formula
DECOMP	RYIELD	Converts nonconventional solids into conventional	Coal → char + volatiles
BURN	RGIBBS	Burns coal with air	Char + volatiles + O2 → CO2 + H2O
CARBONAT	RSTOIC	Carbonation	CaO + CO2 → CaCO3
CALCINER	RSTOIC	Calcination	CaCO3 → CaO + CO2

Table 6.

Process models used for calcium looping with postcombustion capture setup in Aspen Plus.

The heat of the combustion process (without CaL) is obtained by adding the heat from the decomposer, burner, and heat exchangers for ash and flue gas. The heat from the carbonator, and calciner is the heat of the CaL process. These values of heat and the corresponding CO₂ fraction in the final outlet flow are indicative of the performance of the CaL system with postcombustion capture.

Figure 13.

Aspen Plus flow sheet for calcium looping with postcombustion capture.

4.1.2. Results and discussion

The air stream is simulated as a mixture of 79% N₂ and 21% O₂. The other components of air such as Ar and CO₂ are present in such small fractions that they can be ignored with negligible effect on the results. The optimization module in Aspen Plus is employed to find the flow rates of the O₂ and N₂ flow rates such that the burner heat is maximized. The optimized flow rates of O₂ and N₂ are 4.2 kmol/s and 15.8 kmol/s, respectively for 50 kg/s of Illinois #6 coal. The sum of the flow rates of around 20 kmol/s indicates the total air flow in the burner.

The energy penalty of CaL is the sum of the net gain and loss from the carbonation and calcination processes. For each CaO conversion fraction, there is a corresponding CO₂ capture efficiency from the work of Sivalingam [38]. For a defined CaO conversion fraction in Aspen Plus, the experimental value of the CO₂ capture efficiency is matched by adjusting the CaO flow rate into the carbonator, which changes the heat duty of both the carbonator and calciner. Hence, each experimental data point corresponds to one heat duty value for the carbonator and calciner. Additional data points are obtained by considering the extrapolated data from Figure 12 as well, providing to a range of CO₂ capture efficiencies for each CaO conversion fraction due to the varying CaO flow rates instead of just one value that matches the result of Sivalingam [38]. The heat duty of the carbonator and calciner are plotted in Figure 14. It is noted that the heat duty of the calciner from Aspen Plus is negative since the calcination reaction is endothermic. Figure 14 plots the absolute value of the calciner heat absorption for direct comparison with the heat production in the carbonator. For each capture efficiency, the heat absorbed by the calciner is greater than the heat produced in the carbonator, which confirms that there is a net energy penalty associated with the calcium looping process.

Figure 14.

Heat duty of carbonator and calciner for original experimental data and extrapolated data.

Figure 15.

Plot of total energy output vs. CO₂ capture efficiency without CaL and with postcombustion CaL.

It can be observed from Figure 14 that the calcination results fall on a straight line. This is expected since the calculation is based on a stoichiometric relation: the heat produced is proportional to the inflow rate of the reactant, CaCO₃. The CaO conversion fraction does not affect the heat absorbed by the calciner because the calciner has the same temperature for both the inlet and outlet streams (150°C). Any excess CaO that passes through the carbonator and into the calciner has no effect on the heat duty of the calciner. For the carbonator, the extrapolated data around each of the three experimental data points are linear but the lines segments do not coincide. Each straight line segment corresponding to a range of extrapolated data has a reduced y-intercept compared with the previous section and has a more gradual slope compared with the calciner. From the modeling point of view, the only difference between these two reactors, beside the chemical reaction, is the inlet stream temperature. Since the CaO stream entering the carbonator is an external input, the inlet stream has a temperature of 25°C compared with 150°C for the calciner (internal input from the carbonator outlet stream). Thus, some heat is consumed in the carbonator for heating up the inlet stream to the temperature of the outlet. The heat production of the carbonator decreases as the CO₂ capture efficiency is increased since more heat is consumed to heat the higher CaO flow that is required for the increased capture.

Figure 16.

Plot of energy penalty vs. CO₂ capture efficiency for CaL with postcombustion capture.

For a coal feed of 50 kg/s, the heat of combustion is calculated to be 1168 MW without calcium looping. When the net heat from the carbonator and calciner is added, the total heat of the power plant with calcium looping ranges from 1060 to 1130 MW, as shown in Figure 15.

The energy penalty for CaL refers to the fraction of energy produced by a power station that must be dedicated to the carbonation and calcination process in order to capture CO₂. The energy penalty can be defined as

$Energypenalty=\frac{\left|Q_{\text{looping}}\right|}{Q_{\text{total}}-\left|Q_{\text{looping}}\right|}.$E27

where Q_total is the total heat produced by the power plant and Q_looping is the heat for the CaL process. The energy penalty for CaL with postcombustion capture calculated using Eq. (27) is shown in Figure 16. From the figure, the energy penalty ranges from 3.5% to 9.0% over the corresponding range of CO₂ capture efficiencies from 65% to 99%. These results are in line with the work of Cormos and Petrescu [41] that found energy penalties ranging from 5to 10% for capture efficiencies between 92% and 93% for various power plants. More recently, a study by BP Alternative Energy [42] showed that the marginal energy penalty associated with increasing carbon capture efficiency past 98% increases drastically, which further validates the results of the current work. It is noted that only the heat output of the combustion process and the calcium looping setup is assessed in the results presented above; the steam cycle has not been considered. The energy penalty calculated in this work should be considered as a lower bound for any investigation on calcium looping.

5. Summary and conclusion

This chapter presents numerical simulations of the chemical looping combustion and calcium looping processes used for carbon capture from fossil fuel power plants. A three-dimensional CFD simulation for a complete circulating dual fluidized bed system is developed for chemical looping combustion based on the 300 W laboratory-scale experiments of Abad et al. [17]. The oxygen carrier is modeled as an Eulerian fluid phase based on the kinetic theory of granular flow. The results of this study highlight the importance of capturing the diffusion of gases in a CLC reactor in ensuring that accurate results are obtained for the chemical reactions; the results of the three-dimensional simulation are a better match for the outlet concentrations of the gases recorded in the experiment than a two-dimensional simulation previously conducted. It is expected that the simulation accuracy can be increased further if empirical reaction rate data becomes available for the specific oxygen carrier and fuel used in the experiment.

The detailed particle hydrodynamics in a complete circulating CLC system for solid fuels is investigated using the Lagrangian particle-tracking approach known as DEM based on the CD-CLC reactor at NETL [18]. The development of particle flow is investigated as well as the pressure contours in the reactor. Although the simulation time was not sufficient to see the complete particle recirculation in the system, the static pressure contours and gradients showed evidence of favorable conditions for particle recirculation. The cold flow simulation in this work provides a foundation for more detailed simulations of CD-CLC systems in future work where the coal particles and the associated chemical reactions can be taken into account.

The energy penalty associated with the calcium looping process for postcombustion capture of CO₂ is investigated using Aspen Plus; the results indicate that the energy penalty depends on the flow rate of CaO into the carbonator. An important finding from this work is that for CO₂ capture efficiencies above 90%, the marginal energy penalty associated with any further increase in efficiency increases drastically. This suggests that there is a limit on the maximum CO₂ efficiency possible from calcium looping beyond which the process becomes impractical due to the energy consumption.