The Starship Engineer's Notebook: Saha Equation and Simultaneous Equilibrium Dissociation-Ionization Reactions

Abstract

Reaction coordinates derived from the Saha equation and stoichiometry are used to describe the dissociation of diatomic molecules and the ionization of dissociation products. This results in a system of 2 coupled nonlinear equations. In previous literature [1] this was resolved by assuming two different, independent regimes, which works provided that dissociation energy and ionization energy are sufficiently different. For substances like nitrogen, however, these two values are close and this assumption may not be valid. For this reason we solve the system numerically and use this result to test the two regime assumption.

Introduction

The Saha equation predicts the concentration of reactants and products in equilibrium based on temperature. In this paper this equation is first discussed in the simple dissociation/ionization reaction:

$AB \leftrightarrow A + B ~~\Delta E$

Where this basic reaction could represent the dissociation of molecule AB or the ionization of an atom A, where B is an electron. For most systems of practical interest only single ionization is likely, so only the first ionization energy is relevant.

Equation 1 becomes the more complex system of dissociation and ionization by combining both possibilities. In this paper we assume a diatomic molecule, so in the first reaction A and B are the same atom and in the second reaction B is an electron:

$A_2 \leftrightarrow 2A ~~~ \Delta E_D \newline A \leftrightarrow A^* + e^- ~~ \Delta E_I$

In prior literature this system was explored only in the context where dissociation and ionization energy are sufficiently removed that this system can be explored as two separate problems, as in the case of hydrogen dissociation and ionization, with a large gap between the two energy values:

$H_2 \leftrightarrow 2H ~~ \Delta E = 4.476eV \newline H \leftrightarrow H^* + e^- ~\Delta E = 13.6eV$

In this paper we will consider the situation where the two energies are close and this assumption may not be valid, as in the case of nitrogen plasma formation:

$N_2 \leftrightarrow 2N ~~ \Delta E = 9.79eV \newline N \leftrightarrow N^* + e^- ~\Delta E = 14.53eV$

Background

Saha Equation

The Saha equation predicts equilibrium concentrations of products and reactants based on temperature [2]. As derived from the partition function it has the form:

${{[n_A][n_B]}\over{[n_{AB}]}} = ({{2 \pi k T}\over{h^2}} {{M_A M_B}\over{M_{AB}}})^{3/2} ({{Z(i)_A Z(i)_B}\over{Z(i)_{AB}}})$

Where Z(i)_A gives the portion of the partition function related exclusively to internal degrees of freedom. By assuming that AB has no internal structure and only a single bound state, and that A and B have no internal structure, this portion reduces to just e^-^βΔ^E. A similar result from statistical mechanics results by assuming the ratio of statistical weights is near unity [2].

From the Saha equation a reaction coordinate, x, is derived as follows:

$ln({{x^2}\over{1-x}}) = \beta \Delta E - {3 \over 2} ln(\beta \Delta E) + {3 \over 2} ln({{2 \pi \Delta E}\over{n_0^{3/2} h^2}} {{M_A M_B}\over{M_{AB}}})$

In the simple dissociation reaction defined by equation 1 the reaction coordinate x refers to the percentage of AB that has been broken down, and thus the amount of A and B formed:

$[n_{AB}] = (1-x)[n_0] \newline ~~~~~ [n_A] = [n_B] = x[n_0]$

In the more complicated set of dissociation and ionization reactions we must also take into account the stoichiometric coefficients.

Stoichiometry

Stoichiometric coefficients describe the creation and destruction of substances in a sequence of chemical reactions. Consider the dissociation/ionization system described by equation 2 above. In this simple case there are four terms involved: A2, A, A+, and e-. In most practical applications we're not concerned with the formation of electrons so we ignore the e- term:

Reaction Coord.	A₂	A	A⁺
x	-1	+2	0
y	0	-1	+1

Table 1: Stoichiometric coefficients [1]

We modify 7 with a second coordinate, y, to account for the second reaction and the stoichiometric coefficients to reflect the increase in particle count:

$[n_{A2}] = (1-x)[n_0] \newline ~~~~~~[n_A] = 2(x-y)[n_0] \newline ~~~~~~[n_{A+}] = 2y[n_0]$

