A Multiscale Approach for Modeling Oxygen Production by Adsorption.

— A Multiscale Approach for Modeling Oxygen Production by Adsorption — Oxygen production processes using adsorbents for application to CCS technologies (Carbon Capture and Storage) offer potential cost beneﬁts over classical cryogenics. In order to model adsorption

processes an approach using three size scales has been developed. This work is being conducted in the framework of the DECARBit European research project. The first scale is at the size of the oxygen adsorption bed to be modelled as a vertical cylinder filled with pellets. Its length is 0.2 m (scale 10 À1 m). The bed is homogeneous in the transversal direction so that the problem is 1D (independent variables t, x). The physics in the process include gas species (C bk (t, x)) convection and dispersion, thermal convection and conduction (T(t, x)) and hydrodynamics (v(t, x)). The gas constituents involved are N 2 , O 2 , CO 2 and H 2 O. The second scale is at the size of the pellets that fill the adsorber and which are assumed to be of spherical shape with a typical radius of 5 mm (scale 10 À3 m). The independent variable for the pellets is the radius "r p ". At a certain height (x) down in the adsorber all the pellets are the same and are surrounded by the same gas composition but inside the pellets the concentrations may vary. The state variables for the inner part of the pellets are the gas concentrations C pk (t, x, r p ). The pellets are so small that they are assumed to have a uniform temperature. This leads to a 2D transient model for the pellets linked to the 1D transient model for the bulk. The third scale looks into the detailed structure of the pellets that are made of perovskite crystallites. The latter are assumed to be spherical. Oxygen adsorption occurs in the crystallites which have a radius of about 0.5 lm (scale 10 À7 m). All the crystallites at the same radius in a pellet are supposed to behave the same and because they are spherical, the only independent variable for a crystallite located at (x, r p ) is its radius "r c ". The state variables for the crystallites are then the adsorbed oxygen concentration C c2 (t, x, r p , r c ). The crystallites are so small that they are assumed to have a uniform temperature. This leads to a third transient model that is 3D for the crystallite and is linked to the 2D transient model for the pellets which is itself linked to the 1D transient models for the bulk. From the larger to the lower scales, the links between the three models are the following: the bulk concentration and temperature give the boundary conditions surrounding the pellets. The pellet concentration gives the boundary conditions for the crystallites. We chose to solve this multiscale approach that requires the coupling of models of different dimensions in Comsol Multiphysics. The simulator was built to gain knowledge from laboratory experiments in order to estimate whether oxygen separation from air is realistic or not.

