WorkBookSM. 11 Chemical Equilibrium for the dissociation of a homonuclear diatomic

We study the equilibrium properties of the reaction:

	(Subscript[A, 2])  ⇌  2 A						(1)
	
The equilibrium constant is the product of a translational contribution (Subscript[K, t]), a vibrational contribution (Subscript[K, v]), a rotational contribution (Subscript[K, r]) and an electronic contribution (Subscript[K, e]): 

	K =(Subscript[K, t]) (Subscript[K, v]) (Subscript[K, r]) (Subscript[K, e]) 					(2)
	
The translational contribution is:

	(Subscript[K, t])= ((Z(A))^2)/Z((Subscript[A, 2]))					(3)
with
	Z(A) = 
        (RT/(((Λ(A))^3) (Subscript[N, A]) p))					(4)
	
	Z((Subscript[A, 2])) = 
        (RT/(((Λ((Subscript[A, 2])))^3) (Subscript[N, A]) p))					(5)
	
Here \(TraditionalForm\`N\_A = \(\(6.022\)\(\ \)\(×\)\(\ \)\(10\^23\)\(\ \)\)\)is  Avogadro's number,  R is the gas constant, p is the pressure, T is temperature, and Λ is the de Broglie wave length (see Eq. 6).  The partition functions Z(A) and Z((Subscript[A, 2])) are dimensionless quantities.  Since pressure is usually given in atmosphere, I will use R = 82 atm (cm^3)/mol K and calculate Λ in (cm^3).  The formula for Λ is:

	Λ(
          X) =  (Sqrt[((h^2)/(2  π m(X) (Subscript[m, 0]) (Subscript[k, B]) T))]) 				(6)
	
where  m(X)  is the mass of compound X in atomic units and (Subscript[m, 0]) = 1.6605 (10^(-24))  converts atomic mass units into grams. The old Planck constant is h=6.626247224951591`*^-27erg sec and the Boltzmann constant (Subscript[k, B])= 1.3806*(10^(-16))erg/K. Using these values in Eq. 6 gives Λ in cm.
	The vibrational contribution: 
	
	(Subscript[K, v])=1/(Subscript[q, v  ])						(7)
	
where (Subscript[q, v])is the vibrational partition function.

	(Subscript[q, v ])=exp[(Subscript[D, 0])/(Subscript[k, B])  T]/(1 - exp[-(Subscript[T, v])/T])		(8)
	
The vibrational temperature is 

	(Subscript[T, v])  =  ℏω/(Subscript[k, B])						(9)
	
To calculate it I need the vibrational energy ℏω which is provided by spectroscopy. It is likely to be given in (cm^(-1)).  The conversion to erg is 1 (cm^(-1))=1.9862×(10^(-16))erg
	The rotational contribution is:
	
	(Subscript[K, r])  =  (1/(Subscript[q, r])) = σ((Subscript[T, r])/T)					(10)
	
where σ is the symmetry factor (equal to two, for a homonuclear diatomic and to one for a heteronuclear diatomic) and (Subscript[T, r]) is the rotational temperature

	(Subscript[T, r])=B/(Subscript[k, B])						(11)
	
where B is the rotational constant. B is given in units of (cm^(-1))

I will use these equations to calculate the equilibrium constant for three reaction:

	(Subscript[H, 2])⇌2 H
	(Subscript[I, 2])⇌2 I
and 
	(Subscript[K, 2])  ⇌ 2 K
By comparing these reactions to each other we can better understand how the molecular parameters of the molecules control chemical equilibrium.

WorkBookSM.11.1. The data

Off[Remove::rmnsm]
Remove["Global`*"]
Off[General::spell, General::spell1]

The data for the diatomic molecules is given through replacement rules. The replacement list contains   
	m = mass of the atom A in (Subscript[A, 2])(e.g. 1 for (Subscript[H, 2]))
	ϵv = the vibrational energy ℏω in units of (cm^(-1))
	D0 = the dissociation energy in eV
	B = the rotational constant in (cm^(-1))

dataO2 -> {m -> 15.994, ϵv -> 
            1580.19   ,  
          D0 -> 5.115   , 
          B -> 
            1.44563   };
dataH2 = {m -> 1.0080  , ϵv -> 
            4401.21   , D0 -> 4.4781   , 
          B -> 60.853   };
dataI2 = { 
          m -> 126.9045   , ϵv -> 214.50  , 
          D0 -> 1.54238   , 
          B -> 0.0373  };
dataK2 = {m -> 39.098,  ϵv -> 92.3,  D0 -> 0.6849, 
          B -> 0.0561};

The only data I need for atoms are the masses, except in cases when I want to take into account the excited electronic states (for O and perhaps I).
To keep life simple, I will not include the effect of electronic state in the calculations performed here. As an exercise you can include electronic contributions when you study (Subscript[O, 2]) dissociation.

dataI = {m -> 126.9045};
dataH = {m -> 1.008};
dataK = {m -> 39.098};
dataO = {m -> 15.994};

WorkBookSM.11.2.   The translational contribution

I will use the equations (3)-(6)

Cell 1: Calculate the de Broglie wave lengths Λ, from Eq. 6.

The de Broglie wave length for an atom:

\(\(\(Clear[Λa]

    Λa[T_, data_List] := 
      With[{kB = 
            1.3806\ 10^\((\(-16\))\)\  (*erg/
              K*) , 
			\ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ m0 = 
            1.6605\ 10^\((\(-24\))\), 
							h = 
            2\ Pi\ 1.0546\ 10^\((\(-27\))\)}, 
		Sqrt[
            h^2/\((\ 2\ Pi\ \ m\ m0\ \ kB\ T)\)]\  /. 
          data\  (*gives\ Λ\ in\ cm\^\(-1\)*) ]\)\(\ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \)\)\)

Look at the mass and temperature dependence

Λa[T, {}]

If you put in above formula the mass in amu and T in degrees Kelvin, you get the de Broglie wave length in cm

Apply this function to hydrogen and iodine

Λa[T, dataI]

Λa[T, dataH]

Next, calculate the de Broglie wave length for a diatomic molecule.  I give as data the mass of the diatomic. The diatomic is homonuclear. The variable data is a list of replacement rules. The data contains the mass of the atom (in the homonuclear diatomic molecule)

\(\(\(Clear[Λd]

    Λd[T_, data_List] := 
      With[{kB = 
            1.3806\ 10^\((\(-16\))\)\  (*erg/
              K*) , 
			\ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ m0 = 
            1.6605\ 10^\((\(-24\))\), 
							h = 
            2\ Pi\ 1.0546\ 10^\((\(-27\))\)}, 
		Sqrt[
            h^2/\((\ 2\ Pi\ 2\ \ m\ m0\ \ kB\ \ T)\)]\  /. 
          data\  (*gives\ Λ\ in\ cm\^\(-1\)*) ]\)\(\ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \)\)\)

Look at the mass and T dependence

Λd[T, {}]

If you put in above formula the mass of the atom in the homonuclear diatomic molecule in amu, and T in degrees Kelvin, you get the de Broglie wave length in cm

Λd[T, dataH2]

Cell 2:  Calculate the partition functions Z(A) and Z((Subscript[A, 2])) from Eqs. (4) and (5)

The partition function Za for an atom:

Clear[Za]
\(Za[T_, p_ (*atm*) , 
        dataAtom_List] := 
Module[{\ 
          R = 82.0578\  (*atm\ cm^3/
                K\ mol*) , 
\ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ NA = 6.022\ 10^23\ }, 
			R\ T/\((NA\ Λa[T, 
                    dataAtom]^3\ p)\)]\)

The function Zd for a diatomic molecule:

Clear[Zd]
\(Zd[T_, p_ (*atm*) , 
        dataDiatomic_List] := 
Module[{\ 
          R = 82.0578\  (*atm\ cm^3/
                K\ mol*) , 
\ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ NA = 6.022\ 10^23\ }, 
			R\ T/\((NA\ Λd[T, 
                    dataDiatomic]^3\ p)\)]\)

Look at the dependence on T, p and the mass m

Zd[T, p, {}] // PowerExpand

If you put in this formula p in atm, T in degrees Kelvin and the mass m ( of the atom A in (Subscript[A, 2]) ) amu, you get a dimensionless  Z for a diatomic

Zd[298, 1, dataH2]

As a test I check if the function Za for an atom becomes identical to that for a diatomic if I replace in the atom's Z the mass m with 2 m

Za[T, p, {m -> 2 m}] // PowerExpand

The test comes out alright.

Cell 3:  Calculate the translational contribution to the equilibrium constant for a reaction of the type

	(Subscript[A, 2])  ⇌ 2 A
	
Use Eq. (3)

Clear[Kt]
\(Kt[T_, p_ (*in\ Atm*) , dataAtom_List, dataDiatom_List] := 
	
      Module[{\ \ }, 
			\ Za[T, p, dataAtom]^2/
          Zd[T, p, dataDiatom]]\)

Kt[T, p, {}, {}] // PowerExpand

Look at the value at room temperature

Kt[300, 1, dataH, dataH2]

Cell 4: Make a table with  the translational contributions to the equilibrium constant for the dissociation of (Subscript[H, 2]),  (Subscript[I, 2]),   and (Subscript[K, 2])

\(TableForm[
      Table[{T, Kt[T, 1, dataH, dataH2], Kt[T, 1, dataI, dataI2\ ], 
          Kt[T, 1, dataK, dataK2\ ]}, {T, 100, 2100, 200}], 
	
      TableHeadings -> {None, {"T K", \*""(Subscript[K, t])  (Subscript[H, 2])"", \*""(Subscript[K, t]) (Subscript[I, 2])"", \ \*""(Subscript[K, t]) (Subscript[K, 2])""}}]\)

The translational contribution to the equilibrium constant is very large. If translation was the only kind of motion that mattered, than the system at equilibrium would have a large concentration of atoms and very few diatomics. The difference between (Subscript[H, 2]) and (Subscript[I, 2])is rather large even though in (Subscript[K, t]) the only parameter differentiating the molecules is the mass, which is larger for (Subscript[I, 2]). (Subscript[K, t ])goes like (m^(3/2))so the larger the mass the larger the propensity to dissociate

WorkBookSM.11.3. The vibrational  contribution

Cell 1:  Calculate the vibrational temperatures for (Subscript[H, 2]),  (Subscript[I, 2]),  (Subscript[K, 2])

I use the equation
	(Subscript[T, v])  = ℏ(ω/k)
k is Boltzmann constant, ω is the vibrational frequency.  I am given as data the value of ℏω in (cm^(-1))

Clear[vibTemp]
vibTemp[data_]:=Module[{kB=0.0000861771 (*eV*),
conv=1.2398*10^(-4)  (*convert cm^-1 into eV *)},
ϵv*conv/kB/.data (*ϵv is the vibrational energy in reciprocal centimeters*)]

vibTemp[dataH2]

vibTemp[dataH2]

vibTemp[dataI2]

Cell 2:  Calculate (Subscript[K, v]) from:

	(Subscript[K, v])  = 1/(Subscript[q, v])
	
	(Subscript[q, v])  =  
        exp[(Subscript[D, 0])/(Subscript[k, B])  T]/(1 - exp[-(Subscript[T, v])/T ] )
	
	(Subscript[T, v]) = (Subscript[ϵ, v])/(Subscript[k, B]) = ℏω/(Subscript[k, B])

Clear[qv]
qv[T_,data_List]:=
Module[{kB=0.0000861771 (*eV*),Tv},
Tv=vibTemp[data];
Exp[D0/(kB*T)] /(1-Exp[-Tv/T])]/.data (*D0 is the dissociation energy; the zero point is already subtracted from D0*)

qv[T, {}]

qv[300, dataH2]

The vibrational contribution to the equilibrium constant for the dissociation of (Subscript[H, 2]), (Subscript[I, 2]), and (Subscript[K, 2])

\(TableForm[
	
      Table[{T, 1/qv[T, dataH2], 1/qv[T, dataI2], 1/qv[T, dataK2]}, 
		{T, 
          100, 2000, 200}], 
	
      TableHeadings -> {None, {"T", \*""(Subscript[K, V]) for (Subscript[H, 2])"", \ \*""(Subscript[K, v]) for (Subscript[I, 2])"", \*""(Subscript[K, v]) for (Subscript[K, 2])""}}]\)

At room temperature (Subscript[K, v]) is very small.  To understand why, you must remember that 
	(Subscript[K, v])  = 1/
            (Subscript[q, v])  =  (1 - exp[-(Subscript[T, v])/T ] )exp[-(Subscript[D, 0])/(Subscript[k, B])  T].  
	
Here (Subscript[D, 0]) is the binding energy of the diatomic molecule.  When (Subscript[D, 0]) >> (Subscript[k, B]) T  ,  (Subscript[K, v]) is very small.  The dissociation of (Subscript[I, 2]) occurs at lower temperature than that of (Subscript[H, 2]) because  (Subscript[D, 0 ])for (Subscript[I, 2]) is so much smaller than that for (Subscript[H, 2]). It takes much less energy to dissociate (Subscript[I, 2]) and therefore it will happen at lower temperatures. The dissociation of (Subscript[K, 2]) occurs at even lower temperatures, because (Subscript[D, 0])for this molecule is smaller.
This result is consistent with the simplified model discussed in the lecture notes.

WorkBookSM.11.4  Rotational contribution to the equilibrium constants for  (Subscript[H, 2]),  (Subscript[I, 2]),  (Subscript[K, 2])dissociation.

The rotational partition function is 

	 (Subscript[q, r]) = T/(Subscript[T, r]) 
	 
for heteronuclear diatomics and  

	(Subscript[q, v]) = T/(2(Subscript[T, r])) 
	
for homonuclear diatomics. The rotational temperature is

	 (Subscript[T, r])=B/(Subscript[k, B])
	 
where B is the rotational constant.  
	The contribution to the equilibrium constant is (Subscript[K, r]) = 2 (Subscript[T, r])/T (for homonuclear diatomics) and is less than one at all temperatures we normally work at.

Cell 1:  Calculate the rotational partition function

Clear[qr]
\(	qr[T_, dataDiatomic_] := 
      Module[
		{kB = 
            0.6951\  (*Boltzmann\ in\ cm^\(-1\)*) }, 
		T/\((2\ B/kB)\) /. 
          dataDiatomic
	 (*for\ homonuclear\ only*) ]\)

qr[T, {}]

This gives the rotational partition function (dimensionless) if you insert B in (cm^(-1))and T in degrees K.

Cell 2: Calculate and display (Subscript[K, r])for (Subscript[H, 2,])  (Subscript[I, 2]), and (Subscript[K, 2]).

\(TableForm[
	
      Table[{T, 1/qr[T, dataH2], 1/qr[T, dataI2], 1/qr[T, dataK2]}, 
		{T, 
          100, 2000, 200}], 
	
      TableHeadings -> {None, {"T", \*""(Subscript[K, r]) for (Subscript[H, 2])"", \ \*""(Subscript[K, r]) for (Subscript[I, 2])"", \*""(Subscript[K, r]) for (Subscript[K, 2])""}}]\)

The rotational contribution works against dissociation. We see this mathematically because (Subscript[K, r])=  2(Subscript[T, r])/T is much smaller than one at the high temperatures at which we would be interested in dissociation. 
	A very curious thing is the large effect of the symmetry number. It multiplies  (Subscript[K, r]),  hence K,  with 2. This has a large effect on the equilibrium concentration.  For example, the equilibrium constant for dissociation of oxygen with one atom of mass 16 and the other with mass 17 is half that for the dissociation of an (Subscript[O, 2]) molecule in which both atoms are of mass 16.  The factor of 2 comes from σ. The  vibrational frequency and the B constant of \(TraditionalForm\`\(\^16\)O \(\^16\)O\) and \(TraditionalForm\`\(\^16\)O \(\^17\)O\)  are fairly close to each other, so they do not make a big difference in the equilibrium constant.

WorkBookSM.11.5. Put together the translation, rotation and vibrational contributions to calculate the equilibrium constants for (Subscript[H, 2]),  (Subscript[I, 2]),  and (Subscript[K, 2])dissociation.

Cell 1:   Compare the contributions from various kinds of motions to the equilibrium constants of (Subscript[H, 2])

\(TableForm[
	
      Table[{T, ScientificForm[NumberForm[Kt[T, \ 1, dataH, dataH2], 3]], 
          ScientificForm[NumberForm[1/qv[T, dataH2], 3]], 
          ScientificForm[
            NumberForm[1/qr[T, dataH2], 3]], 
		ScientificForm[
            NumberForm[
              Kt[T, \ 1, dataH, dataH2]/\((qv[T, dataH2]*qr[T, dataH2])\), 
              3]]}, 
		{T, 100, 2200, 200}], 
	
      TableHeadings -> {None, {"T", \*""(Subscript[K, t]) for (Subscript[H, 2])"", \ \*""(Subscript[K, v]) for (Subscript[H, 2])"", \*""(Subscript[K, r]) for (Subscript[H, 2])"", \*""K for (Subscript[H, 2])""}}]\)

As you can see the translational motion favors dissociation. This is mainly due to the fact that the translational entropy of two atoms is much larger than the translational entropy of the diatomic. The vibrational contribution opposes dissociation, mostly through the term exp[-(Subscript[D, 0])/(Subscript[k, B])  T]  which says that is rather improbable that the collisions between the molecules in the gas are likely to break the chemical bond. The rotational motion works against dissociation.

Cell 2:   Compare the contributions from various kinds of motions to the equilibrium constants of (Subscript[I, 2])

\(TableForm[
	
      Table[{T, ScientificForm[NumberForm[Kt[T, \ 1, dataI, dataI2], 3]], 
          ScientificForm[NumberForm[1/qv[T, dataI2], 3]], 
          ScientificForm[
            NumberForm[1/qr[T, dataI2], 3]], 
		ScientificForm[
            NumberForm[
              Kt[T, \ 1, dataI, dataI2]/\((qv[T, dataI2]*qr[T, dataI2])\), 
              3]]}, 
		{T, 100, 2200, 200}], 
	
      TableHeadings -> {None, {"T", \*""(Subscript[K, t]) for (Subscript[I, 2])"", \ \*""(Subscript[K, v]) for (Subscript[I, 2])"", \*""(Subscript[K, r]) for (Subscript[I, 2])"", \*""K for (Subscript[I, 2])""}}]\)

We see here that because the binding energy of (Subscript[I, 2]) is less than that of (Subscript[H, 2]), the (Subscript[I, 2]) molecule dissociates more easily.  The comments made for (Subscript[H, 2]) in the cell above apply to (Subscript[I, 2 ])also.

Cell 3:   Compare the contributions from various kinds of motions to the equilibrium constants of (Subscript[K, 2])

\(TableForm[
	
      Table[{T, ScientificForm[NumberForm[Kt[T, \ 1, dataK, dataK2], 3]], 
          ScientificForm[NumberForm[1/qv[T, dataK2], 3]], 
          ScientificForm[
            NumberForm[1/qr[T, dataK2], 3]], 
		ScientificForm[
            NumberForm[
              Kt[T, \ 1, dataK, dataK2]/\((qv[T, dataK2]*qr[T, dataK2])\), 
              3]]}, 
		{T, 100, 2200, 200}], 
	
      TableHeadings -> {{}, {"T", \*""(Subscript[K, t]) for  (Subscript[K, 2])"", \ \*""(Subscript[K, v]) for  (Subscript[K, 2])"", \*""(Subscript[K, r]) for  (Subscript[K, 2])"", \*""K for  (Subscript[K, 2])""}}]\)

Cell 4:  Compare next the equilibrium constants for (Subscript[H, 2]), (Subscript[I, 2]) and (Subscript[K, 2]):

\(TableForm[
	
        Table[{T, \ 
		Kt[T, \ 1, dataH, 
                dataH2]/\((qv[T, dataH2]*qr[T, dataH2])\), 
		Kt[T, \ 1, 
                dataI, dataI2]/\((qv[T, dataI2]*qv[T, dataI2])\), 
		Kt[
                T, \ 1, dataK, 
                dataK2]/\((qv[T, dataK2]*qr[T, dataK2])\)}, 
		{T, 100, 
            2200, 200}], 
	
        TableHeadings -> {{}, {"T", \*""K for (Subscript[H, 2])"", \ \*""K for (Subscript[I, 2])"", \*""K for (Subscript[K, 2])""\ }}] // 
NumberForm[#, 3] &\)

The enormous difference between these values is mainly due to the difference in the dissociation energy. There are other factors though. The translational contribution is  larger for the heavy molecules, however that is countered, to some extent, by the fact that the rotational contribution is smaller.