Oxygen Transport Membranes for Efficient Glass Melting

Glass manufacturing is an energy-intensive process in which oxy-fuel combustion can offer advantages over the traditional air-blown approach. Examples include the reduction of NOx and particulate emissions, improved furnace operations and enhanced heat transfer. This paper presents a one-dimensional mathematical model solving mass, momentum and energy balances for a planar oxygen transport membrane module. The main modelling parameters describing the surface oxygen kinetics and the microstructure morphology of the support are calibrated on experimental data obtained for a 30 μm thick dense La0.6Sr0.4Co0.2Fe0.8O3-δ (LSCF) membrane layer, supported on a 0.7 mm porous LSCF structure. The model is then used to design and evaluate the performance of an oxygen transport membrane module integrated in a glass melting furnace. Three different oxy-fuel glass furnaces based on oxygen transport membrane and vacuum swing adsorption systems are compared to a reference air-blown unit. The analysis shows that the most efficient membrane-based oxyfuel furnace cuts the energy demand by ~22% as compared to the benchmark air-blown case. A preliminary economic assessment shows that membranes can reduce the overall glass production costs compared to oxyfuel plants based on vacuum swing adsorption technology.


Introduction
The term glass generically refers to high-viscosity substances, which do not possess enough mobility during solidification. Thus, they do not have sufficient time to form a regular crystal lattice, resulting in an intermediate metastable stage marked by the formation of an amorphous solid, which lacks the geometrical order typical of crystal lattices [1]. Common raw materials for the glass production are silica (SiO 2 ), sodium carbonate (Na 2 CO 3 ), limestone (CaCO 3 ) and, dolomite (CaMg(CO 3 ) 2 ). However, the actual composition varies significantly, depending on the glass type. These raw materials are heated in a furnace to high temperatures to form molten glass. The molten glass undergoes further processing to originate the final product [2].
A furnace operates glass melting at very high temperature (>1400-1500 • C) in a confined space surrounded by refractory material. The heat input is generally supplied by burning oil or natural gas in air. Electric booster systems can be installed to supply additional energy to increase the melting capacity. The conventional arrangement of a glass furnace includes a regenerative heat exchanger consisting of two refractory chambers. One chamber is used to store heat from the flue gas, at the same time the other chamber is used to preheat the combustion air. The paths of flue gas and combustion air through the regenerator chambers are interchanged every 15-20 min. The raw material is continuously introduced in the furnace through a feeding port, while molten glass is continuously removed from the furnace [3,4].
The theoretical energy requirement consists of the endothermic heat for the glass reaction, the sensible heat for glass heating and the sensible heat for batch gases (gases from glass reaction). According to the theoretical heat requirement reported by Sardeshpande et al. [5], the soda-lime glass (composed mainly of silicon dioxide (71-75%), sodium oxide (12-16%) and calcium oxide (10-15%)) needs 2671 kJ kg −1 of glass, if no cullet (i.e., crushed glass ready to be re-melted) is used. On the other hand, the theoretical minimum is 1886 kJ kg −1 if the furnace operates on 100% cullet. Such a reduction is due to the absence of the reaction heat of the virgin raw material (487 kJ kg −1 ) and the sensible heat of batch gases (298 kJ kg −1 ).
The use of cullet is beneficial since it is cheaper than virgin raw material and allows energy savings. Moreover, the substitution of raw materials by cullet reduces CO 2 process emissions, which are released by the decomposition of carbonates in raw materials. In the industrial practice, external cullet is extensively used only in the container glass and glass wool production [2].
The above mentioned theoretical heat requirements correspond to the lowest thermodynamic limits [5] but even well insulated glass furnaces will be affected by some structural heat losses. In addition, fossil fuel firing is always associated with flue gas production. Although the flue gas heat content can be partly recovered by combustion air and batch preheating, flue gas is released at temperature of at least 200-250 • C. Thus, the minimum achievable energy consumption level is unavoidably greater than the one reported by Sardeshpande et al. [5]. Energy benchmarking for glass furnaces has been thoroughly discussed by Beerkens [6][7][8]. Referring to the most recent paper [6], Figure 1 shows the ranking of the energy efficiency of 168 container glass furnaces, starting from the lowest specific energy consumption (most energy efficient) up to the highest specific energy consumption. The primary energy values are normalized for the case of 50% recycling cullet in the batch and account for the contribution energy required for electric power generation in case of electric boosting and oxygen production. According to the data in Figure 1, the minimum achievable energy consumption level for normal batch would be around 3.8 MJ kg −1 , which reduces to around 2.5 MJ kg −1 for 100% cullet.
Membranes 2020, 10, 2 of 32 material is continuously introduced in the furnace through a feeding port, while molten glass is continuously removed from the furnace [3,4]. The theoretical energy requirement consists of the endothermic heat for the glass reaction, the sensible heat for glass heating and the sensible heat for batch gases (gases from glass reaction). According to the theoretical heat requirement reported by Sardeshpande et al. [5], the soda-lime glass (composed mainly of silicon dioxide (71-75%), sodium oxide (12-16%) and calcium oxide (10-15%)) needs 2671 kJ kg −1 of glass, if no cullet (i.e., crushed glass ready to be re-melted) is used. On the other hand, the theoretical minimum is 1886 kJ kg −1 if the furnace operates on 100% cullet. Such a reduction is due to the absence of the reaction heat of the virgin raw material (487 kJ kg −1 ) and the sensible heat of batch gases (298 kJ kg −1 ).
The use of cullet is beneficial since it is cheaper than virgin raw material and allows energy savings. Moreover, the substitution of raw materials by cullet reduces CO2 process emissions, which are released by the decomposition of carbonates in raw materials. In the industrial practice, external cullet is extensively used only in the container glass and glass wool production [2].
The above mentioned theoretical heat requirements correspond to the lowest thermodynamic limits [5] but even well insulated glass furnaces will be affected by some structural heat losses. In addition, fossil fuel firing is always associated with flue gas production. Although the flue gas heat content can be partly recovered by combustion air and batch preheating, flue gas is released at temperature of at least 200-250 °C. Thus, the minimum achievable energy consumption level is unavoidably greater than the one reported by Sardeshpande et al. [5]. Energy benchmarking for glass furnaces has been thoroughly discussed by Beerkens [6][7][8]. Referring to the most recent paper [6], Figure 1 shows the ranking of the energy efficiency of 168 container glass furnaces, starting from the lowest specific energy consumption (most energy efficient) up to the highest specific energy consumption. The primary energy values are normalized for the case of 50% recycling cullet in the batch and account for the contribution energy required for electric power generation in case of electric boosting and oxygen production. According to the data in Figure 1, the minimum achievable energy consumption level for normal batch would be around 3.8 MJ kg −1 , which reduces to around 2.5 MJ kg −1 for 100% cullet. Example of ranking of energy efficiency from most efficient to lowest efficient container glass furnace of a set of 168 furnaces. Each point refers to one existing container glass furnace. Specific energy data are normalized to 50% cullet and the primary energy consumption takes electricity and oxygen production into account (adapted from Reference [6]).
Oxygen-firing is an effective solution for saving energy in glass furnaces. In the Beerkens' benchmark study [8], the average specific energy consumption (4.3 MJ kg −1 ) of oxygen-fired container glass furnaces is lower than that of new end-port fired regenerative furnaces (4.6 MJ kg −1 including the energy demand for oxygen production) using 50% cullet. Oxygen-fired furnaces with cullet and/or batch preheating also showed values in the range 3.4-3.6 MJ kg −1 .
Different solutions can be implemented to supply oxygen to a glass furnace. In the conventional scenario, oxygen is delivered by a pressure swing adsorption (PSA)/vacuum swing adsorption (VSA) Figure 1. Example of ranking of energy efficiency from most efficient to lowest efficient container glass furnace of a set of 168 furnaces. Each point refers to one existing container glass furnace. Specific energy data are normalized to 50% cullet and the primary energy consumption takes electricity and oxygen production into account (adapted from Reference [6]).
Oxygen-firing is an effective solution for saving energy in glass furnaces. In the Beerkens' benchmark study [8], the average specific energy consumption (4.3 MJ kg −1 ) of oxygen-fired container glass furnaces is lower than that of new end-port fired regenerative furnaces (4.6 MJ kg −1 including the energy demand for oxygen production) using 50% cullet. Oxygen-fired furnaces with cullet and/or batch preheating also showed values in the range 3.4-3.6 MJ kg −1 .
Different solutions can be implemented to supply oxygen to a glass furnace. In the conventional scenario, oxygen is delivered by a pressure swing adsorption (PSA)/vacuum swing adsorption (VSA) Membranes 2020, 10, 442 3 of 32 unit (which represents the most convenient technology in relation to oxygen demand and purity) or cryogenic distillation processes. An emerging alternative for oxygen separation is represented by membranes based on mixed ionic and electronic conducting (MIEC) materials that promise to reduce the energy consumption of oxygen separation compared to conventional technologies. Since they operate at temperature around 1000 • C, an effective thermal integration in high-temperature processes is a key aspect to achieve an efficient oxygen separation. Glass melting, that takes advantage from oxygen combustion and it is operated at temperature significantly above 1000 • C, is one of the most interesting industry for application of these membranes. This paper investigates this topic by an extended approach ranging from implementation and calibration of a newly developed model for performance evaluation of membrane modules, to heat and mass balance simulation of glass melting furnaces based on alternative technologies, to a comparative economic assessment.

