An integrated model of electrochemical cells with co-ionic solid

12 hours ago - One-step Heck Reaction Generates Nonimmunosuppressive FK506 .... PDF (366 KB) ... We consider electrochemical cells with co-ionic solid...
0 downloads 0 Views 377KB Size
Subscriber access provided by Macquarie University

Kinetics, Catalysis, and Reaction Engineering

An integrated model of electrochemical cells with co-ionic solid electrolyte membranes: Coupling of membrane charge-carrier transport and multiple reactions at the triple-phase boundaries Nefeli Charitou, Lydia I Kolitsi, Michael Stoukides, and Stergios Yiantsios Ind. Eng. Chem. Res., Just Accepted Manuscript • DOI: 10.1021/acs.iecr.9b03427 • Publication Date (Web): 29 Aug 2019 Downloaded from pubs.acs.org on August 30, 2019

Just Accepted “Just Accepted” manuscripts have been peer-reviewed and accepted for publication. They are posted online prior to technical editing, formatting for publication and author proofing. The American Chemical Society provides “Just Accepted” as a service to the research community to expedite the dissemination of scientific material as soon as possible after acceptance. “Just Accepted” manuscripts appear in full in PDF format accompanied by an HTML abstract. “Just Accepted” manuscripts have been fully peer reviewed, but should not be considered the official version of record. They are citable by the Digital Object Identifier (DOI®). “Just Accepted” is an optional service offered to authors. Therefore, the “Just Accepted” Web site may not include all articles that will be published in the journal. After a manuscript is technically edited and formatted, it will be removed from the “Just Accepted” Web site and published as an ASAP article. Note that technical editing may introduce minor changes to the manuscript text and/or graphics which could affect content, and all legal disclaimers and ethical guidelines that apply to the journal pertain. ACS cannot be held responsible for errors or consequences arising from the use of information contained in these “Just Accepted” manuscripts.

is published by the American Chemical Society. 1155 Sixteenth Street N.W., Washington, DC 20036 Published by American Chemical Society. Copyright © American Chemical Society. However, no copyright claim is made to original U.S. Government works, or works produced by employees of any Commonwealth realm Crown government in the course of their duties.

Page 1 of 41 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research

An integrated model of electrochemical cells with co-ionic solid electrolyte membranes: Coupling of membrane charge-carrier transport and multiple reactions at the triple-phase boundaries

Nefeli Charitou, Lydia I. Kolitsi, Michael Stoukides, and Stergios G. Yiantsios* Department of Chemical Engineering, Aristotle University of Thessaloniki, Univ. Box 453, GR 541 24, Thessaloniki, Greece.

ABSTRACT: We consider electrochemical cells with co-ionic solid electrolyte membranes and analyze charge transport in the membrane and multiple electrochemical reactions at the triplephase boundaries. Charged-species concentrations at the membrane boundaries are dictated by fast non-electrochemical reactions with the gas phase oxygen, hydrogen and water vapor, which endow the membrane with its co-ionic conductivity properties. The electrochemical reactions are described in terms of the Butler-Volmer equation. Charge conservation is invoked to couple the rates of transport within the membrane and the rates of the electrochemical reactions and close the model system of equations. The electrical potentials of electrodes and membranes at their boundary are unique and, hence, the overpotentials of different reactions cannot in general be simultaneously zero. For a single pair of half-cell reactions the model behavior is as expected, regarding the direction of reactions, the effect of applied voltage and the effect of kinetic constants on overpotentials. When, multiple reactions take place, depending on their relative speeds they may collaborate or compete to determine the net charge flux. The fastest reactions determine the local membrane electrical potential and, hence, the sign and magnitude of the overpotentials of other reactions, and thus their direction.

ACS Paragon Plus Environment

Industrial & Engineering Chemistry Research 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

1. INTRODUCTION Fuel cells convert chemical energy directly into electricity and have the potential to be used in stationary and mobile applications and to replace combustion engines in vehicles1,2. They have attracted research interest for many years because they might offer solutions to concerns on environmental impact and energy efficiency. Of the most common and promising fuel cell types are the solid oxide fuel cells (SOFC’s). At the same time, solid electrolyte cells and, in particular, oxygen ion conductors have been used for several years in the study of catalytic reaction rates and selectivities3,4. Apart from oxygen ion conductors, materials that exhibit protonic conductivity in the solid state have been discovered and introduced into heterogeneous catalytic systems5. These proton conductors are particularly attractive in catalysis because of the temperature range they operate, which favors many industrial hydrogenation and dehydrogenation reactions. Carbon dioxide conversion into useful chemicals is a related issue with obvious environmental implications. Several comprehensive reviews exist on mathematical modeling of SOFCs with oxygen ion conductors, examining three-dimensional or reduced-order models6-8. Detailed models focusing on the workings of microscopic cell sections can be found, as well as more macroscopic ones, describing complete reactors of simpler or more complex configuration9-12. On the other hand, the number of studies on mathematical modeling of SOFCs or electrolytic cells with alternative types of oxides (proton conducting or co-ionic conducting) is relatively limited. As interest in the exploitation of these alternative oxides in electricity production or electrocatalysis grows, there comes the need also to theoretically describe and analyze their behavior.

ACS Paragon Plus Environment

Page 2 of 41

Page 3 of 41 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research

Several studies of fundamental character focus on the transport of multiple ionic species in solid oxide membranes13-21. There have also been a number of studies on electrochemical transport and performance of oxides exhibiting mixed oxygen-ion and electronic conductivity for utilization in fuel cells or batteries and sensors22-26. However, for fuel or electrolytic cells with co-ionic conductors, only very few attempts have been made to comprehensively integrate aspects such as ion transport, electrochemical reaction kinetics, as well as mass, heat and momentum transfer in the flowing reactant and product streams27,28. Moreover, despite, or perhaps due to, the complexity of electrolytic reaction mechanisms, the majority of models focus on simple and single reactions taking place at each of the two sides of the conductor, such as carbon or hydrogen combustion or water dissociation. In particular with electrolytic cells, issues of selectivity necessitate a more detailed and comprehensive consideration and allowance of multiple electrochemical reactions at each electrode interface. Therefore, the scope of the present work is to synthesize some of the key aforementioned aspects towards the development of an integrated model, which, in the framework of continuum concepts, accounts in detail for the transport of several charge-carrying species within the solid electrolyte membranes, as well as for the kinetics of several electrochemical reactions that may take place simultaneously and cooperate or compete. The latter are placed in the framework of the Butler-Volmer equation, as formulated in White et al.29 for water electrolyte solutions and elsewhere30. A key point pertaining to the description of the electrochemical reaction kinetics is that the electrical potentials of electrodes and membranes at their boundary may be different but are unique. Thus, the overpotentials of different reactions cannot in general be simultaneously zero, no matter how large the reaction constants are.

ACS Paragon Plus Environment

Industrial & Engineering Chemistry Research 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The structure of the presentation is as follows: In section 2 the modeling problem is stated and its various aspects are analyzed, starting from the transport of charged species in the membrane, going to their concentrations at the membrane boundaries and closing with the description of the electrochemical reactions at the anode and cathode. The section also contains key aspects on the mathematical and numerical treatment of the system of algebraic and differential equations comprising the model. In section 3 the electrochemical reaction equilibria are considered and then the kinetics and the model behavior are taken up. First, the case is discussed when one pair of half-cell reactions takes place on the anode and cathode. Lastly, the case of more than one reaction at a membrane side is discussed in the light of the effect of kinetic parameters and the developed overpotentials. Finally, a summary and closing comments are presented in section 4.

