© 2022 Horia Metiu

WorkBook T18 . Equilibrium

WorkBook T18.1 The equilibrium constant and the equilibrium composition

The extent of reaction, the initial number of moles, molar fractions and equilibrium constant are related to each other.  In this workbook I examine how to calculate all these quantities if I know one of them.

Cell 1:  Explain the procedure

The reaction is 
	(Subscript[N, 2])(g) + 3 (Subscript[H, 2])(g)->2(Subscript[NH, 3])
with the initial composition 
	(n^0)((Subscript[N, 2]))= 1 ,   (n^0)((Subscript[H, 2])) = 3   and   (n^0)((Subscript[NH, 3]))=0.1
In all that follows I assume that the initial number of moles is the same.  
I will examine the following situations
• Know ξ, calculate everything else
• Know one equilibrium concentration (here x((Subscript[NH, 3])) ) , calculate everything else
• Know equilibrium constant K, calculate everything else. 
The equation used in these calculations are 
	n(i)=(n^0)(i)+ξ ν(i) 									(1)
	x(i) = n(i)/(∑_(i=1)^c)n(i)  									(2)
	 K=((x((Subscript[N, 2]))((x((Subscript[H, 2])))^3))/((x((Subscript[NH, 3])))^2))										(3)

Cell 2.  Input the quantities you need

Remove["Global`*"] ;
Needs["Utilities`CleanSlate`"];
CleanSlate[];

Make a list with compound names, one with stoichiometric coefficients and one with the initial number of moles.

compounds={(Subscript[N, 2]),(Subscript[H, 2]),(Subscript[NH, 3])};
stoCoefs={-1,-3,2};
nInitial={1,3,0.1};

Cell 3. I know  (Subscript[n, i]^0), i = 1, 2, 3 and ξ:  calculate (Subscript[n, i]) , (Subscript[x, i]) and K

Calculate the number of mole as a function of ξ using
	n(i)=(n^0)(i)+ξ ν(i)

nbMoles[stoCoefs_,nInitial_,ξ_  ]:=Module[{},nInitial+stoCoefs ξ]

n=nbMoles[stoCoefs,nInitial,ξ]

Calculate the molar fractions
	x(i) = n(i) / (∑_(i=1)^c)n(i)

molarFractions[stoCoefs_,nInitial_,ξ_]:=Module[{nbMoles},
 nbMoles=Plus@@(nInitial+stoCoefs*ξ);
		(nInitial+stoCoefs*ξ)/nbMoles]

x=molarFractions[stoCoefs,nInitial, ξ]

The molar fraction of (Subscript[NH, 3]) is the third in the list because it is the third in the list compounds

compounds

xNH3=mf[[3]]