Review of Oxygen Transport Membranes
The main purpose of Oxygen Transport Membranes (OTM) is to separate O 2 throughout a non-galvanic electrochemical mechanism occurring within an O 2 -ion transport dense membrane, that is, without the application of an external voltage and the creation of an electric current in an external circuit. The driving force of the process is solely dependent on the difference of chemical potential of O 2 between the two sides of the membrane (i.e., O 2 partial pressure driven OTM) [9]. Sunarso et al. 2008 [10] reviews dense ceramic membranes for O 2 transport. The study identifies the main characteristics of fast oxygen-ion conducting material: (i) fast oxygen-ions exchange at the relevant gas/solid interfaces (i.e., high ionic conductivity in the temperature range of interest for the final application); (ii) high catalytic activity for the Oxygen Reduction Reaction (ORR); (iii) adequate microstructure morphology (i.e., particle and pore size, porosity and tortuosity) to promote O 2 -ions diffusion throughout the material pores; (iv) chemical and physical stability in oxidizing and reducing atmosphere during fabrication and operation. Contrary to galvanic systems in which an external electrical circuit is present, OTMs must ensure a good electrons transport within the dense active material to respect the overall charge neutrality criteria (i.e., simultaneous flux of all charge carriers-ions and electrons-involved in the electrochemical reactions) [11].
In Mixed Ionic-Electronic Conductor (MIEC)-based membranes, charge transport processes are characterized by an increased mobility of electron holes (h • ), other than only ion vacancies (V •• O ), that is, ionic transference number t i < 1. Therefore, MIEC materials are characterized by both high oxide-ion and electron conductivities. This can either be achieved by employing single-phase materials with mixed conducting properties, such as the oxygen-deficient perovskite structure La 1−x Sr x Co 1−y Fe y O 3-δ (LSCF) [12] or by dual-phase membranes, which combine-in a percolated network-materials that show good ionic and electronic conduction properties [13]. In the former case, charged particle transport occurs in the dense membrane bulk without the necessity of a dense Triple Phase Boundary (TPB), since the active sites are distributed across the thickness of the membrane. Consequently, it is possible to speak of Double Phase Boundary (DPB)-i.e., point of contact between the active dense phase and the gas phase-when dealing with MIEC materials [14]. In case a MIEC material is employed, the following requirements should be verified simultaneously: (i) high O 2 -ion conductivity to supply the membrane surface with oxygen vacancies as soon as they are occupied by fresh new O X 0 ; (ii) high electronic conductivity to ensure that electron holes diffusion does not become a rate limiting step in the ion transport process, as also reported by Chen et al., 1995 [15].
From a manufacturing perspective, it is crucial to be able to lay extremely thin dense membrane layers (i.e., 10-50 µm-scale) as this would increase the permeation performance of the membrane proportionally to the decrease in its thickness [16]. This requirement contrasts with the thermo-mechanical resistance capabilities of the membrane when this is subjected to pressure and temperature gradients. Typical industrial applications-such as oxy-combustion [17,18], carbon capture [19] and O 2 separation [20]-require the operation at high pressure and temperature; therefore, membrane's degradation and thermo-mechanical resistance represent important features for a successful process integration. For this reason, ceramic or metallic porous support structures are usually employed to operate as structural resistance, possibly without impeding or hindering the diffusive transport of gaseous species [21].
This work aims to study the integration of OTM for production of pure oxygen and enabling oxy-fuel combustion into the glass industry. Herein we propose to study the thermodynamic performance and economic viability of glass melting furnaces with oxy-fuel combustion via OTM to achieve lower specific energy consumption and higher energy density compared to conventional glass melting furnaces. The first section of this study focuses on the developing of a new 1-D finite volume model of a co-flow and counter-flow planar oxygen transport membrane. The model is then used to design complete OTM modules for oxygen separation and integration into different glass furnace configurations. The second half of this work describes the system thermodynamic analysis and economic assessment of OTM-based glass furnaces. The model is in-house developed at Politecnico di Milano in a Fortran 90 environment. A previous version of the model has already been employed for designing reactive OTMs for small scale H 2 production from natural gas [22]. The model is calibrated against experimental data obtained from the FP7 European Project GREEN-CC using a 30 µm thick dense LSCF membrane layer, supported on a 0.7 mm porous LSCF structure. The model solves the mass, momentum, and energy balances in the membrane's channels, which are subdivided into finite control volumes by a computational grid. The main outcome of the model is the amount of permeated oxygen given the user-defined geometrical (e.g., number of channels, membrane layer and support thicknesses) and thermochemical (e.g., temperature, pressure, composition) inputs in the membrane reactor.
The novelties of this study are: • A new 1-D OTM model is developed to predict the thermodynamic performance and design of planar membrane modules; • The OTM model is calibrated using newly published experimental data from a 30 µm thick dense LSCF membrane; • The designed OTM modules are integrated into three new glass furnace layouts with oxy-fuel combustion; • The proposed layouts are assessed from an energy and economic viewpoint to select the optimal configuration.

Modelling of Supported OTM
The membranes considered in this work are characterized by planar channels with rectangular cross-sectional area. A single membrane channel couple (i.e., feed and sweep sides) is considered in the simulations. Each channel is discretized into control volumes by a uniformly spaced computational grid.

OTM Module Model
The permeation model, the fluid-dynamic model and the heat transfer model are included into a one-dimensional (1D) OTM module model, which aims at designing a complete membrane stack. The OTM model calculates the total surface area required to supply a given oxygen demand and then subdivides the area into layers, stacks, and modules. Various parallel channels constitute each membrane layer, which then form a stack unit. Additionally, depending on the type of application and its design requirements (e.g., the acceptable pressure losses), the integration in series or in parallel of different stacks makes an OTM module.
The development of the OTM model is based on the following assumptions: • One-dimensional analysis: numerical discretization in finite volumes along the axial directory of the streams. In each finite volume, the permeation model, the fluid-dynamic models are solved • Steady state analysis • Ideal gases and ideal mixtures The OTM module model is then used to design membrane modules for glass furnaces and evaluate their performance in various operating conditions-i.e., with and without sweep gas, atmospheric and pressurized operation.

Fluid-Dynamic and Heat Transfer Model
The fluid-dynamics model estimates the viscous friction, which determines the pressure losses within the channels in the feed and permeate compartments. The model accounts for both laminar and completely developed turbulent flow regimes; the former regime is modelled using an interpolation function of tabulated data taken from Reference [23], valid for different channel cross sectional shapes and height/width ratios. The latter is modelled with the Petukhov's correlation: The pressure losses for viscous effects are evaluated via the Darcy-Weisbach equation: where, L is the total length of the considered finite volume, whilst u and ρ are considered at the volume inlet, assuming them constant within each node. Finally, the steady state one-dimensional mass, momentum and energy balances are evaluated per each control volume for both feed and permeate sides.
where, i represents either the feed, the permeate, the membrane or the support and j represents the potential chemical reactions. Since no chemical reactions are present r j terms are null. Note that the heat exchange mechanisms amongst the different parts of the membrane are included as boundary conditions. As far as the heat transfer modelling is concerned, no thermal conduction occurs in the axial direction of the membrane reactor; in other words, the Péclet's number (Pe = Re Pr D h ) is always higher than 1000 for turbulent flow conditions. On the other hand, thermal conduction and convection between the membrane and the feed and permeate streams are considered. The convection heat transfer coefficient is calculated by using the Gnielinski correlation for the Nusselt number in fully developed turbulent and transition flow conditions [23]: where, 0.5 ≤ Pr ≤ 2000 3000 ≤ Re D h ≤ 10 6 .
Laminar conditions are simulated using the heat and mass transport analogy.

Permeation Model
Each control volume is characterized by an O 2 partial pressure profile and permeation process, as the one depicted in Figure 2. The main steps of the permeation process and a qualitative profile of the O 2 partial pressure gradient are highlighted in the z-axis direction in a supported oxygen-ion conducting membrane. The permeation mechanism is sustained by means of an oxygen partial pressure gradient between the feed and permeate compartments. Each control volume is characterized by an O2 partial pressure profile and permeation process, as the one depicted in Figure 2. The main steps of the permeation process and a qualitative profile of the O2 partial pressure gradient are highlighted in the z-axis direction in a supported oxygen-ion conducting membrane. The permeation mechanism is sustained by means of an oxygen partial pressure gradient between the feed and permeate compartments. With reference to the left sketch of Figure 2, four main oxygen transport steps are distinguished in a supported planar membrane: 1. Gas phase diffusion of molecular oxygen from the bulk phase to the membrane surface on the feed side. 2. Diffusion across the ceramic membrane layer; this comprises of the following sub-processes: i.
Oxygen reduction (ORR), which involves: (i) surface exchange reaction on the feed sidechemisorption of the oxygen molecule onto a metal oxide at an oxide-ion vacancy; (ii) dissociation and reduction of O2 with the consequent production of oxygen anions; (iii) incorporation of O 2-ion into the lattice structure [24]. ii.
Ionic transport throughout the dense membrane material-single-phase (MIEC) or composite. Simultaneous bulk-diffusion of charged species and small polarons in the bulk phase. iii.
Surface-exchange reaction on the permeate interface-oxygen anions recombination to form molecular oxygen and oxygen vacancies in the solid structure 3. Oxygen transport in the porous support structure. 4. Gas diffusion and counter-diffusion of chemical species from the gas bulk to the membrane surface and vice-versa.