2. MODEL DEVELOPMENT We consider an elementary segment of an electrochemical cell, comprised of an electrolytic membrane in contact with porous electrodes on both sides and gaseous phases of known composition, as depicted in the Figure 1. Here, the selection of species is only indicative and what is discussed applies to general types of electrolytic or fuel cells. A known voltage is assumed to be applied between the electrodes and, as a result, current flows in the form of electrons through the electrodes from the anode to the cathode and in the form of other charge carrying species within the membrane. For this to be possible, electrochemical reactions take place, which involve the participation of electrons from the electrodes, ionic species from the electrolyte and gaseous species from the gas phase. Thus, these reactions take place at the boundary between electrodes, membrane and the gaseous phase or the triple phase boundary

ACS Paragon Plus Environment

Page 4 of 41

Page 5 of 41 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research

(TPB). Moreover, gaseous species, in particular humidity, react with the membrane at their common interface without the participation of electrons or the electrodes and these nonelectrochemical processes affect the performance and the charge carrier concentration and conductivity of the membrane.

V

e

H2O, O2

H+ O=

CO2, CH4, CO, H2O

Figure 1. Schematic of an elementary electrochemical cell showing the gas phases, the electrodes and the membrane. All these issues must be considered in order to develop a model that can determine the reaction rates and the charge fluxes, given an applied voltage and the compositions of the gaseous phases on either side. Provided that such a model and appropriate relevant data exist, then, at a second step, a macroscopic reactor configuration may be selected and on the basis of mass and energy balances its performance predicted and optimized in terms of conversions, selectivities, etc. 2.1. Charge-carrier transport in the membrane Suppose we have a solid electrolytic membrane, which, in general, may be a co-ionic conductor that allows the transport of protons and oxygen ions. The latter is possible due to the motion of oxygen vacancies. It will be assumed that, in general, four types of charged species are capable

ACS Paragon Plus Environment

Industrial & Engineering Chemistry Research 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 6 of 41

of being transported, namely protons, oxygen vacancies, electron holes and electrons (H+, Vo, p, n). Moreover, it is assumed that electroneutrality holds everywhere, which imposes a constraint on the above concentrations: [H+] + 2[Vo] + [p] – [n] = constant,

(1)

where the constant will be set related to the concentration of the dopant. Finally, it is assumed that a fast reaction is taking place and internal electronic equilibrium always holds: p + n  nil,

(2)

with equilibrium constant Ke17. Therefore, these two species can be dealt with together since [n] = Ke/[p].

(3)

Under steady-state conditions, conservation equations may be written for each of the species in the form

.ji  ri ,

(4)

where ji is the molar flux of species i, and ri is the molar reaction rate. The latter is relevant to electron holes and electrons, while it is zero for protons and oxygen holes. Then the transport of total charge is j q   z i ji .

(5)

i

It follows that .j q  0 ,

(6)

ACS Paragon Plus Environment

Page 7 of 41 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research

i.e. in one dimension there is a constant charge flux from one electrode to the other. The fluxes are given by a generalized form of Fick’s law

ji  

Di ci  i , RT

(7)

where Di is the species diffusion coefficient, ci the concentration and μi is the electrochemical potential

 i   i0  RT ln ai  z i F .

(8)

Here, ai is the activity, zi the species charge, F the Faraday constant, and  the electrical potential. We are going to neglect non-idealities and set the activities equal to concentrations. Moreover, since we assume electroneutrality the electrical potential is governed by

 2  0 .

(9)

Thus, since changes are primarily taking place across the thin dimension of the membrane, the above equation is equivalent to a constant electrical potential gradient and a linear field, which is completely defined by its values at the boundaries, a and c. If we substitute eq (8) into eq (7) we obtain the following system of equations:  d  d [ H  ] F d   [ H  ]   0 , dx  dx RT dx 

(10)

d  d [Vo ] 2 F d   [Vo ]   0 ,  dx  dx RT dx 

(11)

ACS Paragon Plus Environment

Industrial & Engineering Chemistry Research 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

d  d ([ p ]  [n]) F d   ([ p ]  [n])   0  dx  dx RT dx 

Page 8 of 41

(12)

Eq (12) is derived on the basis of equations (4), (7) and (8) separately for the electrons and electron holes by subtracting and using the fact the reaction rates for those two species are equal (eq 2). Note that eq (12) entails the assumption that the diffusivities of electron holes and electrons are equal. If this is not correct, a modified version needs to be written, where one of the two concentrations is multiplied by the ratio of diffusivities. The above equations (10), (11), (12) must be solved together with the equilibrium condition (eq 3) and boundary conditions that are discussed below. 2.2. Charge-carrier concentrations at the membrane/gas interfaces On both sides the membrane is in contact with a porous electrode and with a gaseous phase occupying the electrode pores. Electrochemical reactions may be happening at the points where the three phases meet, called the triple phase boundary (TPB). These are discussed in the next subsection. But in the presence of H2O, H2, O2 in the gas phase, other reactions take place as well, which do not involve the transport of charges across the interface or the presence of electrodes. According to the literature13-18, it is assumed that the following surface reactions take place: 1/2O2 + Vo → O= + 2p

(13)

H2O + Vo + O= → 2OH-

(14)

H2 + 2p → 2H+

(15)

and

ACS Paragon Plus Environment

Page 9 of 41 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research

Furthermore, these reactions will be assumed to be fast enough so that equilibrium is always established17,26,28. Hence, if appropriate equilibrium constants, K1, K2, K3, are available, given the gaseous species concentrations and assuming that O= is fixed within the membrane, the equilibrium relations together with the electroneutrality condition are sufficient to determine the concentrations of protons, oxygen vacancies and electron holes at the interfaces. 2.3. A model of the electrochemical reactions and closure of the system of equations 2.3.1. A single reaction pair To fix ideas, let us assume that a single hypothetical electrochemical reaction is taking place on either side of the membrane. For example, let’s say that at the anode we have H2O → 2H+ + ½ O2 + 2e-,

(16)

where the protons go into the membrane, the electrons go into the electrode and the oxygen goes to the gaseous phase. Similarly, for the sake of discussion, let us assume that at the cathode we have CO2 + 2H+ + 2e- → CO + H2O.

(17)

The rate of the anodic reaction is assumed to be given by the Butler-Volmer equation31, 32:

ra  I 0,a sinh( 

a F 2 RT

(18)

),

where ηa is the overpotential at the anode, defined as

 a   eq ,a   m ,a   e,a .

(19)

ACS Paragon Plus Environment

Industrial & Engineering Chemistry Research 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 10 of 41

Here, the equilibrium potential difference between electrode and membrane is given by the Nernst equation:

 eq ,a  E 0,a 

RT  [O2 ][ H  ] 2 ln 2 F  [ H 2 O ]

  . 

(20)

The equilibrium exchange current is given by



I 0,a  k a [ H 2 O] a [O2 ]1 / 2 [ H  ] 2



(1 a )

(21)