toTable=MapThread[{#1,#2,#3,#4 }&,{compounds, nInitial,n,x }]

TableForm[toTable,TableHeadings->{{},{"comp", "(Subscript[n, 0])","n(ξ)","x(ξ)"}} ]

equilibriumK[stoCoefs_,initialNbMoles_,ξ_]:=Module[{},
Times@@(molarFractions[stoCoefs,initialNbMoles,ξ]^stoCoefs)]

K = equilibriumK[stoCoefs,nInitial,ξ]

Once I know that ξ = 0.88 and the initial number of moles  I can calculate all molar fractions and the equilibrium constant.

molarFractions[stoCoefs,nInitial,0.88]

nbMoles[stoCoefs,nInitial,0.88]

compounds

equilibriumK[stoCoefs,nInitial,0.88]

Cell 4. I know (Subscript[n, i]^0), i = 1, 2, 3  and measured x((Subscript[NH, 3])) : calculate (Subscript[n, i]), (Subscript[x, i]) and K

I apply this equation to ammonia
	(Subscript[x, i]) = (((Subscript[n, i]^0) + (Subscript[ν, i]) ξ)/((∑_( i=1)^( c))(Subscript[n, i]^0) + (Subscript[ν, i]) ξ))
I know x((Subscript[NH, 3]))  so I calculate ξ from this equation. Once I have ξ I can calculate the other molar fractions, the number of moles and the equilibrium constant.

nInitial

nInitial

nInitial={1,3,0}

nInitial[[3]] + stoCoefs[[3]] ξ

Plus@@(nInitial + ξ stoCoefs)

(nInitial[[3]] + stoCoefs[[3]] ξ)/(Plus@@(nInitial + ξ stoCoefs))

Clear[xNH3]
eq = xNH3==(nInitial[[3]] + stoCoefs[[3]] ξ)/(Plus@@(nInitial + ξ stoCoefs))

ruleξ=Solve[eq,ξ]//Flatten//Expand

eq/.ruleξ//FullSimplify//Chop

molarFractions[stoCoefs, nInitial,ξ]

molarFractions[stoCoefs, nInitial,ξ]/.ruleξ//FullSimplify

Plot[molarFractions[stoCoefs, nInitial,ξ]/.ruleξ//FullSimplify,{xNH3,0,1}]

nbMoles[stoCoefs,nInitial,ξ]/.ruleξ

Cell 5. I know  (Subscript[n, i]^0), i = 1, 2, 3 and K = 1: calculate (Subscript[n, i]), (Subscript[x, i]) and ξ

• I have
	n(i)=(n^0)(i)+ξ ν(i) 									(1)
	x(i) = n(i)/(∑_(i=1)^c)n(i) = (((n^0)(i)+ξ ν(i) )/((∑_(i=1)^c)(n^0)(i) + ξ ν(i) ))  						(2)
	 K= (Subscript[Π, i])((x(i))^(ν(i)))											(3)
• When I insert x(i) given by Eq. (2) in Eq. (3) I get a formula connecting K to the extent of reaction.  
• Solve that to get ξ
• Use ξ in Eq. (1) to get n(i)
• Use ξ in Eq. (2) to get the molar fractions x(i)

I defined in the earlier cell all the functions I need
	nbMoles[stoCoefs, nInitial, ξ]   gives n(i) = (n^0)(i)+ξ ν(i)
	molarFractions[stoCoefs, nInitial, ξ]  gives x(i) = (((n^0)(i)+ξ ν(i) )/((∑_(i=1)^c)(n^0)(i) + ξ ν(i) ))
	equilibriumK[stoCoefs, nInitial, ξ] gives K= (Subscript[Π, i])((x(i))^(ν(i)))

For the reaction 
	(Subscript[N, 2]) + 3 (Subscript[H, 2]) = 2 (Subscript[NH, 3])
• The initial number of moles are 
	(Subscript[n, 0]) ((Subscript[N, 2])) = 1
	(Subscript[n, 0]) ((Subscript[H, 2])) = 3
	(Subscript[n, 0])((Subscript[NH, 3])) = 0.1
•  The equilibrium constant is
	K = ((x((Subscript[N, 2]))((x((Subscript[H, 2])))^3))/((x((Subscript[NH, 3])))^2))

Input the information you have

compounds={(Subscript[N, 2]),(Subscript[H, 2]),(Subscript[NH, 3])};
stoCoefs={-1,-3,2};
nInitial={1.0,3.0,0.1};
K = 1.0 ;

The function reaction will verify that my input is OK.  See the previous problem for details.

Solve the equation 
	K =(((((4.1` - 2 ξ))^2) (((0.1` + 2 ξ))^2))/((((3 - 3 ξ))^3) (1 - ξ))) = 1
to get ξ

Clear[K,ξ]
K=1;
eq=K==equilibriumK[stoCoefs,nInitial, ξ]

The equation is fourth order and will have four solutions. I eliminate negative or complex solutions and I get two real solutions.

ruleξ=Solve[{eq,ξ>0},ξ,Reals]//Flatten

Only one represents my system. Which one? In this example I can see that 
	ξ = 1.692
is not possible because I have initially only one mole of (Subscript[N, 2]) so ξ ought to be less than 1.  However, I want to rule out one of the ξ-s based on some computational rules.  For the correct ξ I must have
• All n(i) are positive
• All x(i) must satisfy 0 ≤ x(i) ≤ 1

Try the  first solution:

nbMoles[stoCoefs,nInitial, ξ]/.ruleξ[[1]]

Try the second solution:

nbMoles[stoCoefs,nInitial, ξ]/.ruleξ[[2]]

The second solution is unphysical.  Check that the first solution gives reasonable results for x(i)

molarFractions[stoCoefs,nInitial, ξ]/.ruleξ[[1]]

Plus@@%

molarFractions[stoCoefs,nInitial, ξ]/.ruleξ[[2]]

compounds

Obviously the first solution gives acceptable results. Make a table with the results

toTable=MapThread[{#1,#2,#3 ,#4}&,{compounds,nInitial,nbMoles[stoCoefs,nInitial, ξ]/.ruleξ[[1]],molarFractions[stoCoefs,nInitial, ξ]/.ruleξ[[1]]}]

TableForm[toTable,TableHeadings->{None,{ "comp i","(Subscript[n, 0])(i)","n(i)","x(i)"}}]

WorkBookT.18.2. The equilibrium properties of the reaction (Subscript[F, 2])(g) = 2 F(g)

I will study this reaction to answer two questions
1. The dependence of the equilibrium composition of the equilibrium constant.
2. The dependence of equilibrium composition on initial number of moles.

Cell 1. Make input. Define the list of compounds, of stoichiometric coefficients and of initial number of moles.

I will use two initial number of moles to see how changing them affects equilibrium molar fractions.

Remove["Global`*"] ;
Needs["Utilities`CleanSlate`"];
CleanSlate[];

compounds={"(Subscript[F, 2])","F"};
ν={-1,2};
n01={10,1};
n02={5,1};

Cell 2.  Set up the equation connecting the equilibrium constant to ξ.

Use 
	n(i)=(n^0)(i)+ ξ ν(i)
and 
	x(i)=n(i)/(∑_(j=1)^c)n(j)
where c=2 is the number of compounds participating in the reaction, to calculate how  the number of moles n(i) and the molar fractions x(i) depend on the extent of reaction ξ.  The number of moles in the case 1 is stored in the list n1, the number of moles in the case 2 is in n2, the molar fractions for case 1 are in x1 and the molar fractions for case 2 are in x2.

n1=n01 +ν ξ  ;
n2=n02+ν ξ ;
x1=n1/(Plus@@n1) ;
x2=n2/(Plus@@n2)  ;

The expressions for the equilibrium constants for the two cases are (use K=(((x(F))^2)/(x((Subscript[F, 2]))))):

Clear[K1]
K1[y_]:=Times@@(x1^ν)/.ξ->y
Clear[K2];
K2[y_]:=Times@@(x2^ν)/.ξ->y

Display the results:

For case 1 ( (n^0)((Subscript[F, 2])) =10, (n^0)(F)=1 ) I have:

Print["the number of moles and the molar fractions in case 1"];TableForm[MapThread[{#1,#2}&,{n1,x1}],TableHeadings->{compounds,{"n ","x"}}]

WorkBookT.18.3. The equilibrium composition of the reaction  (Subscript[CH, 4])+(Subscript[H, 2])O  =  CO + 3 (Subscript[H, 2])

The equilibrium constant, at a temperature of 1273 K and pressure of 1 bar, has been calculated in two ways.  An approximate calculations gives K = 26642.4 and a rather precise one gives K = 3200.16.   Calculate the equilibrium composition corresponding to the two values of the equilibrium constant, for the case when the initial number of moles are 
	1 mole (Subscript[CH, 4]), 
	2 mole (Subscript[H, 2])O, 
	0.1 mole (Subscript[H, 2]) 
	0.1mole of CO. 
Do an error analysis: given the large error in K, how large will the error in ξ be?

Plan:  
1.	The equilibrium constant is:
		K=((x(CO) ((x((Subscript[H, 2])))^3))/(x((Subscript[CH, 4])) x((Subscript[H, 2])O))) 						(2)
 2.	 I express the molar fractions in terms of the extent of reaction. 
		x(i) = [ (n^0)(i) + ξ ν(i) ]/ (∑_(1=1)^c)[ (n^0)(i) + ξ ν(i) ]			(3)
	By inserting these expressions in the formula (2) for the equilibrium constant I obtain an equation from which I calculate ξ. 
3.	I obtain multiple solutions and I chose the one that is physically meaningful.   
4.	Knowing ξ I can calculate the number of moles at equilibrium from
		n(i)=(n^0)(i)+ν(i) ξ						(4)
	and the molar fractions from Eq. (3).
I perform these calculations for the approximate and the precise values of K to see how the error made in calculating K affects the calculated equilibrium composition.

Cell 1. Input data:

compounds={"CH4","H2O","H2","CO"};

Stoichiometric coefficients (the order in this list is the same as in the list compounds: the first value -1 corresponds to (Subscript[CH, 4]),etc.)

nu={-1,-1,3,1};

The initial number of moles.

n0={1.,2.,0.1,0.1};

The approximate and the "exact" equilibrium constants Ka and Ke, respectively, are:

Ka=26642.4;
Ke=3200.4;

Cell 2.  Calculate the dependence of the number of moles and the molar fractions on ξ and the initial molar fractions.

The number of moles as a function of the extent of reaction and the initial number of moles (use n(i)=(n^0)(i)+ξ ν(i)):

n=n0+ ξ nu

The molar fractions as functions of initial number of moles and of the extent of reaction (use Eq. 3)

x=(n0+ ξ nu)/(Plus@@(n0+ ξ nu))

The initial molar fractions (not needed for calculations)

xInit=x/.ξ->0

Cell 3: Display the results:

toTable=MapThread[{#1,#2,#3,#4}&,{n0,n,x,xInit}];
TableForm[toTable,TableHeadings->{compounds,{"(n^0)(i)","n(i)","x(i)","(Subscript[x, init])"}}]

Cell 4: Define the equation for the equilibrium constant in terms of the extent of reaction.  Calculate the value of ξ for the two equilibrium constants:

The expression of the equilibrium constant, as a function of ξ, is (use Eq.2 with the mole fractions given by Eq. (3) ). The mole fractions are stored in the list x and the stoichiometric coefficients in the list nu:

eqConst=Times@@( x^nu)

The two equations to be solved for ξ are:

eqa=Ka==eqConst
eqe=Ke==eqConst

Solve the equation eqa for the approximate value of K

sola=Solve[eqa,ξ]

Solve the equation eqe for the "exact" value of K

sole=Solve[eqe,ξ]

Calculate the molar fractions corresponding to the four solutions for ξ, for the approximate equilibrium constant

x/.sola

Only the third set of molar fractions are physically acceptable. Therefore, the third value of ξ in the list sole, which is (Subscript[ξ, a]) =0.999955 is the correct one. The molar fractions corresponding to that value of ξ are:

xa=(x/.sola)[[3]]

Do the same calculations with the solution of the equation eqe, which uses the "exact" equilibrium constant:

x/.sole

Only the third set of molar fractions are physically acceptable (must be positive, real and between zero and one). The extent of reaction for the equilibrium constant (Subscript[K, e])is therefore (Subscript[ξ, e]) =0.999622.  The corresponding molar fractions are:

xe=(x/.sole)[[3]]

Also calculate the number of moles of each compound at equilibrium, for (Subscript[ξ, a])and (Subscript[ξ, e]):

ne=(n/.sole)[[3]]

na=(n/.sola)[[3]]

Note that the number of moles of methane is close to zero. This "coincidence" will play an important role in what follows.

Display the results:

toTable=MapThread[{#1,#2,#3,#4}&,{xa,xe,na,ne}];TableForm[toTable,TableHeadings->{compounds,{"(Subscript[x, a])(i)","(Subscript[x, e])(i)","(Subscript[n, a])","(Subscript[n, e])"}}]

This appears very strange: even though the numerical values of the equilibrium constants (Subscript[K, a]) and (Subscript[K, e]) are very different, the corresponding extents of reactions ((Subscript[ξ, a]) =0.9999545202701221` and (Subscript[ξ, e]) = 0.999622) are very close to each other.  They differ at the fourth decimal! The molar fractions derived from them are also practically the same. We can easily understand why this happens: in both cases the equilibrium constant are very large and the reaction goes practically to completion (until it runs out to methane). This makes the equilibrium extent of reaction very close to 1 (when ξ is 1 there is no methane left, as you can see from n((Subscript[CH, 4]))=(n^0)((Subscript[CH, 4])) -ξ and the fact that (n^0)((Subscript[CH, 4]))=1)

WorkBookT.18.4. Determine the direction of the reaction C+2(Subscript[H, 2])O = (Subscript[CO, 2]) +2(Subscript[H, 2])

The reaction is performed with 
	5 moles of graphite, 
	three moles of water vapor, 
	three moles of (Subscript[CO, 2]), 
	two moles of (Subscript[H, 2]). 
Two experiments are performed:  in one the compounds are mixed at 
	T = 500K and p = 1bar 
and in the other at 
	T = 1000K and p = 1bar.  
I want to know whether the reaction will form more (Subscript[H, 2]) or more water (i.e. I want to know the direction of the reaction).  I calculated that  
	 ((ΔG)^0) = 43.227 kJ at T = 500 K and p = 1 bar 
and 
	((ΔG)^0) = -10.402 kJ, at T = 1000K and p = 1 bar 
	(this calculation is performed in Chapter 19).

To determine the direction of the reaction I have to calculate the sign of chemical affinity A.  I use 
	A=((ΔG)^0) +R T ln K							(1)
with 
	K=((((x((Subscript[H, 2])))^2)x((Subscript[CO, 2])))/((x((Subscript[H, 2])O))^2)) 			          				(2)
	
The molar fraction of graphite does not appear in Eq. 2 because the graphite is solid and does not form a solution with the gases. 
The molar fraction are calculated from their definition
	x(i)=((n(i))/((∑_(j=1)^c)n(i)))							(3)
where the sum runs over the number of moles of gases.

Plan for solution:
1. Make input: define the arrays compounds, ⋁ and n0 containing the names of the gaseous compounds, the stoichiometric coefficients of the gases and their intial number of moles, respectively.  The lists are ordered: the second element in ν and n0 correspond to the second compound in the list compounds.
2. Calculate the molar fractions (use Eq. 3).
3. Calculate the equilibrium constant K (use  Eq. 2).
4. Calculate the affinity A (use Eq. 1).  If A < 0 the reaction goes from left to right. If A > 0 it goes from right to left. If A=0 nothing happens; the initial composition is the equilibrium composition.

Cell 1:  Make input: define the arrays: compounds, ⋁ and n0 containing the names of the gaseous compounds, the stoichiometric coefficients of the gases and their intial number of moles, respectively.  The lists are ordered: the second element in ν and n0 correspond to the second compound in the list compounds.

compounds={"(Subscript[H, 2])O","(Subscript[CO, 2])","(Subscript[H, 2])"}; (*gases only*)
ν={-2,1,2};(*stoichiometric coefficients*)
n0={3,3,2};(*initial number of moles*)

Cell2: Calculate the molar fractions (use Eq. 3)

x=n0/Plus@@n0//N

Cell 3: Calculate the equilibrium constant at the initial composition (note that the system is not at equilibrium and the name equilibrium constant is misleading;   I'll use it to denote the right hand side of Eq. 2).

K=Times@@(x^ν)

Cell 4: Calculate the affinity for the two values of ((ΔG)^0) and the two temperatures (use 
	A=((ΔG)^0) + R T ln[K] )

ΔG0={43.227,-10.402} (*kJ*);
T={500,1000} (*K*);
A=MapThread[(#1+R*#2*Log[K])&, {ΔG0,T}]/.R->8.134*10^-3  (*use R in kJ/molK*)

Display the results in a table

toTable=Transpose[{T,ΔG0,A}];
TableForm[toTable, TableHeadings->{None,{"T K","((ΔG)^0)kJ","A kJ"}}]

At 500K and p = 1 bar the affinity is positive so the reaction goes from right to left;  (Subscript[H, 2]) and (Subscript[CO, 2]) will react to make graphite.  Whether you see this happen in an experiment depends on the rate of the reaction.  Thermodynamics has nothing to say about the rate. If the rate is low it may take millennia for the system to make more graphite than what was initially present. This is what happens when the experiment is done. At T = 1000 and p = 1 bar the affinity is negative. The reaction proceeds to the right and will produce (Subscript[H, 2]) and (Subscript[CO, 2]).  Note that the change of sign occurred because ((ΔG)^0)is negative at 1000 K ( RT ln(K) is negative because K < 1).