Research: Science and Education
Fraction of Molecules Exceeding a Given Energy William McInerny Department of Chemistry, City College of San Francisco, San Francisco, CA 94112
When the beginning chemistry student is introduced to the Maxwell–Boltzmann energy distribution in the context of molecular kinetic theory, several average values are claimed but not proved—for example, Eav = 3⁄2 RT. The determination of these average values is usually postponed until the course in physical chemistry, when the student is expected to evaluate the pertinent improper integrals over the interval (0,∞) (1). However, certain fractional quantities of significance are not readily expressed in the form of easily determined integrals, and these quantities must be evaluated by more sophisticated analytical methods or numerical or graphical means. One such quantity is the fraction of molecules that exceed some specified energy E * in the distribution. This fraction is of interest in assessing the number of sufficiently energetic molecules that can overcome some energy barrier to reaction (2). At City College of San Francisco, the Chemistry 101A course has long included an experiment that requires beginning students to evaluate such fractions by a simple graphical method based on weighing a sheet of paper of uniform density. Although the graphical approach is effective in determining results with 2–3-figure accuracy, there is generally no standard analytical result available for faculty and students to use in judging the overall quality of the graphical results for a broad range of energies and more than 3 significant figures. The expression of the energy fraction by analytical methods is not particularly difficult when the problem is recast in terms of the well-characterized error function or the F0(x) of quantum chemistry. In this paper the fraction is expressed in terms of the function F0(x), and a computer program in the Basic language is provided in a form suitable for implementation on a PC for the numerical determination of the fraction of molecules exceeding a given energy E *. The analytical expression for the normalized Maxwell–Boltzmann distribution density function with the kinetic energy (E ) as the independent variable is given by f (E ) = 2π(1/πRT )3/2E 1/2exp(E /RT ) The fraction of molecules with energies that exceed a specified energy value E * is expressed as the ratio of the two integrals indicated below. The integral in the denominator yields 1 upon integration from 0 to ∞, by virtue of the normalization. Thus, the fractional integral is determined by the numerator integral alone (3). ∞
f (E )dE Fraction(E ≥ E *) =
∞
E*
=
∞
2π E*
f (E )dE
1 πRT
3/2
E 1/2exp E dE RT
0
The ultimate integral is easily expressed in terms of the integral 1
F0 x =
0
exp xt 2 dt
by the substitutions E /RT = u 2 and z = (E */RT )1/ 2 and the attendant differential relation dE = 2uRTdu:
Fraction(E ≥ E *) = 4 π
∞
z
u 2exp u 2 du
This last integral can be simply separated into two integrals by dividing the u axis into two parts, (0,z) and (z,∞). One gets ∞
0
u 2exp u 2 du =
z
u 2exp u 2 du +
0
∞
z
u 2exp u 2 du
and – the integral to the left of the equal sign has the value √ π /4. The energy fractional integral now can be written in the form
Fraction(E ≥ E *) = 1 – 4 π
z
0
u 2exp u 2 du
A last integration by parts will yield the result
z Fraction(E ≥ E *) = 1 – 4 exp z 2 + 1 π 2 2
z
0
exp u 2 du
This result can be rewritten with one last change of variable zy = u to give the ultimate expression. Thus,
Fraction(E ≥ E *) = 1 – 2 z π 1 – 2 E* π RT
1
1/2
exp 0
1
0
exp z 2 y 2 dy – exp z 2 =
E *y 2 E * dy – exp RT RT
The final integral is the well-studied function F0(x) from quantum chemistry, where it is needed in the evaluation of two-electron integrals when employing Gaussian orbitals. The evaluation of the above fraction can be effected by standard functions available in all forms of Basic, and, except for the determination of the value of the F0(x), the programming required is easy and brief. The problem of evaluating F0(x) has been of interest among quantum chemists since the time of Boys’s first publication (4 ). The standard recursion formula for the determination of F0(x) is derived and explained in Shavitt’s article (5): that is,
Fn =
2x F n+1 + exp x 2n + 1
Further results by Cristofferson and Shipman in the 1970s advanced the evaluation schemes to 16-digit results that were valuable in determining the accuracy of the error function evaluation programs themselves (6 ). Recent research has produced low-order rational polynomial methods that are capable of 18-digit accuracy for the full range of the variable x (7 ). It should be noted that rational polynomial methods require the use of a computation method and a programming language to determine the specific coefficients employed.
JChemEd.chem.wisc.edu • Vol. 78 No. 6 June 2001 • Journal of Chemical Education
801
Research: Science and Education
Table 1. Fraction of Molecules with Energy ≥ E * E*/(J mol 1)
T/K
Fraction with E ≥ E*
11,000
300 500
0.8490329 0.9230260
15,000
300 500
0.2604860 0.4926220
10,000
300 500
0.0456383 0.1861804
There are, however, classical methods that are applicable to the integral evaluation of F0(x), and their implementation is relatively easy on the PC. For example, a direct 20-point Gaussian integration can produce 16 places of decimals but at the expense of large numbers of exponential evaluations (8; private communication with Harris, F. E., Mar 1992). The determination of F0(x) employed here uses two separate approaches for two different intervals of x. For 0 ≤ x ≤ 13.2, the downward recursion method explained by Shavitt is employed (5). Inasmuch as 7 digits of decimals are sought here, the recursion method starts at n = 44 and proceeds downward in unit steps to n = 0. The formula is Fn =
2
56 70
2x F n+1 + exp x 2n + 1
and here F44 is set to 0 while n is decreased in unit steps until n = 0 and F0(x) is obtained. This is the standard approach and it is effective for all x, but it requires larger initial values of n for increasing values of x. For all x > 13.2, the evaluation of F0(x) can be performed using a limiting functional form that produces 7 or more digits. The formula is F0(x) = (π/4x)1/2 (8). The above formulas suffice to complete the exact determination of the energy fraction expression for all E * ≥ 0. It is to be noted that even as the fraction of energies exceeding some specified energy E * is usually expressed in terms of the error function, it can equally well be expressed in terms of the F0(x) of quantum chemistry, and the resulting values are computed to 7 places of decimals. The effort to find an accurate and efficient method for an expression from the area of molecular kinetic theory can be fruitfully expressed in terms of well-researched analysis from quantum chemistry. Thus, investigations from one area of chemistry cross-fertilize those from another. The results given in Table 1 are indicative of the results required of Chemistry 101A students in our experiment on the Maxwell–Boltzmann distribution. Since six fractional results are needed, students must actually perform eight graphical integrations because the denominators (areas of whole graphs) in the fraction definition must be evaluated graphically for each temperature specified. A simple Basic program is shown below for the general determination of the above fraction for any finite cutoff energy
802
E *. The implementation of the program is now carried out using the HP 350 MHz computers recently acquired by the City College Department of Chemistry. Students in our Chemistry 101A laboratories can now calculate the required fractions and confirm their graphically determined results in minutes using the departmental computers, rather than in days, as in the past. The program follows.
99
CLS:PRINT:PRINT:PRINT:PRINT R=8.314472 INPUT ”INPUT ENERGY IN JOULES = “; E INPUT ”INPUT TEMP IN KELVIN = “; T XX =E/R/T PRINT ”E/RT = “; XX PI = 3.14159265 F = 0 A = EXP(-XX) IF XX > 13.2 THEN 56 B = 2*XX FOR I = 44 TO 0 STEP -1 C = 2*I + 1 F = (B*F +A)/C NEXT GOTO 70 F = SQR(PI/4/XX) P = SQR(PI) PRINT” F0 = “;F X1 = SQR(XX) FRACT = 1 - 2*X1*(F-A)/P PRINT”FRACTION = “; FRACT PRINT INPUT”DO YOU WISH TO STOP(Y,N)”; A$ Y$ = LEFT$(A$,1) IF Y$ = “Y” THEN 99 GOTO 2 END
This short program can be used to determine the function F0(x) itself, if one is so inclined. One simply sets the values of R and T to 1.0, and the variable E is assigned the value of x. The student can restart the calculation by entering any value at the stop question except the character Y. Compilation will improve performance, but it may not be useful when the student wants to modify the program for other uses. Literature Cited 1. Castellan, G. W. Physical Chemistry, 3rd ed.; Addison-Wesley: Reading, MA, 1983; pp 68–69. 2. Castellan, G. W. Ibid., p 70. 3. Castellan, G. W. Ibid., p 71. 4. Boys, S. F. Proc. R. Soc. London 1950, A200, 542. 5. Shavitt, I. In Methods of Computational Physics; Alder, B.; Fernbach, S.; Rotenberg, M., Eds.; Academic: New York, 1963; Vol. 2, p 8. 6. Shipman, L. L.; Christofferson, R. E. Comput. Phys. Commun. 1971, 2, 201–206. 7. Ralston, A.; Rabinowitz, P. A First Course in Numerical Analysis; McGraw-Hill: New York, 1978; p 98. 8. Handbook of Mathematical Functions; Abramowitz, M.; Stegun, I., Eds; National Bureau of Standards Applied Mathematics Series 55; U.S. Government Printing Office: Washington, DC, 1970; p 916.
Journal of Chemical Education • Vol. 78 No. 6 June 2001 • JChemEd.chem.wisc.edu