where a is the symmetry factor, taken to be ½. In a similar fashion, at the cathode we have

rc  I 0,c sinh(

c F 2 RT

(22)

),

where ηc is the overpotential at the cathode, defined as

 c   eq ,c   m ,c   e ,c .

(23)

The equilibrium potential difference between electrode and membrane is given again

 eq ,c  E 0,c 

RT  [CO2 ][ H  ] 2 ln 2 F  [ H 2 O ][CO ]

  . 

(24)

The equilibrium exchange current is given by





I 0,c  k c [ H 2 O][CO ] [ H  ] 2 [CO2 ] a

(1 a )

where again we assume a = ½ .

ACS Paragon Plus Environment

(25)

Page 11 of 41 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research

Provided that the compositions of the gaseous phases are known and reaction rate constants are available, we have essentially four unknowns, which are the potentials of the electrodes and the membrane at the two interfaces, namely e,a, e,c m,a, m,c. However, we may suppose that the potential difference V between the two electrodes is given and, moreover, we may set arbitrarily one of them equal to zero since only differences matter. Therefore, we may write

 e , a  V ,  e ,c  0 .

(26)

Finally, the two left unknowns m,a and  m,c, may be determined by equating the charge fluxes at the two interfaces, expressed in terms of membrane transport and reaction rates, i.e.

2ra   ji z i ,

(27)

2rc   ji z i ,

(28)

i

i

which provide the last two equations to close the system. 2.3.2. Multiple electrochemical reactions It is more realistic to consider that several reactions may be taking place simultaneously at the electrodes. For example, at the anode we may have ra1: H2O → 2H+ + ½ O2 + 2e-,

(29)

ra2: O2-→ ½ O2 + 2e- .

(30)

and

Similarly, at the cathode suppose that we have

ACS Paragon Plus Environment

Industrial & Engineering Chemistry Research 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 12 of 41

rc1: CO2 + 2e- → CO + O2-

(31)

rc2: CO2 + 2H+ + 2e- → CO + H2O

(32)

rc3: CO2 + 4H+ + 8e- → CH4 + 2O2-

(33)

Provided that kinetic constants are available, each reaction rate is expressed in a similar fashion as before, through the Butler-Volmer equation and its own equilibrium potential and overpotential, though all containing the unique values of membrane and electrode potentials. What now changes is the charge flux continuity conditions at the two interfaces, which become

2ra1  2ra 2   ji z i ,

(34)

2rc1  2rc 2  8rc 3   ji z i ,

(35)

i

i

to provide closure for the system of equations again. 2.4. Analytical and numerical methods and model solution The procedure to solve the system of algebraic and differential equations comprising the models involves several steps outlined below. 2.4.1. Charge-carrier anodic and cathodic concentrations We start with the determination of charge-carrier concentrations at the membrane boundary with the aid of equilibrium relations for the reactions given by eqs (13), (14), (15). The equilibrium relations are listed below.

[ p]2 , K1  [Vo]PO12/ 2

(36)

ACS Paragon Plus Environment

Page 13 of 41 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research

K2 

[ H  ]2 , [Vo]PH 2O

(37)

[ H  ]2 , K3  [ p ] 2 PH 2

(38)

Here the fixed oxygen ion concentration in the membrane has been incorporated in the equilibrium constants. The above are used together with the equilibrium condition (eq 3) and the neutrality condition (eq 1). The method to solve the system is iterative and proceeds in the following way. From eq (37) we have [Vo]  [ H  ] 2 / PH 2O / K 2 . Substituting into the electroneutrality condition we obtain a quadratic equation for [H+], namely:

(

2 )[ H  ] 2  [ H  ]  [ p]  [n]  [dopant ]  0 , K 2 PH 2O

(39)

Given guesses for [p] and [n] this is solved and the positive root provides an estimate for [H+]. Then on the basis of this estimate, [Vo] is obtained from eq (37). Lastly, improved guesses [p] and [n] are obtained from eqs (36) and (3) and the cycle is repeated until convergence, which is very fast. 2.4.2. Equilibrium potentials of electrochemical reactions Having obtained surface concentrations for all the species involved, the second step is to find the equilibrium potentials for the electrochemical reactions through the Nernst equations (eqs (20) and (24)). The standard potentials at the temperature under consideration are evaluated from the equation

ACS Paragon Plus Environment

Industrial & Engineering Chemistry Research 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

T T  C P H 0  TS 0 1  0  E   E 25   C dT  T dT P 25 T  , nF nF  25 0

Page 14 of 41

(40)

where the specific heats are expressed as33:

C P (T )  a  bT 1 / 4  cT 1 / 2  dT 3 / 4 .

(41)

The integrals may be evaluated analytically. 2.4.3. Charge-carrier concentration distributions and charge flux The next step is to evaluate the charge flux through the membrane and this involves all the charge carriers. For this task we need the potential drop across the membrane, which at this stage is tentative, by introducing guesses for the membrane electrical potentials at the anode and cathode, m,a and m,c. Thus, the constant potential gradient is /L, where L is the membrane thickness. Eqs (10) and (11) can be integrated analytically to give

[ H  ]( x)  [ H  ] a 

(1  e



(1  e

F  x RT L F  RT

)

([ H  ]c  [ H  ] a ) ,

(42)

)

and

[VO ]( x)  [VO ] a 

(1  e



(1  e

2 F  x RT L



2 F RT

)

([VO ]c  [VO ] a ) .

(43)

)

Eq (12) together with the equilibrium condition (eq 3) for electron holes and electrons is solved numerically with the help of a shooting technique and the Newton-Raphson iterative method. Integrating once we have

ACS Paragon Plus Environment

Page 15 of 41 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research

d ([ p ]  [n]) F   ([ p ]  [n])  j pn . dx RT L

(44)

where jpn is an unknown constant, proportional to the relevant flux. We set [p] - [n] = y and with the aid of the equilibrium condition (eq 3) we replace [p] + [n] in terms of y to obtain dy F   dx RT L

y 2  4K e 2

 j pn .

(45)

Having obtained the values of [p] and [n] at the anode, and hence of y also, the above nonlinear equation can be integrated towards the cathode to obtain a tentative boundary value, yc. For this task we need to provide a guess for the constant jpn. Then, through a Newton-Raphson procedure the error between yc and the known boundary value at the cathode provides a correction, jpn, for the constant jpn :

[ p ] c  [ n] c  y c 

dy c j pn . dj pn

(46)

Here, the derivative is evaluated numerically by perturbing jpn by a small fixed value during each iteration. 2.4.4. Electrochemical reaction rates and system closure Knowing the concentration of all species, the equilibrium potentials for all the electrochemical reactions and having assumed the membrane potentials, one can calculate the reaction rates from the corresponding Butler-Volmer equations. Then the last step of the numerical procedure is to update the guesses for the membrane potentials by enforcing the charge conservation equations at the anode and cathode. This step is implemented again through a Newton-Raphson procedure. More specifically, we define the charge conservation errors at the anode and cathode as

ACS Paragon Plus Environment

Industrial & Engineering Chemistry Research 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 16 of 41

f a   ji z i   a k rka ,

(47)

f c   ji z i   a k rkc .

(48)

i

i

k

k

Here, ak is the electron stoichiometric coefficient in the respective reactions. Then, since these errors must be driven to zero we write