Gas Phase Bulk Diffusion
Consider the absolute flux of species in a mixture along the direction normal to the membrane wall for diffusion in a non-stationary medium: With reference to the left sketch of Figure 2, four main oxygen transport steps are distinguished in a supported planar membrane:

1.
Gas phase diffusion of molecular oxygen from the bulk phase to the membrane surface on the feed side.

2.
Diffusion across the ceramic membrane layer; this comprises of the following sub-processes: i. Oxygen reduction (ORR), which involves: (i) surface exchange reaction on the feed side-chemisorption of the oxygen molecule onto a metal oxide at an oxide-ion vacancy; (ii) dissociation and reduction of O 2 with the consequent production of oxygen anions; (iii) incorporation of O 2− ion into the lattice structure [24]. ii.
Ionic transport throughout the dense membrane material-single-phase (MIEC) or composite. Simultaneous bulk-diffusion of charged species and small polarons in the bulk phase. iii.
Surface-exchange reaction on the permeate interface-oxygen anions recombination to form molecular oxygen and oxygen vacancies in the solid structure

3.
Oxygen transport in the porous support structure.

4.
Gas diffusion and counter-diffusion of chemical species from the gas bulk to the membrane surface and vice-versa.

Gas Phase Bulk Diffusion
Consider the absolute flux of species i in a mixture along the direction normal to the membrane wall for diffusion in a non-stationary medium: where, D ik is the mutual diffusion coefficient of species i into k and x i is the molar fraction of species i. For the case of an OTM, i = O 2 and k = mix in air. Considering the conditions on the feed side, the system is characterized by the conservation of the number of moles throughout the whole gas phase, due to the absence of any chemical reaction; therefore, the following holds: Moreover, the concentration profile for the oxygen features a decreasing slope (∂C O 2 /∂z > 0) from the bulk to the membrane wall; on the other hand, the flux of all the other species in the feed mixture is zero, since they are not subjected to any driving force (N mix = 0). Following these considerations, Equation (8) becomes [23]: which can be integrated between the bulk phase (z = δ m ) and the membrane wall (z = 0) through the diffusion layer thickness (δ m ), which is taken equal to the hydraulic diameter of the membrane channels [25]. We assume perfect mixing outside the diffusion boundary layer thickness that is, The integrated form results: where, b and w represent the conditions of points E, A and D, B respectively on the feed and permeate side. Moreover, the ratio D O 2 ,mix /D h is defined as the mass transfer coefficient (h SFM ) in the stagnant film model [26][27][28], which holds if the velocity of the flow normal to the membrane wall is zero across the whole bulk phase. Within the hypotheses of the stagnant film model, the mass transfer coefficient (h SFM ) is rigorously equal to the overall mass transfer coefficient defined as: where, the mutual diffusion coefficients can be estimated with the Fuller's method, as described in the kinetic theory [29]. where, and the molecule diffusion volumes (v i ) can be found in literature [29].
From the definition of the Sherwood number Sh = h m D h /D O 2 ,mix , Equation (10) applied to the O 2 flux becomes: where k represents the feed or permeate side; the binary diffusion coefficient per each component in the mixture is calculated as follows [30]: Moreover, x k O 2 ,b and x k O 2 ,w represent the oxygen molar fraction in the gas bulk phase and at the membrane wall, respectively. For one-dimensional motion in a straight channel in laminar regime, the Sherwood number is solely dependent on the geometry of the channel. Assuming the heat transfer-mass transfer analogy, it is possible to calculate the Sherwood number with the same correlations as for the Nusselt number in rectangular cross-sectional area channels [23]. For channel width-to-height ratios (b/a) between 1 and 4, the following interpolation function may be employed.
For higher values of b/a ratios, a piecewise linear interpolation of the values in Reference [23] is used. A similar relation can be written for the permeate side.

Solid Phase Bulk Diffusion
As far as the oxygen diffusion mechanism in ionic conductors is concerned, the flux density of the charge carriers is expressed by: where, σ k is the conductivity of the charge carrier k (i.e., electrons or oxygen anions), z k is the number of electrical charges of the carrier k (z O 2− = 2) and µ and φ are the chemical and electrostatic potentials across the membrane, respectively. MIECs materials are known to have electronic conductivity σ e − orders of magnitude greater than the ionic one σ O 2− . The coupling of Ohm's law and Nernst-Einstein equation is used to calculate the electrical conductivity of a material with respect to a defined charge carrier; the first one allows defining the current density as a function of the electrostatic potential gradient in an isotropic medium: where, c k is the concentration of mobile charge carriers (i.e., electrons, holes) and v k is the velocity of the charge carrier moving in the lattice structure. The latter shall be replaced with the charge carrier drift mobility defined as u k ≡ v k /∇φ. The Nernst-Einstein equation is used to connect the drift mobility to the electrical conductivity throughout a conductivity-related diffusivity factor D k . hence, where, D k can be expressed in an Arrhenius-like formulation D k = D • k exp(−∆H k /RT). Consequently, considering the definition of the chemical potential for an ideal component, µ k = µ • k + RT ln C k and assuming that the gradient of the electrostatic potential is negligible in comparison to the chemical potential due to the higher electronic conductivity of electrons with respect to oxygen anions, the flux of the charge carrier k can finally be expressed in the form of a Fick's first law: Membranes 2020, 10, 442 9 of 32 Considering the main charge carrier as the oxygen vacancies in the solid lattice (k = V •• O ) and knowing that j V •• O = −2j O 2 , the oxygen flux through the bulk of the solid perovskite lattice structure (step 3) is found and expressed in the preferential direction normal to the membrane wall as:

Superficial Reactions
The overall permeation process and the concentration of the charge carriers on each side of the membrane is also governed by the surface reactions (Equation (23) and Equation (24)), which correspond to steps 2i and 2iii, respectively.
where, V •• O and h • are the charge carriers for the ionic and electronic conduction process respectively and O X O are the oxygen atoms embedded in the lattice structure of the perovskite membrane. Moreover, k f and k r are the forward and reverse reaction rate constants and can be expressed with an Arrhenius-like formulation: Considering the reported reactions as elementary steps in the overall mechanism and applying the law of mass-actions, which states that the np product has to be constant at a given temperature, the charge carriers' concentrations can be related to the oxygen partial pressure, as follows: where, k eq = k f /k r is the equilibrium constant of the surface reaction. Finally, the reaction rate of the two reactions rates, which also correspond to the oxygen flux through the membrane, shall be expressed as: where, the apices and indicate respectively the feed (reduction, surface B in Figure 2) and the permeate (oxidation, surface C) side. The final assumption for the surface reaction model follows the hypothesis of Xu and Thomson 1999 [31], which states that the two surface reactions could be considered as pseudo-zero order, as far as the concentration of oxygen ions and electron vacancies; therefore, the two equations-Equation (29) and Equation (30) [32], does not contribute to a great variation in the oxygen flux, which is only slightly overestimated in the zero-order equation.
The combination of Equation (22), Equation (29) and Equation (30) allows obtaining a simplified equation for the oxygen flux through a perovskite membrane dependent only upon the bulk diffusion, the surface reaction process and the thickness of the membrane itself.
However, the oxygen transport across membranes characterized by a thin dense layer and operating at high temperature are found to be mostly governed by kinetic processes, whilst the ionic transport resistance within the lattice structure results to be negligible. This assumption was already verified by Chen et al., 1997 [33] and reported in Plazaola et al., 2019 [13], whereby the oxygen permeation becomes partially limited by surface oxygen exchange for thicknesses below 1 mm. The rate limiting step is constituted by the recombination of oxygen ions (O 2− ) that permeate through the membrane lattice to form molecular oxygen on the permeate side surface.
In this scenario, the first term of the denominator in Equation (31) can be neglected and the oxygen flux can be calculated by the following equation: where, p O 2 ,B and p O 2, C are the oxygen partial pressure on the opposite surfaces of the membrane while k r is the kinetic constant of the recombination reaction.