Dissociation-Ionization Reaction

For the specific dissociation-ionization reaction we must modify the Saha equation and derived reaction coordinate to account for the additional terms. It is necessary to account for the formation of 2 particles of A for every dissociation of A₂, and to cancel out similar mass terms:

${[n_A]^2 \over [n_{A2}]} = {1 \over e^{\beta E_D}} ({{2 \pi k T}\over h^2}{M_A \over 2})^{3/2} \newline ~~~{[n_{A+}] \over [n_A]} = {{1}\over{e^{\beta E_I}}} ({ {2 \pi k T}\over h^2} M_e)^{3/2}$

We've used the simplifying assumption here that A₂ is simply twice the mass of a single A and that the change in mass from losing a single electron is negligible.

Combine these simplified Saha equations with the coefficients from table 1 to find the reaction coordinates, x and y:

$ln({{[2(x-y)]^2} \over {1-x}}) = \beta E_D - {3 \over 2} ln(\beta E_D) + {3 \over 2} ln({{2 \pi E_D}\over{n_0^{3/2} h^2}} {{M_A}\over{2}}) \newline ~~~~~~~ ln({{(2y)^2}\over{2(x-y)}}) = \beta E_I - {3 \over 2} ln(\beta E_I) + {3 \over 2} ln({{2 \pi E_I}\over{n_0^{3/2} h^2}} M_e)$

Method

The coupled nonlinear system given in equation 10 is analytically solvable, however the solution is unwieldy and, in general, the resulting system need not be solvable. In previous literature [1] this was resolved by assuming two different, independent regimes. In the dissociation regime y is uniformly 0 and in the ionization regime x is uniformly 1. This assumption works provided that dissociation energy and ionization energy are sufficiently different. For substances like nitrogen, however, these two values are close and this assumption may not be valid. Therefore we will attempt to solve the system numerically. This machinery can also be used to test the validity of the two regime assumption for various substances.

The system of equations is very intimidating to look at. However it is possible to reduce them into a relatively simple form by noting that the right hand side of both equations is independent of x and y:

$ln({[2(x-y)]^2 \over(1-x)}) = f_D(\beta) \rightarrow {[2(x-y)]^2 \over (1-x)} = D \newline ~~~~~~~~~ln({(2y)^2 \over {2(x-y)}}) = f_I (\beta) \rightarrow {(2y)^2 \over {2(x-y)}} = I$

Although still coupled nonlinear equations, the set 11 at least looks a bit more manageable. This is a system of 2 nonlinear equations and 2 unknowns, and is solved using Newton's method for nonlinear systems, where f is the vector of functions f(x)=0 and F is its Jacobian[3]:

$\vec x_{n+1} = \vec x_n - F^{-1}(\vec x_n)f(\vec x_n) \newline f=\begin{bmatrix} [2(x-y)]^2 - D(1-x)\\ (2y)^2 - 2I(x-y) \end{bmatrix} = 0 \newline F = \begin{bmatrix} D+8(x-y) & 8(y-x)\\ -2I & 2(I-4y) \end{bmatrix}$

This can then be iteratively solved until x and y converge. Once x and y are found they give the values of concentration for the 3 products of interest at a given temperature. In the case of nitrogen:

Dissociation Energy	Ionization Energy	Atomic Mass
9.79 eV	14.53 eV	14 amu

Table 2: Nitrogen Properties [2]

It is worth noting that Newton's method is extremely sensitive to the initial guess of x₀ and y₀ so care must be taken with the initial estimate to avoid the system diverging or converging to nonphysical stable points[3]. Finding the inverse of the Jacobian analytically can also help avoid oddly singular or non-invertable values during computation.

Results

The process described above was performed in a range of interest specifically 1,500 – 15,000 K, where it is expected that any overlap between ionization and dissociation reactions will occur. Figure 1 shows the results in percent concentration vs temperature in the initial range of interest. Brief observation of the results will show that, in fact, there is trivial overlap between 8,000 – 9,000 K. The assumption of two distinct regimes would result in minimal loss of accuracy. Nitrogen was specifically selected since its dissociation and ionization energies are among the closest of any common gas. This indicates that the two regime model is applicable to any diatomic gas of practical interest.