0  f a ( m ,a ,  m ,c ) 

f a f  a  a  c ,  a  c

(49)

0  f c ( m ,a ,  m ,c ) 

f c f  a  c  c .  a  c

(50)

This forms a linear system for the corrections a, c to the membrane potentials and the whole scheme is repeated until convergence. The partial derivatives in the above equations are evaluated numerically by perturbing the assumed values of the membrane potentials in each iteration. To summarise: Given kinetic, thermodynamic data and the applied voltage, the membrane concentrations are determined first, as explained in section 2.4.1. Thereafter the equilibrium potentials are determined as outlined in section 2.4.2. In the third step, assuming the membrane potentials, charge fluxes are obtained as analyzed in section 2.4.3 and reaction rates are obtained on the basis of the Butler-Volmer equations. Finally, the charge conservation errors are evaluated on the basis of eqs (47-48). If these errors exceed preset tolerances the membrane potentials are updated on the basis of eqs (49-50) and the iterative scheme goes back to the third step.

ACS Paragon Plus Environment

Page 17 of 41 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research

3. RESULTS AND DISCUSSION In the following results, which have only the character of demonstration of the model features, the assumption was made that the gas phase at the anode is composed mainly of water vapor and also some oxygen. The partial pressures were taken to be 1 and 0.001 Atm, respectively. Similarly, the gas phase at the cathode was assumed to be dominated by CO2, and also to contain traces of H2O, H2, CH4. The partial pressures of those species were assumed to be 1, 10-3, 10-6, and 10-6 Atm respectively. The solid electrolyte membrane SrCe0.95Yb0.05O2.975 considered here has been taken from Tan et al.17, who provide data for the composition, species diffusivities and equilibrium constants of the surface reactions with water vapor, oxygen and hydrogen. According to Uchida et al.34, mixed proton–hole conductivity is often observed in perovskite type materials (ΑΒΧ3) and especially in perovskite oxides such as ABO3, at high temperatures. More specifically, SrCe0.95Yb0.05O2.975 (SCY) is a common example of this category. SCY membranes demonstrate hole conductivity (p–type) when surrounded by a water- or hydrogen-free environment. Alternatively, protons become the predominant charged carrier. At the same time, when Ce4+ is replaced by cationic dopants, such as Yb3+, the production of oxygen vacancies or holes is enhanced and proton conductivity is greatly affected. Furthermore, it is stated in literature that mixed proton-hole conductivity in materials can be transformed gradually to oxygen ion-electron conductivity, when the temperature is increased. Consequently, the desired co-ionic behavior can be accomplished at some intemediate temperature range17,34.

ACS Paragon Plus Environment

Industrial & Engineering Chemistry Research 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 18 of 41

Data for the diffusivities of charge carriers in the membrane are shown in Table 1. Table 1. Diffusion coefficient of charged species17. 𝑐𝑚2

diffusion coefficient (

𝑠

protons Η+ oxygen vacancies VO∙∙ electrons n/holes p

)

value

0.122 ∙ 𝑒 24.24 ∙ 𝑒 49.38 ∙ 𝑒

ACS Paragon Plus Environment







7,851 𝑇

23,467 𝑇 12,588 𝑇

Page 19 of 41 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research

Equilibrium data for the non-electrochemical reactions taking place at the membrane gas interface are presented in Table 2. Table 2. Equilibrium constants of anodic and cathodic reactions and concentration of dopant17.

description

value 𝑚𝑜𝑙

equilibrium constant Κ1 (

) (eq 13)

1

𝑐𝑚3 ∙ 𝑎𝑡𝑚2 𝑚𝑜𝑙

equilibrium constant Κ2 (

) (eq 14)

1

𝑐𝑚3 ∙ 𝑎𝑡𝑚2 𝑚𝑜𝑙

equilibrium constant Κ3 (

1,371.3

3.24 ∙ 10

―8

𝛵

∙𝑒

1368.4

2.087 ∙ 10

―2

𝛵

∙𝑒

𝐾𝑤 ∙ 𝐾2

) (eq 15)

1

𝐾1

𝑐𝑚3 ∙ 𝑎𝑡𝑚2

equilibrium constant of water formation 29,589

reaction Κw (

𝑚𝑜𝑙

0.536 ∙ 𝑒

)

1

𝛵

𝑐𝑚3 ∙ 𝑎𝑡𝑚2

equilibrium constant Κe (

𝑚𝑜𝑙

) (eq 2)

1

𝑐𝑚3 ∙ 𝑎𝑡𝑚2 𝑚𝑜𝑙

concentration of ytterbium [ΥbCe΄] (𝑐𝑚3)

―28,520

1.14 ∙ 10 ―6 ∙ 𝑒

𝛵

1.06 ∙ 10 ―3

In the developed model, we tentatively considered six electrochemical reactions in total, anodic and cathodic. These are presented in Table 3. Of course, the above list is not exhaustive nor necessarily realistic for all conditions, membranes and catalysts. However, given the complexity of the system and the multitude of the parameters, what we aim at is an understanding of the model behavior and major qualtitative features, rather than presenting a tool with quantititive predictive capabilities. For this reason, the emphasis has been placed only on some of those

ACS Paragon Plus Environment

Industrial & Engineering Chemistry Research 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 20 of 41

reactions, as will become apparent in the subsequent discussion. Data for the equilibrium potentials at standard conditions and unit activities for the reactions considered are also shown in Table 3. Table 3. Anodic and cathodic electrochemical reactions and their standard equilibrium potentials.

reaction

Ε025 (V)

Ref.

a1

H2O → 2H+ + ½ O2 + 2e-

-1.229

3535

a2

O2-→ ½ O2 + 2e-

-0.688

3636

c1

CO2 + 2e → CO + O2-

-0.645

35,3635,36

c2

CO2 + 2H+ + 2e- → CO + H2O

-0.104

3535

c3

CO2 + 4H+ + 8e- → CH4 + 2O2-

-0.372

35, 3635,36

c4

2H+ + 2e- → H2

0

3535

3.1. Charge-carrier concentrations and membrane fluxes On the basis of eqs (13-15), charge-carrier concentrations at the anode and cathode are presented in Figure 2 as functions of temperature. The membrane thickness was taken to be 0.4 mm. It can be seen that protons and oxygen vacancies dominate. Moreover, they change little over the temperature range selected because the dopant concentration is fixed and the equilibrium constant K2 changes only by a factor of about 2.

ACS Paragon Plus Environment

Page 21 of 41

Anode

Cathode 0.001

(a)

3

3

Concentration (mol/cm )

0.0001

Concentration (mol/cm )

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research

-6

10

-8

10

-10

10

-12

10

-14

10

protons oxygen vacancies electron holes electrons

-16

10

-18

10

600

(b)

-5

10

-7

10

-9

10

-11

10

protons oxygen vacancies electron holes electrons

-13

700

800

900

1000

Temperature (K)

1100

1200

10

600

700

800

900

1000

Temperature (K)

1100

1200