Permeation in Porous Media
The model considers the possibility of simulating supported membranes (step #3 in Figure 2). One of the most employed models for the simulation of mass transport in a porous media is the Dusty Gas Model (DGM), which consists of an application of the Maxwell-Stefan diffusion equation [30]. The DGM considers the parallel contribution of three different diffusion mechanism: (i) bulk diffusion (molecule-molecule collision dominated); (ii) Knudsen diffusion (molecule-wall collision dominated); (iii) surface diffusion/micro-pore diffusion (diffusion of adsorbed species along the wall). The model implicit formulation can be written as: where, c i is its concentration and µ i,T is the chemical potential at a fixed temperature; moreover, D e f f ij and D e f f Kn,i are the effective binary and Knudsen diffusion coefficients corrected by the ratio of the porosity and tortuosity ( /τ) of the porous structure. The binary diffusion coefficient is estimated through the kinetic theory [29] as already reported in Equation (13). The Knudsen coefficient is equal to: where, r pore is the average radius of the pores in the support. The chemical potential can be expressed as a function of the pressure gradient as: Membranes 2020, 10, 442 11 of 32 Finally, introducing the overall flux of the species i expressed by Equation (8) as, where, the velocity in the advective term can be modelled by means of the Darcy's law for viscous flows [34]: Considering that the support bulk structure is not characterized by any chemical or electrochemical reactions processes, the number of moles remains the same throughout the whole of its thickness. This hypothesis supports using the Graham's law which states the absence of a pressure gradient as a consequence of reaction stoichiometry [30].
Finally, considering the methodology reported in Zhu et al. 2005 [35] to develop a direct equation for the DGM, the following simplified oxygen flux formulation is obtained (i = O 2 ).
Therefore, at fixed partial pressures on the feed and permeate side (i.e., p A and p E ), the model employs Equation (14) (for both feed and sweep sides), Equation (32) and Equation (39) and solves the nonlinear equation system employing a modified Powell hybrid algorithm and finite-difference approximation to the Jacobian to obtain j O 2 , p B , p C and p D .

Experimental Set-Up
The permeation of asymmetric La 0.6 Sr 0.4 Co 0.2 Fe 0.8 O 3-δ (LSCF) membranes and bare porous support are measured at RSE laboratories in the frame of the FP7 European Project GREEN-CC. Disk-shaped samples are manufactured and provided by Forschungszentrum Jülich GmbH (Jülich, Germany). Asymmetric membranes consist of a 30 µm thick dense LSCF layer, supported on a 0.7 mm porous LSCF structure. No catalytic active layer was deposited on the tested samples. The membrane has a circular shape with a diameter of 13.7 mm and its active area is 147.4 mm 2 . The support porosity is estimated by a quantitative image analysis by using a Nikon LUCIA software and it is approximately 0.45.
A SEM (Scanning Electron Microscopy) image of the polished cross section of an asymmetric membrane is shown in Figure 3.
LSCF bare porous supports, with a thickness of 1 and 2 mm are tested at room temperature, in pure Helium, by feeding a known feed gas flow rate to the test section in the absence of sweep gas and by measuring the pressure difference across the sample, as depicted in Figure 4. LSCF structure. No catalytic active layer was deposited on the tested samples. The membrane has a circular shape with a diameter of 13.7 mm and its active area is 147.4 mm 2 . The support porosity is estimated by a quantitative image analysis by using a Nikon LUCIA software and it is approximately 0.45.
A SEM (Scanning Electron Microscopy) image of the polished cross section of an asymmetric membrane is shown in Figure 3. LSCF bare porous supports, with a thickness of 1 and 2 mm are tested at room temperature, in pure Helium, by feeding a known feed gas flow rate to the test section in the absence of sweep gas and by measuring the pressure difference across the sample, as depicted in Figure 4. High temperature permeation tests are performed on asymmetric membranes in a glass permeation cell, similar to the one described in Reference [36]. The permeation cell is placed in an electric furnace to heat the system to the operating temperature, monitored by a thermocouple in the sweep gas tube (Figure 5a). The gas tightness between the membrane and the permeation cell is ensured at high temperature (approx. 1050 °C) by two pure gold O-rings, on both side of the sample. Tests are carried out in the 700-950 °C temperature range, by feeding air or pure oxygen as process gas on the membrane side and helium or nitrogen as sweep stream on the porous support side. The A SEM (Scanning Electron Microscopy) image of the polished cross section of an asymmetric membrane is shown in Figure 3. LSCF bare porous supports, with a thickness of 1 and 2 mm are tested at room temperature, in pure Helium, by feeding a known feed gas flow rate to the test section in the absence of sweep gas and by measuring the pressure difference across the sample, as depicted in Figure 4. High temperature permeation tests are performed on asymmetric membranes in a glass permeation cell, similar to the one described in Reference [36]. The permeation cell is placed in an electric furnace to heat the system to the operating temperature, monitored by a thermocouple in the sweep gas tube (Figure 5a). The gas tightness between the membrane and the permeation cell is ensured at high temperature (approx. 1050 °C) by two pure gold O-rings, on both side of the sample. Tests are carried out in the 700-950 °C temperature range, by feeding air or pure oxygen as process gas on the membrane side and helium or nitrogen as sweep stream on the porous support side. The High temperature permeation tests are performed on asymmetric membranes in a glass permeation cell, similar to the one described in Reference [36]. The permeation cell is placed in an electric furnace to heat the system to the operating temperature, monitored by a thermocouple in the sweep gas tube (Figure 5a). The gas tightness between the membrane and the permeation cell is ensured at high temperature (approx. 1050 • C) by two pure gold O-rings, on both side of the sample. Tests are carried out in the 700-950 • C temperature range, by feeding air or pure oxygen as process gas on the membrane side and helium or nitrogen as sweep stream on the porous support side. The test configuration is reported in Figure 5b. Both the process gas and the sweep gas are fed at atmospheric pressure and are mass flow controlled-the process gas flow rate is kept constant at 250 NmL min −1 and 200 NmL min −1 , in the case of air and pure oxygen, respectively, while the sweep gas flow rate is varied between 100 and 300 NmL min −1 .
The oxygen permeation is estimated by measuring the oxygen concentration in the permeate stream by an online micro-Gas Chromatographer (m-GC), while the absence of membrane leakage is ensured by continuously monitoring the nitrogen concentration in the permeate side. The experimental error on the O 2 flux measure results to be at maximum 4%.

Results of Experimental Tests
The helium permeance of six different bare supports, 1 or 2 mm thick, is determined at room temperature as a function of the average pressure measured across the sample, to verify the homogeneity of the samples' porosity. The comparison between results obtained on samples with different thickness, is reported in Figure 6 as a function of the average pressure [37].
Membranes 2020, 10, 12 of 32 test configuration is reported in Figure 5b. Both the process gas and the sweep gas are fed at atmospheric pressure and are mass flow controlled-the process gas flow rate is kept constant at 250 NmL min −1 and 200 NmL min −1 , in the case of air and pure oxygen, respectively, while the sweep gas flow rate is varied between 100 and 300 NmL min −1 . The oxygen permeation is estimated by measuring the oxygen concentration in the permeate stream by an online micro-Gas Chromatographer (m-GC), while the absence of membrane leakage is ensured by continuously monitoring the nitrogen concentration in the permeate side. The experimental error on the O2 flux measure results to be at maximum 4%.

Results of Experimental Tests
The helium permeance of six different bare supports, 1 or 2 mm thick, is determined at room temperature as a function of the average pressure measured across the sample, to verify the homogeneity of the samples' porosity. The comparison between results obtained on samples with different thickness, is reported in Figure 6 as a function of the average pressure [37]. As expected, the He permeance measured on supports with a thickness of 2 mm is about two times lower than the permeance measured on thinner supports (1 mm), therefore the reduction of the He permeance is proportional with the increase of the support thickness.  The oxygen permeation is estimated by measuring the oxygen concentration in the permeate stream by an online micro-Gas Chromatographer (m-GC), while the absence of membrane leakage is ensured by continuously monitoring the nitrogen concentration in the permeate side. The experimental error on the O2 flux measure results to be at maximum 4%.

Results of Experimental Tests
The helium permeance of six different bare supports, 1 or 2 mm thick, is determined at room temperature as a function of the average pressure measured across the sample, to verify the homogeneity of the samples' porosity. The comparison between results obtained on samples with different thickness, is reported in Figure 6 as a function of the average pressure [37]. As expected, the He permeance measured on supports with a thickness of 2 mm is about two times lower than the permeance measured on thinner supports (1 mm), therefore the reduction of the He permeance is proportional with the increase of the support thickness. As expected, the He permeance measured on supports with a thickness of 2 mm is about two times lower than the permeance measured on thinner supports (1 mm), therefore the reduction of the He permeance is proportional with the increase of the support thickness.
Oxygen permeation tests are performed at high temperature on an asymmetric membrane by feeding synthetic air (with a constant flow rate of N 2 = 197.5 NmL min −1 , O 2 = 52.5 NmL min −1 ) or pure oxygen as process gas (with a constant flow rate of 200 NmL min −1 ) and helium or nitrogen as sweep gas, by varying the flow rate in the range 100-300 NmL min −1 . Results are shown in Figure 7, in terms of oxygen flux as a function of the sweep gas flow rate for different temperature [38].
As predicted by Equation (31), the oxygen flux increases with the temperature. By feeding pure oxygen instead of air as process gas, an improvement of the oxygen permeation is observed, due to the increase of the oxygen partial pressure in the process side p O 2 and, consequently, of the process driving force.
The membrane permeation can be also described by an Arrhenius-type correlation. As shown in Figure 8, in the range 950-750 • C a different Arrhenius behavior can be observed, which suggests a variation in the membrane permeation mechanism as a function of the temperature. At high temperature, the bulk diffusion through the membrane dense layer mainly controls the permeation, while at lower temperature the oxygen surface exchange is predominant [12,39]. The estimated activation energy values are in accordance with the results reported by Serra et al. 2013 [12].
the increase of the oxygen partial pressure in the process side ( ) 2 O p ′ and, consequently, of the process driving force.
The membrane permeation can be also described by an Arrhenius-type correlation. As shown in Figure 8, in the range 950-750 °C a different Arrhenius behavior can be observed, which suggests a variation in the membrane permeation mechanism as a function of the temperature. At high temperature, the bulk diffusion through the membrane dense layer mainly controls the permeation, while at lower temperature the oxygen surface exchange is predominant [12,39]. The estimated activation energy values are in accordance with the results reported by Serra et al. 2013 [12].  The oxygen flux through the membrane increases with rising temperature, both in the presence of helium and nitrogen. At high temperature (850 °C-900 °C), oxygen flux is slightly influenced by the increase of sweep gas flow rate due to the higher oxygen dilution in the permeate and by the increase of the process driving force. On the other hand, at lower temperatures permeation is almost The oxygen flux through the membrane increases with rising temperature, both in the presence of helium and nitrogen. At high temperature (850 • C-900 • C), oxygen flux is slightly influenced by the increase of sweep gas flow rate due to the higher oxygen dilution in the permeate and by the increase of the process driving force. On the other hand, at lower temperatures permeation is almost constant in the investigated sweep gas flow rate range. Moreover, the oxygen flux is slightly higher when helium is used as sweep gas, instead of nitrogen.
In summary, the experimental data are available at 950 • C, 900 • C and 850 • C and atmospheric pressure. The feed and permeate streams are changed in the following conditions: I air feed (200 NL min −1 ) and helium as sweep gas (100-300 NL min −1 ), I oxygen feed (200 NL min −1 ) and helium as sweep gas (100-300 NL min −1 ), I oxygen feed (200 NL min −1 ) and nitrogen as sweep gas (100-300 NL min −1 ).

Support Calibration
Given the geometric characteristics of the model presented, it cannot be directly applied to evaluate the permeation of the circular membrane as the one used for the tests. Introducing the actual axisymmetric flow field occurring in the real disk-shaped membrane setup in the presented model requires the application of CFD simulations that is beyond the scope of this study. In order calibrate the main parameters of the model against the experimental test, the permeation model (Equation (1) to Equation (32)) is applied to evaluate the resulting O 2 flux considering the thermodynamic and chemical variables averaged between inlet and outlet on each side of the circular membrane in a single discretization step.
To calibrate the unknown parameters of the support (r pore and τ), the set of permeance data for helium are considered (see Figure 5). At the highest end of the measured flow rates range, the flow regime results turbulent; therefore, Equation (39) is not able to completely describe the mass transfer phenomena occurring in the porous support, because it is developed for laminar diffusion regime. Hence, a mixed laminar-turbulent diffusion semi-empirical model is employed [40].
where, The parameters α and β are calibrated against the support permeance experimental data minimizing the sum of the root mean square of the relative difference between the model and the measured data. The calibrated values in the considered operating range are α = 4.27 × 10 −8 m 2 , β = 3.41 × 10 −9 m. The average error results to be 3.14% with a maximum error comprised between −6.24% and 5.33%. The result of the calibration is represented in Figure 9.
Subsequently, the pressure drop obtained from the semi-empirical model are employed to calibrate the physical model in Equation (39) for the DGM in the laminar flow regime (1-12 NmL min −1 -cm 2 ). Being the experimental data relative to the permeation of a pure species (i.e., He) across the porous support, Equation (39) reduces to a simplified form which depends solely upon the Knudsen diffusivity D Kn and not on the molecular diffusivity in a binary mixture D He,j . Consequently, knowing the geometry of the support-thickness 1.0 and 2.0 mm-and its porosity 0.45, the only unknown characteristics is the ratio between the pore radius and the tortuosity, that is, r pore /τ. The DGM is calibrated in the flow regime of interest for laminar mass transport diffusion in a porous media by finding the r pore /τ ratio which minimises the sum of the root mean square of the relative errors between the model and the measured data. The resulting pore radius-to-tortuosity ratio r pore /τ is 0.581 µm. As far as the value of tortuosity of the porous media, a value equal to τ = 1.765 is taken, as results from applying the , as reported in Kong et al., 2015 [41]; the resulting value of tortuosity is also in line with the simulations performed via microstructure reconstruction by Joos et al., 2011 [42]. The resulting r pore is therefore 1.03 µm. where, • J [m 3 s −1 m −2 ] is the volumetric flux across the porous structure. The parameters α and β are calibrated against the support permeance experimental data minimizing the sum of the root mean square of the relative difference between the model and the measured data. The calibrated values in the considered operating range are α = 4.27 × 10 −8 m 2 , β = 3.41 × 10 −9 m. The average error results to be 3.14% with a maximum error comprised between −6.24% and 5.33%. The result of the calibration is represented in Figure 9. Subsequently, the pressure drop obtained from the semi-empirical model are employed to calibrate the physical model in Equation (39) for the DGM in the laminar flow regime (1-12 NmL min −1 -cm 2 ). Being the experimental data relative to the permeation of a pure species (i.e., He) across

Supported Membrane Calibration
To calibrate the supported membrane model, the kinetic coefficients (k 0 r , E 0 r ) of the Arrhenius-like equation for the kinetic constant in Equation (32) are calibrated by minimizing the mean squared deviation between the experimental and calculated oxygen flux in a temperature range between 850-950 • C. The minimization procedure employs the genetic algorithm solver available with Microsoft Excel. The calibration is performed by considering a planar membrane of the same active area (147.4 mm 2 ) and porous support thickness (i.e., 0.7 mm) of the experiments.
A rectangular geometry with co-current flows is considered. Length is assumed equal to the radius of the circular membrane (6.85 mm) tested. The resulting channel width is 21.5 mm. Since the width-to-height ratio (b/a) for both sweep and feed channels is impossible to define in the experimental set-up, this parameter is treated as the third calibration parameter, constrained in a range 1-∞. Moreover, the assumption of negligible pressure losses is employed for the calibration case.
In summary, the calibration of the supported membrane consists of a constrained optimization with the following objective function, variables, and limits:

Box 1
Objective Function: The calibrated parameters for the supported membrane result: The average root mean square of the error on the experimental data available is 3.99%, thus considered acceptable throughout the whole operating range. Figure 10 reports the results of the calibration procedure in terms of comparison between modelling and experimental data ( Figure 10a) and prediction error (Figure 10b). The values of permeation activation energies support the claim that the main process is governed by oxygen recombination reaction on the permeate side.

Heat and Mass Balance of Glass Melting Furnace
In this work, we apply the OTM technology to a glass melting furnace. The work by Sardeshpande et al. [5] proposed a model to achieve the minimum specific energy consumption in glass furnaces. Based on field analyses and experimentations in the glass industry, the paper provides the design and operating parameters of a case study representing the current state of the art of glass furnaces. The furnace's capacity and draw are fixed at 100 ton day −1 and 90 ton day −1 , respectively, whereas the cullet fraction is set at 40%. Comparing with Beerkens [6], it is possible to appreciate that it is a very efficient benchmark, even considering the actual cullet (40% vs. 50%), as the energy input amounts to 3833 kJ kg −1 .
Starting from the work by Sardeshpande et al. [5], heat and mass balances are estimated with GS, a computer code in-house developed in the past years at the Department of Energy of Politecnico di Milano [43]. The code is a powerful and flexible tool that can accurately predict the performance of a wide variety of chemical processes and power generation systems. GS was originally designed to assess the performance of gas-steam cycles for power production [44] and has been progressively developed and improved to calculate complex systems including coal gasification, chemical reactors, fuel cells and essentially all the processes present in advanced plants for power generation [45][46][47][48]. When addressing the new oxy-fuel furnace solutions, the detailed membrane model, as previously

Heat and Mass Balance of Glass Melting Furnace
In this work, we apply the OTM technology to a glass melting furnace. The work by Sardeshpande et al. [5] proposed a model to achieve the minimum specific energy consumption in glass furnaces. Based on field analyses and experimentations in the glass industry, the paper provides the design and operating parameters of a case study representing the current state of the art of glass furnaces. The furnace's capacity and draw are fixed at 100 ton day −1 and 90 ton day −1 , respectively, whereas the cullet fraction is set at 40%. Comparing with Beerkens [6], it is possible to appreciate that it is a very efficient benchmark, even considering the actual cullet (40% vs. 50%), as the energy input amounts to 3833 kJ kg −1 .
Starting from the work by Sardeshpande et al. [5], heat and mass balances are estimated with GS, a computer code in-house developed in the past years at the Department of Energy of Politecnico di Milano [43]. The code is a powerful and flexible tool that can accurately predict the performance of a wide variety of chemical processes and power generation systems. GS was originally designed to assess the performance of gas-steam cycles for power production [44] and has been progressively developed and improved to calculate complex systems including coal gasification, chemical reactors, fuel cells and essentially all the processes present in advanced plants for power generation [45][46][47][48]. When addressing the new oxy-fuel furnace solutions, the detailed membrane model, as previously described, is employed to propose a reactor design which could fit the process operating conditions derived from GS.
This section presents four cases of glass melting furnaces: Case A: air-blown reference furnace; Case B: oxy-fuel furnace with heat recovery; Case C: oxy-fuel furnace with integrated 3-end OTM module; Case D: oxy-fuel furnace with integrated 4-end OTM module. The results of mass and energy balances in terms of mass flow rates, temperatures and compositions are reported in Appendix A, along with the calculation assumptions for the fluid machinery. The fuel input to the furnace is always natural gas with a lower heating value (LHV) of 46.48 MJ kg −1 (see Appendix A for the composition).

Air-Blown Reference Case
The first case considered is the benchmark, according to the work by Sardeshpande et al. [5]. Figure 11 shows a simplified lay-out of the air-blown system and Table A1 in Appendix A reports the main streams' thermochemical conditions as resulting from the mass and energy balances.
The combustion air flow rate is set to obtain an O 2 content of 2% for the in the gas exiting the furnace [5].
along with the calculation assumptions for the fluid machinery. The fuel input to the furnace is always natural gas with a lower heating value (LHV) of 46.48 MJ kg −1 (see Appendix A for the composition).

Air-Blown Reference Case
The first case considered is the benchmark, according to the work by Sardeshpande et al. [5]. Figure 11 shows a simplified lay-out of the air-blown system and Table A1 in Appendix A reports the main streams' thermochemical conditions as resulting from the mass and energy balances. The combustion air flow rate is set to obtain an O2 content of 2% for the in the gas exiting the furnace [5]. Figure 11. Schematic of the air-blown reference furnace.
A larger mass flow rate can be observed at station #7 compared to station #6 because of air infiltration, which is calculated according to an O2 content at the regenerator output of 8% [5]. In the case of no air infiltration, a temperature of 610.3 °C would be calculated instead of 361.8 °C, considering the heat losses from the regenerator. The after-treatment station cools down the gas at the furnace exit to 150 °C by diluting with ambient air. The resulting stream is de-dusted and scrubbed before being blown to the chimney (stream #8 is significantly larger than stream #7).

Oxyfuel Case with VSA Unit
Compared to a similar air-blown case, Beerkens and Muysenberg [49] reported that an efficient oxy-fuel glass melting furnace requires ~16.7% less fuel; however, such saving does not take into account the corresponding energy consumption required for oxygen production. Therefore, with reference to the calculations in Case A, the input fuel flow rate is tuned to maintain constant the thermal heat loss of the furnace for a fixed production capacity (90 ton day −1 ) [49]. Without varying the temperature of the gas leaving the furnace, the heat loss is reduced by 23.5%, compared to Case A.
The overall energy consumption for an oxy-fuel glass melting furnace is approximately equal to the one in case of air combustion [7]. Beerkens [7] assumed 0.375-0.4 kWh for the electricity consumption per Nm 3 of pure oxygen generation. Based on an electric plant efficiency of 40%, ~3400 kJ Nm −3 are assumed for the primary energy consumption. The new oxy-fuel case results in almost equal total fuel consumption (exactly 1.1% higher) compared to Case A, fully consistent with the technical literature [7].
Non-regenerative furnaces may be advantageous from an energy viewpoint compared to conventional regenerative furnaces because the exit temperatures are higher and allow for a better heat recovery section. Inspired by the HotOxyGlass project [50], Figure 12 shows the schematic of an A larger mass flow rate can be observed at station #7 compared to station #6 because of air infiltration, which is calculated according to an O 2 content at the regenerator output of 8% [5]. In the case of no air infiltration, a temperature of 610.3 • C would be calculated instead of 361.8 • C, considering the heat losses from the regenerator. The after-treatment station cools down the gas at the furnace exit to 150 • C by diluting with ambient air. The resulting stream is de-dusted and scrubbed before being blown to the chimney (stream #8 is significantly larger than stream #7).

Oxyfuel Case with VSA Unit
Compared to a similar air-blown case, Beerkens and Muysenberg [49] reported that an efficient oxy-fuel glass melting furnace requires~16.7% less fuel; however, such saving does not take into account the corresponding energy consumption required for oxygen production. Therefore, with reference to the calculations in Case A, the input fuel flow rate is tuned to maintain constant the thermal heat loss of the furnace for a fixed production capacity (90 ton day −1 ) [49]. Without varying the temperature of the gas leaving the furnace, the heat loss is reduced by 23.5%, compared to Case A.
The overall energy consumption for an oxy-fuel glass melting furnace is approximately equal to the one in case of air combustion [7]. Beerkens [7] assumed 0.375-0.4 kWh for the electricity consumption per Nm 3 of pure oxygen generation. Based on an electric plant efficiency of 40%,~3400 kJ Nm −3 are assumed for the primary energy consumption. The new oxy-fuel case results in almost equal total fuel consumption (exactly 1.1% higher) compared to Case A, fully consistent with the technical literature [7].
Non-regenerative furnaces may be advantageous from an energy viewpoint compared to conventional regenerative furnaces because the exit temperatures are higher and allow for a better heat recovery section. Inspired by the HotOxyGlass project [50], Figure 12 shows the schematic of an oxy-fuel glass furnace where air (#11) is pre-heated in a recuperator by the high-temperature gas exiting the furnace (#6). This hot air stream (#12) is then directed to other heat exchangers for pre-heating both the oxidizer resulting from a VSA unit and the natural gas. In this new case a further reduction in natural gas consumption is achieved, with an overall fuel saving of 23.8% compared to Case A. Although such a fuel saving does not include the primary energy related to the oxygen production, this result is promising since the HotOxyGlass project [50] declared a saving of 25% ± 2%. After including the necessary energy contribution for oxygen production, the final saving will be lower.
As reported in Table A2 and differently from the previous case, streams #6 and #7 are equal in flow rate. However, a larger air dilution in the after-treatment station is required in this case since the temperature of stream #7 is higher compared Case A.
2%. After including the necessary energy contribution for oxygen production, the final saving will be lower.
As reported in Table A2 and differently from the previous case, streams #6 and #7 are equal in flow rate. However, a larger air dilution in the after-treatment station is required in this case since the temperature of stream #7 is higher compared Case A.

Oxyfuel Cases with Integrated OTM Systems
Two innovative cases are investigated in this paper, both considering a membrane-based system for the separation of oxygen from a hot air stream. The layouts are designed to self-produce the oxygen required for the combustion in the glass melting furnace and the electricity for all auxiliary loads. The proposed cases consider 3-and 4-end membrane systems, without and with sweep stream fed to permeate side.
A first layout for the advanced glass melting furnace is schematized in Figure 13 (Case C). Details of the main streams are reported in Table A3 in the Appendix A. Air (#7) is delivered by a compressor and pre-heated up to 600 °C through a recuperative heat exchanger by the hot gas from the turbine T1 (#11) and the high-temperature oxygen from the membrane. Subsequently, air is heated up to ~700°C in another heat exchanger by the high-temperature gas (#4) at the furnace exit. After the combustion with natural gas to rise the temperature up to 950 °C, the stream (#9) enters the membrane system at 4 bar at high (>18 mol.%) O2 content. The OTM separates 50 wt.% of the oxygen in the feed stream on the permeate side at 0.25 bar, which is then cooled down to 40 °C (#12). Ultimately, the oxidizer is delivered to the furnace by a vacuum pump (C1) through a recuperative heat exchanger, where it is pre-heated (up to 550 °C) together with the natural gas required at the burners.

Oxyfuel Cases with Integrated OTM Systems
Two innovative cases are investigated in this paper, both considering a membrane-based system for the separation of oxygen from a hot air stream. The layouts are designed to self-produce the oxygen required for the combustion in the glass melting furnace and the electricity for all auxiliary loads. The proposed cases consider 3-and 4-end membrane systems, without and with sweep stream fed to permeate side.
A first layout for the advanced glass melting furnace is schematized in Figure 13 (Case C). Details of the main streams are reported in Table A3 in the Appendix A. Air (#7) is delivered by a compressor and pre-heated up to 600 • C through a recuperative heat exchanger by the hot gas from the turbine T1 (#11) and the high-temperature oxygen from the membrane. Subsequently, air is heated up to~700 • C in another heat exchanger by the high-temperature gas (#4) at the furnace exit. After the combustion with natural gas to rise the temperature up to 950 • C, the stream (#9) enters the membrane system at 4 bar at high (>18 mol.%) O 2 content. The OTM separates 50 wt.% of the oxygen in the feed stream on the permeate side at 0.25 bar, which is then cooled down to 40 • C (#12). Ultimately, the oxidizer is delivered to the furnace by a vacuum pump (C1) through a recuperative heat exchanger, where it is pre-heated (up to 550 • C) together with the natural gas required at the burners. The results reported in Table A3 are calculated with the hypothesis of no air infiltration in the sub-atmospheric pressure zones of the plant. Air infiltration in permeate side of the membrane would affect the temperature profiles of the multi-flow heat exchanger downstream turbine T1, as well as the temperatures at stations #5 and #6. However, a minimum ΔTapproach ≈ 50 °C is always verified from the calculations, hence guaranteeing an ample operative margin in case some air infiltrations occur. The difference in flow rates #6 and #13 is justified by air dilution in the after-treatment station, as in the previous cases.
It is found that the new oxy-combustion furnace consumes 291.6 kg h −1 of fuel. This amount is greater than the fuel demand of Case B due to the fuel required to rise the temperature of the feed stream to the membrane-based system. Nevertheless, the fuel demand is reduced by ~6% compared to Case A.
Another layout for the advanced glass melting furnace integrated with the OTM technology is The results reported in Table A3 are calculated with the hypothesis of no air infiltration in the sub-atmospheric pressure zones of the plant. Air infiltration in permeate side of the membrane would affect the temperature profiles of the multi-flow heat exchanger downstream turbine T1, as well as the temperatures at stations #5 and #6. However, a minimum ∆T approach ≈ 50 • C is always verified from the calculations, hence guaranteeing an ample operative margin in case some air infiltrations occur. The difference in flow rates #6 and #13 is justified by air dilution in the after-treatment station, as in the previous cases.
It is found that the new oxy-combustion furnace consumes 291.6 kg h −1 of fuel. This amount is greater than the fuel demand of Case B due to the fuel required to rise the temperature of the feed stream to the membrane-based system. Nevertheless, the fuel demand is reduced by~6% compared to Case A.
Another layout for the advanced glass melting furnace integrated with the OTM technology is proposed in Figure 14 (Case D). Here we use of a 4-end membrane system with sweep stream fed to the permeate compartment. Differently from Case C, a sweep stream is provided to the membrane by partially recycling the gas exiting the after-treatment station. It is worth noting that the sweep gas (stream #13 in Figure 14) is characterized by a high CO 2 content. In these conditions, some authors report an increased degradation and performance loss of the LSCF membrane due to segregation of Sr in the form of SrCO 3 [51]. Other experimental studies suggest the possibility of increased oxygen permeation flux for increasing CO 2 concentrations in the sweep gas at operating temperatures of 1000 • C [12]. Within the framework of this study, we neglect any effect that CO 2 might have on the permeation performance of the membrane, confident that the methodology proposed here is equally valid for other types of materials, like NiFe 2 O 4 -Ce 0.8 Tb 0.2 O 2-δ (NFO-CTO) that are less susceptible to CO 2 influence. the temperatures at stations #5 and #6. However, a minimum ΔTapproach ≈ 50 °C is always verified from the calculations, hence guaranteeing an ample operative margin in case some air infiltrations occur. The difference in flow rates #6 and #13 is justified by air dilution in the after-treatment station, as in the previous cases.
It is found that the new oxy-combustion furnace consumes 291.6 kg h −1 of fuel. This amount is greater than the fuel demand of Case B due to the fuel required to rise the temperature of the feed stream to the membrane-based system. Nevertheless, the fuel demand is reduced by ~6% compared to Case A.
Another layout for the advanced glass melting furnace integrated with the OTM technology is proposed in Figure 14 (Case D). Here we use of a 4-end membrane system with sweep stream fed to the permeate compartment. Differently from Case C, a sweep stream is provided to the membrane by partially recycling the gas exiting the after-treatment station. It is worth noting that the sweep gas (stream #13 in Figure 14) is characterized by a high CO2 content. In these conditions, some authors report an increased degradation and performance loss of the LSCF membrane due to segregation of Sr in the form of SrCO3 [51]. Other experimental studies suggest the possibility of increased oxygen permeation flux for increasing CO2 concentrations in the sweep gas at operating temperatures of 1000 °C [12]. Within the framework of this study, we neglect any effect that CO2 might have on the permeation performance of the membrane, confident that the methodology proposed here is equally valid for other types of materials, like NiFe2O4-Ce0.8Tb0.2O2-δ (NFO-CTO) that are less susceptible to CO2 influence.
The 4-end case reproduces a solution that is intermediate between the air-blown case and the one with oxygen only. With reference to Figure 14, Case D results much more complex compared to Case C.  Similarly to Case C, it is necessary to raise the temperature of stream #19 up to 950 • C at the inlet of the OTM system, as seen in Table A4.
A sweep stream (#13) enters the permeate side of the membrane system and exits with a fixed O 2 content (50 mol.%). The gas exiting the glass melting furnace, downstream of the after-treatment station, is partly vented and partly recirculated to the membrane, after mixing with air (#9). Consequently, the O 2 content in stream #11 is sufficient for raising the temperature of the sweep stream (#13) up to 950 • C by catalytic combustion. Based on the results in Table A4, air flow rate #9 is negligible compared to the exhaust gas recirculation (#3).

Energy Assessment
An overview of the cases presented in the previous paragraphs is reported in Table 1, which details the heat fluxes entering and exiting the glass furnace. The heat in the glass is always the same a fixed production of 90 ton day −1 (1628.1 kW), so the fuel reduction for cases B to D is very interesting compared to the reference case. The heat recovery for cases B to D is clearly lower compared to Case A but heat losses are reduced. According to the previous description, Case C, more than Case D, improves Case B from an energy efficiency viewpoint. Heat losses in the flue gases are significantly reduced by transferring heat to the compressed air stream. On the other hand, the heat recovered and recycled into the furnace is lower in Case C compared to Case B, because the oxygen stream delivered by the vacuum pump is already at~250 • C. When considering Case D, additional fuel is introduced as stream #12 in Figure 14, necessary to raise the temperature of the sweep stream up to 950 • C. Moreover, the heat supplied by O 2 separation is higher for Case D compared to Case C, where the oxidant is cooled (stream #12 in Figure 13) prior to be delivered by the vacuum pump. As detailed in Table 1, the heat transferred to the micro-gas turbine cycle in Case D is higher than in Case C. An overview of the energy consumption results of the four cases is reported in Table 2. Details about calculation assumptions for all balance of plant components are reported in Appendix A. The net power from the micro-gas turbine cycle is a system output ad is reported as negative. The overall electric power consumption and production is converted into an equivalent fuel figure, considering an electric conversion of 50% and a LHV of the fuel equal to 46.48 MJ kg −1 . With reference to the primary fuel results, it is possible to realize the superiority of Case C, from an energy saving point of view. Indeed, Case D does not seem to be more attractive than Case B, owing to the complexity of the deeply integrated system.

Design and Performance of the Integrated OTMs
The design of OTM module for the system configuration with a 3-end OTM and 4-end OTM membrane (i.e., cases C and D) are summarized. For the sweep and non-sweep membrane modules, the flow-field configuration is respectively planar counter-flow and planar co-flow. Each channel geometry is characterized by a fixed width (i.e., 21.5 mm) and a fixed height/width ratio on both feed and sweep side, b/a = 1.0 and b/a = 3.01, respectively. Moreover, a mechanical support structure of semi-thickness equal to 2.5 mm on every side of each channel forms the flow-field geometry. This structure does not influence the fluid-dynamic performance of the channels, but it is considered for the calculation of the final module volume. The membrane and support structure and transport parameters are equal to the calibration case. Moreover, each membrane layer is characterized by 20 equal adjacent channels and each stack is made up of 50 layers. This fixes the width and the height of each stack to 0.48 m and 1.59 m, respectively.
The number and length of the membrane channels must be calculated to find the overall dimensions of the membrane module (i.e., membrane area and volume). These two design parameters influence the fluid-dynamic performance, especially in terms of pressure losses. The membrane area is the main input data for the economic analysis.
A Monte Carlo simulation is used to explicit the functional dependency of the membrane area and the outlet pressure to the number and length of the membrane channels. For a fixed value of oxygen molar flow rate, there are multiple combinations of channel length and number of channels in a membrane module that provide the design constraint, as depicted in Figure 15a. These couples of values are characterized by different values of membrane area and overall pressure loss. Therefore, the objective function of the design optimization is the minimization of the product of the membrane area and the pressure loss. Figure 15b,c represent the functional dependence of the total membrane area and the pressure loss with the two design parameters (channel length and number of channels). Both the membrane area and the overall pressure loss find their minimum with the minimum channel length.  To find a physical solution, a minimum channel length is fixed equal to at least the stack width (0.48 m). With these assumptions, the best design is found, as depicted in Figure 15d for the specific conditions of a 3-end membrane operating at 4 bar and with a SF = 40%.
The summary of the geometric characteristics for all considered scenarios is reported in Table 3  Table 3. Design characteristic summary for the membrane modules in Case C and Case D.

Economic Assessment
This section summaries the economic benefits of oxy-fuel glass melting furnaces when oxygen production occurs by means of OTM technology instead of conventional solutions. We focus the comparison between Case B and Case C.
An investment costs (C) of 5 M€ and 1 M€ are assumed for a plant capacity (Q) of 2000 Nm 3 /h and 200 Nm 3 /h, respectively. Assuming an exponential trend in investment cost as a function of the capacity ratio, that is, C = C 0 ·(Q/Q 0 ) f , an exponential factor f = 0.7 can be inferred from the previous values. Considering that Case B requires around 630 Nm 3 h −1 of pure oxygen, it is possible to estimate 2228.6 k€ as the total cost of the plant for oxygen production via conventional technology. In addition, an energy consumption of 0.4 kWh per Nm 3 of pure oxygen is assumed.
Referring to Case C, component costs are more thoroughly broken down. We make reference to the methodology proposed by Galanti and Massardo [52] for a quick estimation of these components costs. The direct material costs for the micro-gas turbine system amounts to 491 k€ (as a reference, a cost of 210 k€ for a 200 kW micro-gas turbine is considered). This cost is increased by 30% to take customization into account and further 65.6 k€ are considered for additional recuperative heat exchangers and oxygen vacuum pumping equipment. The sum of such direct material costs is then increased to include EPC services (32%), construction costs (20%) and other costs (8%), resulting iñ 1126 k€.
A specific cost of 1900 €/m 2 is assumed for the OTM module [53]. Based on the calculated membrane area (484.16 m 2 ), the cost for the membranes results~920 k€, plus 10% as the estimated cost for the containing vessel. The resulting direct material costs for the membrane system is then increased to include EPC services (16%), construction (10%) and other costs (8%). The total cost for the membrane system is~1356 k€. Ultimately the overall Case C system cost results 2482.1 k€.
Additional 93.8 k€ are considered in Case B for recuperative heat exchangers and blowers, not present in Case C's layout. Such direct material costs ultimately come to 150.1 k€ as the total costs for additional plant engineering.
The total plant costs of the two compared technologies, excepting the glass melting furnace and the after-treatment station, do not seem to differ significantly.
The additional following economic assumptions are made: 10-year plant life, 8% discount rate, energy valued at the prices of the Italian market (0.79 € Nm −3 for natural gas, 0.159 € kWh −1 for electricity) and 8000 h of plant operation per year. Economic results are reported in Table 4 for the two cases, considering that the discounted net present value of 10-unit payments equals 6.71 according to the assumed rates. No maintenance costs are considered for the sake of simplicity, even though the membranes could be subject to fouling as well as occlusion. In general, the more complex the system, the greater the maintenance costs. The energy superiority highlighted in Table 2 for the OTM-based oxy-fuel furnace is further remarked by the results in Table 4. Compared to Case B, the solution proposed in this paper is~1.8 M€ cheaper on a total cost basis. Further saving can be calculated in case of larger scale production. In this case, after reducing the cost for the membrane down to 1050 € m −2 [53], the cost of plant in Table 2 can be revised as 1875.5 k€. This revision leads to possible economic saving of~2.4 M€ compared to Case B. Indeed, these figures are based on a long-term stability of the membranes, which has not been secured. Therefore, the discounted total costs in Table 4 could be revised in case one or even two membrane stack replacements are considered. Based on such event, Case C would result in a final cost of 17.11 M€ or 18.47 M€, if one or two membrane system replacements are necessary, respectively. Ultimately, switching to the scenario with lower cost of membrane (1050 € m −2 ), the total discounted costs for Case C result 15.91 M€ and 16.66 M€, highlighting again the promising economic results of the OTM technology.

Conclusions
Oxygen transport membranes are studied in this work for application into the glass production industrial sector. Oxygen membranes are identified as a potential solution for low specific energy consumption and high energy density solution for oxygen separation, in competition with other more conventional systems, such as pressure or vacuum swing adsorption and cryogenic separation.
A detailed mathematical modelling of the physical and chemical processes in an oxygen transport membrane is presented to predict its thermodynamic performance when integrated into a more complex system. Starting from an air-blown reference benchmark case, the mass and energy balances of glass melting furnaces are reproduced for three oxy-fuel furnaces: the first one with oxygen produced by a vacuum swing adsorption unit and the other two where the oxygen is separated from a hot air stream in 3-and 4-end membrane-based systems integrated with a micro-gas turbine.
The assessment of the primary fuel demand of the furnaces highlights the superiority of the oxy-fuel technology compared to the air-blown reference case. The energy demand results to be reduced respectively by~12% for the VSA based plant and~22% for the 3-end membrane system.
The preliminary economic analysis shows the latter one is more convenient in terms of discounted value than the oxy-fuel solution based on a VSA system. On the other hand, calculations of a 4-end membrane system do not suggest similar conclusions, owing to a more complex and capital-intensive system.

Appendix A
This section reports for tables with details of the streams as numbered in Figures 11-14, as well as details of some calculation assumptions for the fluid machinery. Table A1. Mass flow rates, temperatures, pressures and compositions (mol.%) at the main stations in Figure 11 for a glass melting furnace with production of 90 ton day −1 . Two blowers are present in the layout of the reference case A. Auxiliary energy consumptions are calculated based on:

Streamṁ
• a pressure ratio of 1.06 for the blower C1 (air to the burners); • a pressure ratio of 1.1 for the blower C2 (extraction of the gas from the after-treatment station).
Isentropic and mechanical-electrical efficiency values equal to 0.75 and 0.9, respectively, are assumed. Thus, the calculated energy consumptions result in 11.2 kW for the blower C1 and 93.9 kW for the blower C2.
Three main auxiliaries are present in Case B. According to: • a pressure ratio of 1.06 for the blower C1 (oxygen to pre-heater then to burners); • a pressure ratio of 1.1 for the blower C2 (extraction of the gas from the after-treatment station); • a pressure ratio of 1.1 for the blower C3 (delivering the hot carrier through the regenerative heat exchanger for oxygen and fuel pre-heating); Electric energy consumptions of 1.8 kW, 41.3 kW and 4.3 kW are calculated respectively for the three blowers. The auxiliary consumptions total 47.3 kW. In addition, the energy demand related to the oxygen production must be included. Based on the calculations, slightly less than 700 Nm 3 h −1 are required, so the oxygen can be produced by pressure or vacuum swing adsorption technology instead of a small cryogenic distillation process. Assuming an electric consumption of 0.4 kWh per Nm 3 of oxygen, along with an electric conversion efficiency of 50%, the related electric consumption will result equal to 275.2 kW.
As regards the advanced oxy-fuel case with OTM technology in Figure 13, four main auxiliaries are present. The following assumptions have been made: • pressure ratios of 4.43 and 1.1 for the vacuum pump C1 and the blower C2, respectively, with the previously considered efficiency values; • a pressure ratio of 4.2 with an isentropic efficiency of 83% for the compressor C3; • a pressure ratio of 3.88 with an isentropic efficiency of 87% for the turbine T1; • a mechanical-electrical efficiency of 0.9 for the micro-gas turbine including the compressor C3 and the turbine T1. Table A2. Mass flow rates, temperatures, pressures and compositions (mol.%) at the main stations in Figure 12 for a glass melting furnace with production of 90 ton day −1 .  Table A3. Mass flow rates, temperatures, pressures and compositions (mol.%) at the main stations in Figure 13 for a glass melting furnace with production of 90 ton/day. The electric power consumptions have resulted equal to 55.2 kW and 21.9 kW for the vacuum pump C1 and the blower C2, respectively, whereas a net power of 319.6 kW has been calculated for the micro-gas turbine system.

Streamṁ
Ultimately, always four machineries are present in the schematic of the advanced oxy-fuel case with OTM technology in Figure 14. The following assumptions have been made: • pressure ratios of 1.1 for both the blowers C1 and C2, with the previously considered efficiency values; • a pressure ratio of 4.28 (differently from the case in Figure 13, one more passage through heat exchanger is considered here) with an isentropic efficiency of 83% for the compressor C3; • a pressure ratio of 3.88 with an isentropic efficiency of 87% for the turbine T1; • a mechanical-electrical efficiency of 0.9 for the micro-gas turbine including the compressor C3 and the turbine T1.
The electric power consumptions have resulted equal to 3.3 kW and 5.1 kW for the blowers C1 and C2, respectively, whereas a net power of 374.2 kW has been calculated for the micro-gas turbine system. Table A4. Mass flow rates, temperatures, pressures and compositions (mol.%) at the main stations in Figure 14 for a glass melting furnace with production of 90 ton day −1 .