.Figure 2 shows the reaction coefficients, x and y, as functions of temperature. Again it can be seen clearly that there is trivial overlap between the two. At no point where y significantly exceeds 0 is x also significantly less than 1. The divide between the two regions seems to be around 9,000 – 10,000 K for nitrogen. A similar analysis would likely show the threshold for other gasses which, it's expected, would have a much larger gap.

Conclusions

Using this framework the dissociation-ionization equilibrium should be possible to calculate for any diatomic molecule. Since standard atmosphere is almost entirely N₂ and O₂ for atmospheric plasma study the obvious next step is to perform the same technique for oxygen.

Once the breakdown of both oxygen and nitrogen is understood it becomes possible to describe the formation and decomposition of the products of N + O reactions. This equilibrium also follows the Saha equation.

In non-equilibrium plasma formation it is also possible to decouple electron, vibrational, and translational temperatures. Since the decomposition reaction depends on vibration, recombination depends on translation, and ionization depends on electron temperature then all of the processes under study can be described in terms of their individual temperatures rather than global temperature as was done here. Since thermal non-equilibrium is rarely associated with chemical equilibrium use of the Saha equation here would tend to overpredict production. However this analysis could still provide valuable predictions of NO_x formation in atmosphere.

References

2: Robert Gilmore, Saha's Equation and the Partition Function,

2: Alexander Fridman, Plasma Chemistry, 2008

3: Kendall Atkinson, An Introduction to Numerical Analysis, 1989

Appendix 1: Table of Symbols

[n_a] = concentration of species “a” in #/cm³
M_a = mass of species “a” in amu
k = Boltzmann constant 8.617*10^-5 eV/K
h = Planck constant 4.135*10^-15eV-s
E_I, E_D = ionization/dissociation energy
I, D = reaction coordinate energy dependence
β = inverse temperature
x, y = reaction coordinates
f = function vector
F = Jacobian

Appendix 2: SCILAB Code

clear

//universal and chemical constants, mks units
k=1.38*10^-23, h=6.626*10^-34, Mp=1.66*10^-27;
ED=9.79*1.602*10^-19, EI=14.53*1.602*10^-19, MA=14*Mp, Me=9.11*10^-31; 

uD=MA/2, uI=Me;
tol=10^-6, t=1500, dt=10, n0=101325/(k*t), x=0, i=1;

x=0.5, x0=2, y=0.25, y0=2;

while t<15000
    //solve for the temperature and energy dependence
    B=1/(k*t);
    D=exp(-B*ED - 3/2*log(B*ED) + 3/2*log((2*%pi*uD*ED)/(n0^(2/3)*h^2)));
    I=exp(-B*EI - 3/2*log(B*EI) + 3/2*log((2*%pi*uI*EI)/(n0^(2/3)*h^2)));
    
    //solve x by Newton's method
    while abs(x-x0)+abs(y-y0)>tol
        f=[(2*(x-y))^2-D*(1-x);(2*y)^2-2*I*(x-y)];
        a=D*(I+4*y) + 32*y*(x-y);
        F=1/a*[I+4*y , 4*(x-y); I, D/2+4*(x-y)];
        x0=x, y0=y;
        X=[x0;y0] - F*f;
        x=X(1), y=X(2);
    end
    
    Z(i,1:2)=[x,y];
    
    //find concentrations
    na2= (1-x)*n0;
    na = 2*(x-y)*n0;
    nap= 2*y*n0;
    
    //store results and increment
    N(i,1:3)=[na2,na,nap]/(na2+na+nap)*100;
    T(i)=t;
    i=i+1, t=t+dt;
end

clf()
plot2d(T,N)
xtitle("Nitrogen Dissociation and Ionization vs. Temperature","Temperature (K)","Percent Substance")
legend("N2","N","N+",-1)

The Starship Engineer's Notebook

Dropdown Navigation

March 27, 2016

Saha Equation and Simultaneous Equilibrium Dissociation-Ionization Reactions