Figure 2. Charge-carrier concentrations at the anode (a) and the cathode (b), respectively, as functions of temperature. These surface concentrations, which serve as boundary conditions, can be used to obtain the charge carrier concentrations within the membrane and give the individual and total charge fluxes, as discussed in the previous section. These results are shown in Figure 3 for the case when 1 V of potential difference exists between the two sides of the membrane. The boundarylayer type behavior of the concentrations near the cathode is apparent. It may be observed that proton conductivity dominates due to the high mobility of those ions within the crystal structure of the perovskitic membrane.

ACS Paragon Plus Environment

Industrial & Engineering Chemistry Research

(a)

0.001

0.0001

Fluxes (mol/cm *s)

3

Concentration (mol/cm )

-5

10

-7

10

-9

10

protons oxygen vacancies electron holes electrons

-11

10

-6

10

-8

10

-10

10

-12

10

protons oxygen vacancies electrons/electron holes

-14

10

-16

10

-13

10

(b)

2

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 22 of 41

0

0.2

0.4

0.6

0.8

1

600

Dimensionless position across membrane

800

1000

1200

1400

1600

Temperature (K)

Figure 3. (a) Distributions of charge-carrier concentrations within the membrane and (b) total molar fluxes. The potential drop across the membrane is taken to be 1 Volt. 3.2. Electrochemical reactions and equilibrium potentials In Figure 4 the equilibrium potentials of the above half-cell reactions are presented as functions of temperature.

ACS Paragon Plus Environment

Page 23 of 41

1

Electrode potential (V)

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research

0 -1 -2 -3 -4 -5 600

Ea1 Ec2 Ec3 Ec1

800

1000

1200

1400

1600

Temperature (K) Figure 4. Equilibrium electrode potential as a function of temperature for half-cell electrochemical reactions.

The information contained in this diagram is a useful guide as to the behavior of individual pairs of reactions when these are considered in isolation. The difference between the anodic and cathodic equilibrium potential at any given temperature can be interpreted as the open circuit potential that would be measured if such a system were assembled. The sign of the difference also tells whether a given pair of reactions would proceed in the forward or reverse direction if the circuit were closed but without applying an external voltage. Lastly, the sign and the magnitude of the difference suggest the voltage that needs to be applied to reverse the direction of a reaction pair. These features are indeed borne out by the model predictions, as will be discussed in the next section, where we examine the case of single reactions taking place at the anode and cathode. Note here that by forward direction it is meant that protons move within the membrane from the anode to the cathode and electrons are pumped through the electrical circuit

ACS Paragon Plus Environment

Industrial & Engineering Chemistry Research 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

from the anode to the cathode. In other words, both anodic and cathodic reactions proceed to the right as written in Table 3. 3.3. Electrochemical reaction rates for single anodic and cathodic reaction pairs 3.3.1. Spontaneous direction of reactions and effect of the applied voltage By examining the pair of reactions a1 and c1 in Figure 4 we can see that the equilibrium potential difference is negative at low temperature and changes sign above 850 oC. This means that under no applied external potential the reactions would spontaneously proceed in the forward direction at relatively low temperatures and in the reverse direction near and above 900 oC.

Indeed, this is the model prediction as may be seen from Figure 5a, for the anodic and

cathodic reaction rates. These rates are equal since both reactions involve two electrons in their stoichiometry. On the other hand, if we consider the pair a1-c2, the difference in the equilibrium potentials is positive at low temperatures and changes sign at a temperature near 1100 oC. Again, as may be observed from Figure 5b, the model prediction for the direction of the reactions is compatible with the expected behavior.

ACS Paragon Plus Environment

Page 24 of 41

Page 25 of 41

-6

-7

3.5 10

3.5 10

(a)

o

-7

T = 800 C o T = 900 C

3 10

-6

3 10

2.5 10

Reaction rate

-7

2 10

-7

-6

2 10

-6

1.5 10

1.5 10

-7

1 10

-8

-6

1 10

-7

5 10

5 10

0

0

-7

-8

-5 10

(b)

o

T = 800 C o T = 1200 C

-6

-7

2.5 10

Reaction rate

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research

-8

10

-7

10

-6

10

-5

10

0.0001 0.001 0.01

0.1

1

-5 10 0.0001 0.001

0.01

0.1

1

10

100

1000

kc2

kc1

Figure 5. Reaction rates of different pairs when no external current is applied. (a) a1-c1 pair as a function of kc1 at 800 and 900 K. (b) a1-c2 pair as a function of kc2 at 800 and 1200 K. If we now consider the pair a1-c3, then Figure 4 suggests that the reactions proceed spontaneously in the reverse direction and a voltage close to 3 V needs to be applied in order to change their direction. Indeed, as Figure 6 shows, when the applied voltage is 1 V both reaction rates are negative, whereas when the applied voltage is increased to 3 V the reactions proceed to the forward direction. Note here that the two reactions rates are not equal but proportional, since the number of electrons involved in reaction c3 is four times that of a1.

ACS Paragon Plus Environment

Industrial & Engineering Chemistry Research

-7

5 10

-6

4 10

Reaction rate (mol/(cm *s))

(a) 0

2

2

Reaction rate (mol/(cm *s))

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

-7

-5 10

-6

-1 10

-6

-1.5 10

ra1 rc3 -6

-2 10

Page 26 of 41

(b)

ra1 rc3

-6

3 10

-6

2 10

-6

1 10

0

-6

-9

10

-7

10

-5

10

0.001

0.1

10

1000

5

10

-1 10

-9

10

-7

10

kc3

-5

10

0.001

0.1

10

1000

5

10

7

10

kc3

Figure 6. Rates of reactions a1 and c3 as functions of reaction coefficient kc3. (a) Applied voltage V = 1 Volt, (b) V= 3 Volt.

3.3.2. Effect of the reaction kinetic constants We now consider the single pair of reactions a1 and c1 and focus on the effect of the kinetic constants in the reaction rates, which are expressed in terms of the Butler-Volmer equations (eqs 18 and 22). Both constants are varied and the results are shown in Figure 7 in the form of a family of curves representing reaction rates as functions of kc1 for various fixed values of ka1. Calculations for two values of the applied voltage, namely 1 and 2 V, are presented separately. For any specific value of ka1 the reaction rates display the same behavior, namely they first increase to finite values as kc1 is increased and then they reach a plateau. This means that the cathodic reaction step becomes fast enough and eventually not rate-controlling. The higher the value of ka1 the higher is also the plateau. But when ka1 exceeds a certain threshold the curves tend to be superposed and an upper limit is reached. When this happens the reaction at the anode is also sufficiently fast so that the rate-controlling step becomes the transport of charge carriers

ACS Paragon Plus Environment

Page 27 of 41

within the membrane. Indeed, as may be seen from Figure 7b, if the applied voltage is increased to 2 V the limiting reaction rates double as well. -5

-5

2 10

2 10

(a)

-5

ka1=10 -4 ka1=10 -3 ka1=10 -2 ka1=10 -1 ka1=10 ka1=1 ka1=10

-5

1.5 10

-5

1 10

-6

ka1=10 -4 ka1=10 -3 ka1=10 ka1=10-2 -1 ka1=10 ka1=1 ka1=10

-5

1 10

-6

5 10

0 -15 -13 -11 10 10 10

(b)

-5

Reaction rate

-5

1.5 10

Reaction rate

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research

5 10

-9

10

-7

10

-5

10

0.001

0.1

10

0 -15 -13 -11 10 10 10

-9

10

-7

10

-5

10