C bk
Bulk concentration (mol/m 3 ) C c2 Crystallite oxygen concentration (oxygen #2 only) (mol/m 3 ) C c2m Maximum adsorbed oxygen concentration at equilibrium at PO 2 (mol/m 3 ) C c2M Maximum adsorbed oxygen concentration at infinite oxygen pressure (mol/m 3 ) C pg Gas thermal capacity (J/mol/K) C pk Pellet concentration (mol/m 3 ) C s Pellet thermal capacity (J/mol/K) D c2 Diffusion coefficient for the adsorbed oxygen in the crystallites (m 2 /s) D pk Dispersion of species ( Gas velocity (m/s) V p Volume occupied by the pellets in one cubic meter of the porous media (m 3 /m 3 ) x

INTRODUCTION
Processes generally involve porous media. They are natural in the case of reservoir characterization or can be artificial when heterogeneous catalysts are involved. Porous media are decimetric to metric in size for industrial and lab process but knowledge deals with millimetric or nanometric length scale for pore sizes, wettability or even more specific in our case, adsorption ability and diffusibility in perovskite [1]. Numerical simulators are of course representative of the highest scale [2,3].
To use low level information at the highest scale, some authors condensate lower level knowledge in some kind of overall parameters that could be calculated (Thiele modulus) or estimated on experiments or given by known correlations [2,4].
The Thiele modulus is a calculated parameter. It consists in an analytical solution of the diffusion from the bulk to the catalyst pellets. Assumptions make possible to solve analytically the diffusion and to derive a global parameter (the Thiele modulus) that takes diffusionnal limitation into account as a global parameter used at the reactor level.
Artificial millimetric networks [5] or structured catalytic reactors [6] are also designed to be representative of actual porous media. The main advantage is that their geometrics are simple and fully controlled and large scale parameters can be derived from these networks through numerical simulations. Another way to control the geometrics is to finely measure the solid matrix structure in order to reconstruct the porous media. Solid matrix can be described by series of images [4] or by high resolution computed micro tomography [7]. Large scale characteristics such as capillary pressure and relative permeabilities can then be calculated [7]. Fractal geometrics were very popular in the 90's to describe porous media geometry [8]. The difficulty with fractals was to link the large scale parameter to the fractal decimal dimension that characterizes the pattern. Barthe´le´my et al. [9] mixed fractal and statistical descriptions of fractured geological layers to derive hydraulic properties.
Another way to integrate multiscales in porous media is to perform homogenisation techniques on real equations and variables at the pore scale. Marle [10] developed a complete theory to derive generalized twophase flow in porous media based on the homogenisation of the Navier-Stocks equations and variables. Kalaydjian [11] worked on Marle's equations towards two-phase Darcy law in porous media while Pavone [12] completed Marle's equations by the derivation of a macroscopic capillary pressure equation for porous media based on wettability and on solid-liquid areas. Similar researches were performed by Quintard and Whitaker [13] by the same time.
The paper presented hereafter is different in the sense that it solves different length scales simultaneously in a 1D + 2D + 3D numerical simulator. A 2D approach involving Knudsen diffusion in fuel cells is also described in [14]. In this case multiscale is solved by a 2D model for which one dimension is 1.0 9 10 À3 m and the other is 4.0 9 10 À1 m. Following the same idea, Ingram et al. [15] solve a catalytic packed bed reactor using a discretisation scheme of N pellets exchanging with a bulk solved in N meshes. Discretisation allows to keep 1D models for the pellets whereas we solve a 3D model in a continuous approach.
The approach is applied to an oxygen production process using solid sorbents for application to CCS technologies (Carbon Capture and Storage). These processes can offer potential cost benefits over classical cryogenics. The Ceramic Autothermal Recovery (CAR) process is one such oxygen production system investigated by BOC that uses pressure swing adsorption principles.
Here an alternate technology is being researched that implements the sorbents in a rotating bed and uses the oxygen partial pressure differences between air and a sweep gas for collecting the oxygen.

MODEL GEOMETRY
One can have knowledge on adsorption and desorption at deferent size levels. Of course the adsorber can be described at its own level with its length, injection area, flow rates, etc.
At a lower scale, the pellet one, other information can be added like the pellet size, their porosity, internal diffusion coefficient, etc.
And finally, because a pellet is made of a huge amount of perovskite crystallites aggregated together, information like the mean crystal size can be useful. If it occurs that additional information are available for the perovskite crystallites such as the oxygen diffusion parameters and/or the oxygen adsorption isotherm, these should be integrated at this level. This is the reason why we need to describe the adsorber at three different scales, from the bulk to the crystal going through the pellet one as an intermediate scale.

The Bulk Sorbent Bed
The fixed bed adsorber is a vertical cylinder filled with pellets ( Fig. 1). On the absorber scale, the bed is considered globally homogeneous at all axial (vertical) positions such that the radial (horizontal) gradients are zero. The state variables for the bulk gas at the sorbent bed length scale are the gas concentrations C bk t; x ð Þ, the pressure P t; x ð Þ, the temperature T t; x ð Þ and the gas velocity v t; x ð Þ. State variables depend on the time "t" and on the axial position "x". This leads to a 1D transient model for the bulk.

The Pellets
The pellets are actually cylindrical but will be assumed to be spherical so as to model them with a 1D spherical symmetry approach. The independent variable for the pellets is their radius "r p ". The state variables for the pellets are the concentrations C pk t; x; r p À Á in the pellets and the temperature T p t; x; r p À Á . Because the pellets are relatively small compared to the adsorber size, the temperature in the pellet is assumed to be homogeneous and equal to the bulk temperature (T p t; x; r p À Á ¼ T t; x ð Þ). Pellet state variable C pk depends on the time "t", on the vertical position "x" and on the position inside the pellet "r p ". This leads to a 2D transient model for pellets.
There is no adsorption or gas reaction at the pellet scale, only gas transport occurs within the porosity of the pellets.

The Crystallites
As presented in Figure 2, a pellet is made of crystallites that are assumed to be spherical as well. The independent variable for the crystallites is their radius "r c ". The state variables for the crystallites are the absorbed oxygen concentrations C c2 t; x; r p ; r c À Á and the temperature T c t; x; r p ; r c À Á . The "#2" in C c2 stands for oxygen that is referred as #2. Because the crystallites are very small compared to the adsorber size, the temperature in the crystallites is assumed to be homogeneous and equal to the bulk temperature (T c t; x; r p ; r c À Á ¼ T t; x ð Þ). The state variable C c2 depends on the time "t", on the axial position "x", on the crystallite position inside the pellet "r p " and on the position inside the crystallite itself "r c ". This leads to a 3D transient model for crystallites.
x C bk (t,x) C pk (t,x,r p ) r p Reactor Pellets x Figure 1 From the sorbent bed to the pellet geometry.   From pellet geometry to crystal geometry.

State Variables
The absorber is perfectly characterized by its state variables which are time (t) and space (x, r p , r c ) dependent. These state variables are: -species (#k) concentration in the bulk: C bk t; x ð Þ (mol/m 3 ), -species (#k) concentration in the pellet pores: C pk t; x; r p À Á (mol/m 3 ), -adsorbed O 2 (#2) concentration in the crystallites: C c2 t; x; r p ; r c À Á (mol/m 3 ), -temperature in the catalyst: T t; x ð Þ(K), -superficial gas velocity in the bulk: v t; x ð Þ(m/s), -gas pressure: P t; x ð Þ(Pa).

Components
The components are the air components (mainly N 2 and O 2 ) plus a gas for desorption that can be condensed to be separated from oxygen with two possibilities: CO 2 and H 2 O.

Mass Balance Equations for the 1D Bulk
The bulk gas transport in the bed is convective and diffusive with a source term (S bk ) that derives from the gas transfer between the bulk and the pellets: where C bk (t, x) is the (k) species concentration in the bulk gas, e b is the bulk porosity and D bk is the dispersion of species (k). The molar sources in the bulk flow (S bk ) occur at the external surface of the pellets and is due to adsorption/desorption processes within the pellets. The model assumes that the molar flux is proportional to the difference between the concentration of a species in the bulk phase (C bk (t, x)) and at the surface of the pellets (C pk (t, x, R p )): where R p is the radius of the pellet and k pk is the proportionality factor. The bulk source term can be derived and calculated as the flux summed over the surface of the pellets: S bk ¼ Àn #p k pk C bk t; x ð Þ À C pk t; x; where n p (in #/m 3 ) is the number of pellets per cubic meter and is derived from the bulk porosity: Here V p (in m 3 /m 3 ) is the volume occupied by the pellets in one cubic meter of the porous media. By substitution the source term becomes: The boundary conditions are standard Dirichlet-Neumann conditions driven by the flux at the entrance and convective flux at the exit: where Q in,k is the injected mol flux of species (k):

Mass Balance Equations for the 2D Pellets
At a distance "x" from the absorber entrance, it is assumed that the perovskite pellets are spherical and surrounded by an homogeneous bulk gas mixture. Inside the pellets, the species transport is only diffusive. The mass balance equations for the gas phase species are: e p @C pk @t þ 1 r 2 p @ ÀD pk r 2 p @C pk @r p where C pk is the (k) species concentration in the pellet, e p is the pellet porosity and D pk is the dispersion of species (k). S pk is the source terms of species (k) that derives from the oxygen adsorption and desorption in the crystallites. As oxygen is #2, S pk ¼ 0 for k #2: where n #c (in #/m 3 ) is the number of crystals per cubic meter of pellets, F a is the adsorption flux, F d is the desorption flux and R c is the crystallite radius: and The adsorption flux (F a ) is positive in the case of adsorption and null in the case of desorption. The adsorption flux tends towards zero when equilibrium is reached. Inversely, the desorption flux (F d ) is positive in the case of desorption and null in the case of adsorption. The desorption flux tends towards zero when the equilibrium is set as well: where k a and n a are constant, C c2 is the adsorbed oxygen concentration and C c2m is the adsorbed oxygen concentration that would be at equilibrium with the oxygen partial pressure in the pellet (P O 2 ).
C c2M is the maximum adsorbed concentration that the medium can hold under an infinite oxygen partial pressure and K cm is a constant together with n c . For the desorption flux: where K d and n d are constant, P c2m is the oxygen partial pressure that would be in equilibrium with the oxygen adsorbed concentration in the pellet: The boundary condition at the surface of the pellets simulates the conservation of gas flux from the bulk to the pellets. Because the length scales of the bulk and the pellets are not the same, the gas concentration can not be taken as continuous from one scale to the other. To account for that the model introduces a jump in the concentrations and a mol flux proportional to it that corresponds to the source term of the bulk phase.
@C pk @r p r p ¼0 ¼ 0 and ÀD pk @C pk @r p r p ¼R p ¼ k pk C bk t; x ð Þ À C pk t; x; R p À Á À Á

Mass Balance Equations for the 3D Crystallites
It is assumed that a perovskite crystallite within a pellet is spherical and surrounded by a homogeneous gas mixture. Inside the crystallites, adsorbed oxygen can transfer by some kind of ionic diffusion here described with a standard Fick's law. Ionic diffusion between crystallites is not considered. The mass balance equation is the following: where e c is a volume fraction that can be occupied by the adsorbed oxygen and D c2 is the diffusion coefficient for the adsorbed oxygen in the crystallites. The boundary conditions are standard where it is assumed that there are two oxygen fluxes at the crystallite boundary, one for adsorption and one for desorption:

Molar Balance in the Bulk, Pellets and Crystallites
The molar balance Q k can be calculated with: where Q pk is the accumulated mol in one pellet: where Q ck is the accumulated mol in one crystallite: The total accumulated mol for component #k in the three domains is: where S reac is the adsorber surface.

Thermal Equation
As pellets are small, the simulator assumes that the pellets are homogeneous in temperature and that this temperature is the same as for the bulk at the same spacial position. As a result the temperature just depends on two independent variables (x,t): where C pg is the bulk gas thermal capacity, C s is the pellet thermal capacity and D T is the bulk averaged thermal conductivity. Q b is the energy source derived from adsorption and desorption fluxes and the heat of adsorption or reaction. The boundary conditions are:

Hydrodynamics
The hydrodynamics just concern the bulk flow. The model uses the Darcy equation that links the flow velocity with the pressure head: K d is the medium permeability and l g is the gas viscosity.

Transient Results
Material is a LSCF (Lanthanum Strontium Cobalt Ferrite) perovskite investigated experimentally in the framework of the European funded DECARBit project. The experiments allowed to estimate the model parameters for the simulations presented here. The calculations that are presented correspond to cyclic adsorption and desorption phases lasting 45 s each. The presented simulation lasts 45 s for desorption and the same for adsorption. The results are presented in Figure 3 and Figure 4. At first, the medium is saturated with air then at time zero and during 45 s, steam is injected (Fig. 3). As a result, oxygen is desorbed and released from the medium (blue curve in Fig. 4). Then for a new time period of 45 s, oxygen is injected at the same amount of 0.21 mol/mol. Normally 0.79 mol/ mol N 2 should be injected as well but to give a reference, CO 2 is "numerically" injected as a non adsorbed gas, like a tracer, to show what would be the oxygen level in the case of no adsorption. Hence injection is 0.21 O 2 , 0.21 CO 2 (acting as N 2 ) and 0.58 N 2 .
In Figure 4, we can see the adsorption and desorption effects looking at the difference between the blue curve (O 2 ) and the green curve (CO 2 acting as a inactive gas).
Looking at the oxygen amount remaining in the perovskite can give interesting insight into the process. Figure 5 presents the amount of oxygen still trapped in the perovskite compared to the total amount at the beginning of the desorption. In this example, the figure shows that the process deals with less than 15% of the adsorbed oxygen.

Instantaneous Results
It is interesting to look at three instantaneous snap shots. The first one gives the oxygen content in the bulk gas in the adsorber and is 1D according to the model construction (Sect. 2.3). The second one gives the oxygen content in the gas in the pellets and is 2D according to the model construction (Sect. 2.4). Finally the third one presents the amount of oxygen adsorbed in the crystallites and is 3D (Sect. 2.5). This last is a little more complex to analyze than the two others.

1D Oxygen Mol Fraction in the Bulk at the Adsorber Scale
At the highest scale, the bulk concentrations are quite easy to understand because they are 1D. Figure 6 presents oxygen concentrations versus adsorber length at three different times: one at the end of air injection and the two others just after the switching from air to steam. The first curve, the upper one, shows that injected oxygen mol fraction is 0.21 and decreases in the adsorber due to oxygen adsorption. The two lower curves show that no more oxygen is injected but due to oxygen release, the bulk still contains oxygen with a decreasing amount with the time. Figure 7 gives the oxygen mol fraction in a pellet in the adsorber (2D map at bottom). The horizontal axis gives the position of the pellet in the adsorber, from the inlet (x = 0 at left) to the outlet (x = 0.2 at right). The vertical axis gives the oxygen mol fraction inside the pellet located a`x. In the example, the 1D graph gives the oxygen mol fraction in the gas for a pellet located at (x = 0.01) from the centre of the pellet to the boundary.

Adsorbed Oxygen at the Crystallite Scale
The 3D graph presented in Figure 8 is complex. Surprisingly, the adsorber is cylindrical but the graph is cubic! Actually, this graph gives the adsorbed oxygen (in mole): -in the crystallites located, -in the pellets, -in the adsorber.
The vertical dimension is the adsorber length, going from the inlet (x = 0 m) at the bottom to the outlet (x = 0.2 m) at the top. The location of the crystallite in the pellet is given by the horizontal green axis (from R p = 0 to R p = 5 9 10 À3 m). The crystallites located on the pellet surface are presented on the vertical side at the right of the cube. The ones at the centre of the pellets are on the opposite hidden surface. And finally the graph gives also the concentration in the crystallite (from R c = 0 to R c = 0.5 9 10 À6 m) along the blue horizontal axis.
The figure presents the adsorbed oxygen just after the beginning of a new air injection phase. Air is injected Oxygen mol fraction stored in the perovskite. Oxygen in the bulk at the adsorber scale.
upwards at the visible angle (dark red angle in Fig. 8).
We can see on the figure the diffusion of the adsorbed oxygen in the three directions.

CONCLUSIONS
A multi-scale simulator is presented that allows to account for a metric adsorber (2 9 10 À1 m) filled with millimetric pellets (5 9 10 À3 m), each pellet being an agglomerate of micrometric crystallites (5 9 10 À7 m). The simulator meshes and solves these three length scales at a time in a 1D model plus a 2D model plus a 3D model, all linked and coupled together. This simulator has been adapted to match lab experiments and can be used as a tool to test the effects of microscopic parameters such as the crystallite size and diffusivity on oxygen adsorption and desorption.   Oxygen profile in a pellet close to the adsorber entrance.