0.001

0.1

10

kc1

kc1

Figure 7. Rates of reactions a1 and c1 as functions of reaction coefficient kc1 for different values of coefficient ka1. (a) Applied external voltage V= 1 Volt. (b) V = 2 Volt.

The above observations are supported and more completely understood by the calculations of the overpotentials at the anode and cathode, as presented in Figure 8. As kc1 and, hence, the speed of cathodic reaction increase, the overpotential at the anode adjusts by increasing as well. This leads also to a concomitant increase in charge flux since the electrical potential difference within the membrane increases. Since a balance of the three steps must exist, namely the two reaction rates and the charge flux, the increase of kc1 leads also to the reduction of the overpotential at the cathode, as seen in Figure 8b. Thus, a plateau is reached, which becomes lower for increasing values of ka1. Eventually, when both reaction constants are sufficiently high, both overpotentials tend to zero, as shown in the figure and as expected. Then, the membrane potentials have reached fixed values and the charge transport through the membrane has become the ratelimiting step. Thus, the low plateau is cathodic-reaction controlled, the high plateau is anodic-

ACS Paragon Plus Environment

Industrial & Engineering Chemistry Research

reaction controlled, and when the anodic reaction also becomes fast enough charge transport takes over.

1.2

0.5 -5

0.8

ka1=10 -4 ka1=10 -3 ka1=10 -2 ka1=10 -1 ka1=10 ka1=1

-5

cathode overpotential

1

anode overpotential

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 28 of 41

0.6 0.4 0.2 0 -0.2 -15 -13 -11 10 10 10

0

ka1=10 -4 ka1=10 -3 ka1=10 -2 ka1=10 -1 ka1=10 ka1=1

-0.5

-1

(b)

(a) -9

10

-7

10

kc1

-5

10

0.001

0.1

10

-1.5 -15 -13 -11 10 10 10

-9

10

-7

10

-5

10

0.001

0.1

10

kc1

Figure 8. Anodic (a) and cathodic (b) overpotential for the combination of reactions a1 and c1 as a function of reaction coefficient kc1 for different values of the coefficient ka1 and for external voltage equal to 1 Volt. The arrows indicate the direction of increasing ka1.

3.4. Multiple electrochemical reactions at an electrode/membrane/gas phase interface Finally, we examine the case of multiple reactions taking place at an electrode. Now the parameter space increases dramatically and for this reason we aim only at identifying and highlighting some key features of the predicted behavior. Added to the reaction pair a1-c1 discussed in the previous subsection, we now assume that a second reaction, namely c2, can also take place at the cathode. It is also assumed that the anodic reaction is fast enough and not ratelimiting by setting ka1 equal to 1, which was found in the previous subsection to be sufficiently high. The system behavior is exhibited in the set of Figures 9(a-d), where, for different fixed values of kc1 anodic and cathodic reaction rates are presented as functions of the rate constant kc2. The selected fixed values of kc1 are 10-6, 10-5, 10-4, 10-3, thus covering a range where, as

ACS Paragon Plus Environment

Page 29 of 41 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research

discussed in the previous subsection, c1 changes behavior from rate-controlling to sufficiently fast. At low speeds of reaction c2, as seen from Figure 9a, the rates of reactions a1 and c1 are equal and can be found also from Figure 7. They increase from Figure 9a to Figure 9d in accordance with kc1. Now, as reaction c2 comes into play by increasing kc2 from negligible to significant values, its rate increases whereas that of c1 decreases. Still the net outcome is a net increase of charge flux and the anodic reaction, which almost double. Thus, under these conditions the two reactions collaborate, although c2 eventually becomes dominant for large values of kc2. A similar behavior is observed in Figure 9b where kc1 is now taken larger. Again the two reactions collaborate, although now c1 has a larger share in the overall charge flux. The trends change in Figure 9c where kc1 is still larger. Now reaction c2 tends to proceed in the opposite direction and thus compete with c1. Eventually both reaction rates reach a plateau, where they are large and opposite. The overall charge flux, represented by the rate of the anodic reaction, remains almost unchanged. Lastly, in Figure 9d the same trend but with even larger opposite rates is observed if kc1 is further increased.

ACS Paragon Plus Environment

Industrial & Engineering Chemistry Research

-6

7 10

-6

7 10

a1 c2 c1

-6

6 10

-6

6 10

-6

5 10

-6

Reaction rates

Reaction rates

5 10

-6

4 10

-6

3 10

-6

2 10

-6

0

4 10

-6

3 10

-6

2 10

-6

0

(a)

-6

-1 10

a1 c2 c1

-6

1 10

1 10

(b)

-6

-7

10

-5

10

0.001

0.1

10

1000

5

10

7

10

-1 10

9

10

-7

10

-5

10

0.001

0.1

kc2

-5

1000

5

10

7

10

9

10

0.0004 a1 c2 c1

-5

3 10

a1 c2 c1

0.0003 0.0002

-5

Reaction rates

2 10

-5

1 10

0 -5

-1 10

0.0001 0 -0.0001 -0.0002

-5

-2 10

-7

10

-5

10

0.001

0.1

10

1000

5

10

7

10

(d)

-0.0003

(c)

-5

-3 10

10

kc2

4 10

Reaction rates

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 30 of 41

9

10

-0.0004 -7 10

-5

10

0.001

0.1

10

1000

5

10

7

10

9

10

kc2

kc2

Figure 9. Rates of reactions a1, c1 and c2 as a function of reaction coefficient kc2 for external applied voltage V = 1 Volt and fixed values of kc1. (a) kc1 = 10-6, (b) kc1 = 10-5, (c) kc1 = 10-4, (d) kc1 = 10-3.

Suppose that a small-scale experimental reactor with well-mixed gas chambers exists, where one would attempt to materialize the above simulated conditions. In practice, the prediction that one of the competing reactions proceeds in the reverse direction would imply that the corresponding species would be depleted, thus pushing the selectivity towards the species related to the forward proceeding competitor.

ACS Paragon Plus Environment

Page 31 of 41 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research

What has been presented in Figure 9 can be rationalized if it is examined together with the overpotentials of the two cathodic reactions, as shown in Figure 10. It is worth emphasizing that since the electrode and membrane potentials at the triple-phase boundary are unique, the overpotentials of the two reactions cannot in general be simultaneously zero, no matter how large the reaction constants are. When the value of kc2 is small, the membrane potential is determined by the characteristics of reaction c1. In Figure 10a we see the overpotentials of c2 and c1, which initially are both finite and negative. Their difference, which depends only on local concentrations and not on the values of membrane or electrode potential, remains always fixed. As kc2 is increased the membrane potential adjusts so that the overpotential of c2 tends eventually to zero and the one of c1 changes accordingly but always remains negative. Thus, the two reactions always proceed in the same direction and collaborate in the net outcome as regards the charge flux and the reaction rate at the anode. The same behavior is manifested in Figure 10b, where kc1 has been increased by an order of magnitude. Now, with reaction c1 being faster, the overpotentials are initially smaller, but still both negative. When kc1 is further increased in Figure 10c, the overpotential of reaction c1 is initially negative but even smaller. Thus, since the difference is necessarily fixed, now the overpotential of reaction c2 is positive. Hence, if the c2 reaction rate constant is significant it proceeds in the reverse direction and now competes with c1. Eventually, with kc2 being increased, the overpotential of reaction c2 tends to zero, which means that the overpotential of c1 becomes more negative, thus enhancing its rate and also the competition between the two reactions. Lastly, the same behavior is observed when kc1 is further increased as seen in Figure 10d. Now the membrane potential is such that the overpotential of c1 is close to zero and that of

ACS Paragon Plus Environment

Industrial & Engineering Chemistry Research

c2 even higher. Again, the same transition is observed, with the two competing reactions now having even higher rates.

0.1

0.2 c2 c1

-0.2 -0.4 -0.6 -0.8 -1 -7 10

-5

0.001

0.1

10

1000

5

10

7

10

-0.1 -0.2 -0.3 -0.4

(a) 10

c2 c1

0

Overpotential (V)

Overpotential (V)

0

-0.5 -7 10

9

10

(b) -5

10

0.001

0.1

kc2

5

1000

10

1000

10

7

10

9

10

0.4 c2 c1

0.2

0.2

Overpotential (V)

c2 c1

0.1 0 -0.1 -0.2 -0.3

0

-0.2

-0.4

-0.4 -0.5 -7 10

10

kc2

0.3

Overpotential (V)

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 32 of 41

(c) -5

10

0.001

0.1

10

1000

5

10

7

10

9

10

-0.6 -7 10

(d) -5

10

0.001

kc2

0.1

10

5

7

10

9

10

kc2

Figure 10. Overpotential of reactions c1 and c2 as a function of reaction coefficient kc2 for external applied voltage V = 1 Volt and fixed values of kc1. (a) kc1 = 10-6, (b) kc1 = 10-5, (c) kc1 = 10-4, (d) kc1 = 10-3.

Lastly, we examined the same systems with the same parameters but simply increased the applied voltage to 2 V. The results showed exactly the same trends and transitions. Simply, the

ACS Paragon Plus Environment

Page 33 of 41 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research

threshold value of kc2 where reaction c2 changes direction is modified to a slightly larger value and also the net charge flux reflected in the rate of reaction a1 is increased. 4. CONCLUSIONS The development of a generalized integrated model of electrochemical cells with co-ionic solid electrolyte membranes is the major feature of the present work. The model accounts for charged species transport in the membrane, including protons, oxygen vacancies, electrons and electron holes. Their concentrations at the membrane boundaries are dictated by fast non-electrochemical reactions with the gas phase oxygen, hydrogen and water vapor, which endow the membrane with its co-ionic conductivity properties. At the same time, multiple electrochemical reactions may take place at the electrode/membrane/gas-phase boundaries, which are described in terms of respective Butler-Volmer equations. A key point pertaining to the description of the electrochemical reaction kinetics is that the electrical potentials of electrodes and membranes at their boundary may be different but are unique. Thus, the overpotentials of different reactions cannot in general be simultaneously zero, no matter how large the reaction constants are. Charge conservation is invoked to couple the rates of transport within the membrane and the rates of the electrochemical reactions. If the speed of the non-electrochemical reactions is assumed to be finite, the model can be extended by invoking proton and oxygen-ion conservation at the membrane-electrode interfaces, in addition to charge conservation. The most obvious simplification of the proposed approach is when a single pair of half-cell reactions is considered to be dominant, as is usually the case in the modeling literature. Moreover, when the applied overpotentials are not significant membrane resistance becomes the

ACS Paragon Plus Environment

Industrial & Engineering Chemistry Research 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

limiting factor. In all those cases the model behavior is as expected with regard to the direction of the reactions, the effect of the applied voltage and the effect of kinetic constants on the developed overpotentials and the relative significance of membrane resistance to the overall rates. When, multiple reactions take place, depending on their relative speeds they may collaborate or compete to determine the net charge flux. The fastest reactions determine the local membrane electrical potential and, hence, the sign and magnitude of the overpotentials of other reactions, and thus their direction, as well. The prediction that some of the reactions may proceed in the reverse direction in practice would imply that the corresponding species would tend to be depleted, thus pushing the selectivities in favor of the species related to the forward proceeding competitors. Such an integrated model may be helpful in analyzing experimental data and guiding experimental procedures and conditions to explore in the quest of optimal selectivities. It may be also helpful by serving as a module in the development of design models for macroscopic reactor configurations, where mass, energy and momentum balances must be incorporated together with features, such as gas transport through the porous electrodes, etc.

AUTHOR INFORMATION Corresponding Author: *E-mail: [email protected]

ACS Paragon Plus Environment

Page 34 of 41

Page 35 of 41 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research

ACKNOWLEDGMENTS This work was supported by the General Secretariat of Research and Technology, Greece, in the framework of the Bilateral and Multilateral Greek-German Research and Innovation cooperation for the project “Proton and oxygen co-ionic conductors for CO2/H2O co-electrolysis and intermittent RES conversion to methanol and other chemicals towards EU Sustainability”.

NOMENCLATURE c = molar concentration (kmol m-3) CP = specific heat (J kmol-1 K-1) D = diffusivity (m2 s-1)

E = Nernst equilibrium potential (V) F = Faraday’s constant (A s kmol-1) H = molar enthalpy (J kmol-1) I = exchange current density (A m-2) j = molar flux (kmol m-2 s-1) k = reaction rate constant (reaction dependent) K = reaction equilibrium constant (reaction dependent) L = membrane thickness (m)

ACS Paragon Plus Environment

Industrial & Engineering Chemistry Research 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

P = partial pressure (Pa) r = electrochemical reaction rate (kmol m-2 s-1) R = gas constant (J kmol-1 K-1) S = molar entropy (J kmol-1 K-1) T = temperature (K) V = applied voltage (V) z = species ionic number Greek Letters

 = overpotential (V)  = chemical potential (J kmol-1)  = electrical potential (V) Subscripts a = anode c = cathode e = electrode m = membrane q = charge

ACS Paragon Plus Environment

Page 36 of 41

Page 37 of 41 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research

eq = equilibrium

REFERENCES (1) Yamamoto, O. Solid Oxide Fuel Cells: Fundamental Aspects and Prospects. Electrochim. Acta 2000, 45 (15-16), 2423–2435. (2) Cheddie, D.; Munroe, N. Review and Comparison of Approaches to Proton Exchange Membrane Fuel Cell Modeling. J. Power Sources 2005, 147 (1-2), 72–84. (3) Garagounis, I.; Kyriakou, V.; Anagnostou, C.; Bourganis, V.; Papachristou, I.; Stoukides, M. Solid Electrolytes: Applications in Heterogeneous Catalysis and Chemical Cogeneration. Ind. Eng. Chem. Res. 2011, 50, 431–472. (4) Marnellos, G.; Sanopoulou, O.; Rizou, A.; Stoukides, M. The Use of Proton Conducting Solid Electrolytes for Improved Performance of Dehydrogenation Reactors. Solid State Ionics 1997, 97, 375. (5) Stoukides, M. Electrochemical Studies of Methane Activation. J. Applied Electrochem. 1995, 25, 899. (6) Bavarian, M.; Soroush, M.; Kevrekidis, I. G.; Benziger, J. B. Mathematical Modeling, Steady-State and Dynamic Behavior, and Control of Fuel Cells: A review. Ind. Eng. Chem. Res. 2010, 49 (17), 7922−7950. (7) Bhattacharyya, D.; Rengaswamy, R. A Review of Solid Oxide Fuel Cell (SOFC) Dynamic Models. Ind. Eng. Chem. Res. 2009, 48 (13), 6068−6086.

ACS Paragon Plus Environment

Industrial & Engineering Chemistry Research 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

(8) Hajimolana, S. A.; Hussain, M. A.; Daud, W. M. A. W.; Soroush, M.; Shamiri, A. Mathematical Modeling of Solid oxide Fuel Cells: A Review. Renew. Sustain. Energy Rev. 2011, 15 (4), 1893−1917. (9) Debenedetti P. G.; Vayenas C. G. Steady-State Analysis of High Temperature Fuel Cells. Chem. Engng. Sci. 1983, 38 (11), 1817-1829. (10) Vayenas C. G.; Debenedetti, P. G.; Yentekakis, I.; Hegedus, L. L. Cross-Flow Solid-State Electrochemical Reactors: A Steady-State Analysis. Ind. Eng. Chem. Fundam. 1985, 24, 316-324. (11) McKenna, E. A.; Othoneos, A.; Kiratzis, N.; Stoukides, M. Synthesis of HCN in a SolidElectrolyte-Cell Reactor. Ind. Eng. Chem. Res. 1993, 32, 1904-1913. (12) Costamagna, P.; Honegger, K. Modeling of Solid Oxide Heat Exchanger Integrated Stacks and Simulation at High Fuel Utilization. J. Electrochem. Soc 1998, 145, 3995. (13) Norby, T. Proton Conduction in Oxides. Solid State Ionics 1990, 40/41, 857-862. (14) Bonanos, N. Transport Properties and Conduction Mechanisms in High-temperature Protonic Conductors. Solid State Ionics 1992, 53-56, 967-974. (15) Ahlgren, E. O. Thermoelectric Power of Gd-doped Barium Cerate. J. Phys. Chem. Solids 1997, 58 (9), 1475-1480. (16) Glockner, R.; Islam, M. S.; Norby T. Protons and Other Defects in BaCeO: a Computational Study. Solid State Ionics 1999, 122, 145–156. (17) Tan, X.; Liu, S.; Li, K.; Hughes, R. Theoretical Analysis of Ion Permeation through Mixed Conducting Membranes and its Application to Dehydrogenation Reactions. Solid State Ionics 2000, 138, 149–159.

ACS Paragon Plus Environment

Page 38 of 41

Page 39 of 41 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research

(18) Song, S. J.; Wachsman, E. D.; Dorris, S. E.; Balachandran, U. Defect Chemistry Modeling of High-temperature Proton-conducting Cerates. Solid State Ionics 2002, 149, 1 – 10. (19) Kreuer, K.D. Aspects of the Formation and Mobility of Protonic Charge Carriers and the Stability of Perovskite-type Oxides. Solid State Ionics 1999, 125, 285-302. (20) Kreuer, K.D.; Paddison, S. J.; Spohr, E.; Schuster, M. Transport in Proton Conductors for Fuel-Cell Applications: Simulations, Elementary Reactions, and Phenomenology. Chem. Rev. 2004, 104, 4637-4678. (21) Gorbova, E.; Maragou, V.; Medvedev, D.; Demin, A.; Tsiakaras, P. Investigation of the Protonic Conduction in Sm Doped BaCeO3. J. Power Sources 2008, 181, 207–213. (22) Riess, I. Current-Voltage Relation and Charge Distribution in Mixed Ionic Electronic Solid Conductors. J. Phys. Chem. Solids 1986, 47 (2), 129−138. (23) Ross, Jr, P. N.; Benjamin, T. Thermal Efficiency of Solid Electrolyte Fuel Cells with Mixed Conduction. J. Power Sources 1977, 1 (3), 311−321. (24) Yuan, S.; Pal, U. Analytic Solution for Charge Transport and Chemical Potential Variation in Single Layer and Multilayer Devices of Different Mixed Conducting Oxides. J. Electrochem. Soc. 1996, 143, 3214. (25) Riess, I. Mixed Ionic-Electronic Conductors – Material Properties and Applications. Solid State Ionics 2003, 157 (1−4), 1−17. (26) Demin, A.; Gorbova, E.; Tsiakaras, P. High Temperature Electrolyzer Based on Solid Oxide Co-ionic Electrolyte: A Theoretical Model. J. Power Sources 2007, 171, 205–211. (27) Huang, J.; Yuan, J.; Mao, Z.; Sunden, B. Analysis and Modeling of Novel lowTemperature SOFC with a Co-ionic Conducting Ceria Based Composite Electrolyte. J. Fuel Cell Sci. Technol. 2010, 7 (1), 011012.

ACS Paragon Plus Environment

Industrial & Engineering Chemistry Research 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

(28) Bavarian, M.; Kevrekidis, I. G.; Benziger, J. B.; Soroush, M. Modeling and Bifurcation Analysis of a Co-ionic Conducting Solid Oxide Fuel Cell. Ind. Eng. Chem. Res. 2013, 52, 3165−3177. (29) White, R.E.; Lorimer, S.E.; Darby, R. Prediction of the Current Density at an Electrode at which Multiple electrode reactions occur under Potentiostatic Control. J. Electrochem. Soc. 1983, 130 (5), 1123-1126. (30) Yarlagadda, V.; van Nguyen, T. A 1D Mathematical Model of a H2/Br2 Fuel Cell. J. Electrochem. Soc. 2013, 160 (6) F535-F547. (31) Bard, A. J.; Faulkner, L. R. Electrochemical Methods and Applications, 2nd ed.; John Wiley: New York, 2001. (32) Bazant, M. Theory of Chemical Kinetics and Charge Transfer based on Nonequilibrium Thermodynamics. Accounts Chem. Res. 2013, 46 (5), 1144-1160. (33) Gyftopoulos, E.; Beretta, G. Thermodynamics: Foundations and Applications; Dover Publications: Mineola, N.Y, 2005. (34) Uchida, H.; Yoshikawa H.; Iwahara. H. Formation of Protons in SrCeO3-based Proton Conducting Oxides. Part I. Gas Evolution and Absorption in Doped SrCeO3 at high Temperature. Solid State Ionics 1989, 34 (1-2), 103-110. (35) Pegis, M.; Roberts, J.; Wasylenko, D.; Mader, E.; Appel, A.; Mayer, J. Standard Reduction Potentials for Oxygen and Carbon Dioxide Couples in Acetonitrile and N,NDimethylformamide. Inorg. Chem. 2015, 54 (24), 11883-11888. (36) Tsiplakides, D.; Vayenas, C. Electrode Work Function and Absolute Potential Scale in Solid-State Electrochemistry. J. Electrochem. Soc. 2001, 148 (5), E189.

ACS Paragon Plus Environment

Page 40 of 41

Page 41 of 41 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research

For Table of Contents Only

e

H2O, O2

V

H+ O=

CO2, CO, CH4, H2O

H2O → 2H+ + ½ O2 + 2e-

CO2 + 2H+ + 2e- → CO + H2O

O2-→ ½ O2 + 2e-

CO2 + 4H+ + 8e- → CH4 + 2O2-

ACS Paragon Plus Environment