Laplace transformation approach, saddle point approximation method for calculation of number and density of vibrational states

JOURNAL OF SCIENCE OF HNUE Natural Sci., 2008, Vol. 53, N ◦ . 5, pp. 81-88 LAPLACE TRANSFORMATION APPROACH, SADDLE POINT APPROXIMATION METHOD FOR CALCULATION OF NUMBER AND DENSITY OF VIBRATIONAL STATES Tran Vinh Quy and Nguyen Van Hanoi National University of Education To Ba Ha Minitary Technical Academy Nguyen Dinh Do Hanoi University of Mining and Geology Abstract. In the paper the calculation of density of molecular vibrational states has been performed by the Laplace transformation approach, i.e re- duced to the evaluation of inverse Laplace transform of partition function by the aid of saddle point approximation method: The exponential of inter- grand in inverse Laplace integral is presented in the Taylor series, in which at the saddle point the first derivative vanishes, then, the density of vibrational states - the inverse Laplace transform is integrated through the steepest de- scent which is in imaginary direction having the saddle point as abscissa on the real axis. Finally, this process of treament gives the expression of density of vibrational states containing terms from first up to fifth orders of accuracy, which have been easily calculated on computer. 1. Introduction In the reaction rate calculation of the investigation of unimolecular reaction kinetics the number of states and the density of quantum vibrationalstates play an important role. There have been, so far, many methods of calculating the numbers of states and the densities of quantum vibrational states of reacting molecules and activated complexes: Method of direct counting; Whitten-Rabinovith method; Laplace trans- formation method followed by saddle point approximation. However, before, all sci- entific projects using Laplace transformation-saddle point method [1-4] have con- fined to the second order of approximation only, and, in evaluating the unimolecular rate constant by RRKM method, one usually made a mistake of the methodologic unconsistency (using different methods for the calculating the numbers of quantum states: the method of direct counting for activated complexes at low energy, and 81 Tran Vinh Quy, Nguyen Van, To Ba Ha and Nguyen Dinh Do Whitten-Rabinovith method for reacting molecules at high energy), as a result, this unconsistency led to great deviation of theory from experiments at measurements. In this paper we make use of the Laplace transformation method-saddle point ap- proximation with very high accuracy for calculating numbers of states and densities of states of both reacting molecules and activated complex. We have developed the representation of Eyring H., Lin S. H., Lin S. M., [1] in order to obtain results with high accuracy up to fifth, sixth,... orders of saddle point approximation and then we applied to evaluate the rate constant of methylisocyanide isomerization. 2. Laplace transformation method and Saddle point ap- proximation 2.1. Density of quantum vibrational states as the inverse Laplace transform of vibrational statistical sum Statistical sum of harmonic oscillator system has the form Q(β) = ∑ i e−βEi = ∞∫ E0 e−βEW (E)dE, (2.1) where, E0 is the lowest energy level and W(E) is density of states of the system. From the mathematical view point, the integral of (2.1) has the form of Laplace integral, therefore, statistical sum Q(β) is the Laplace transform of real function W(E). Density of states satisfies these conditions W ≥ 0, lim E→∞ W (E)e−αE = 0, (α > 0), (2.2) From (2.1) and (2.2) we get to following properties: 1- If β is real then integral (2.1) converges for β > 0 and Q(β) is bounded. If β is complex, the mentioned above affirmation is right for Reβ > 0 (Reβ is real part value of β). 2- For β > 0 we have ∂2lnQ ∂β2 > 0. (2.3) 3- The inverse formule of (2.1), presenting density of states W(E) via statistical sum Q(β), has the form [4;119 or 3;122] W (E) = 1 2pii ∫ C Q(β)eβEdβ, (2.4) 82 Laplace transformation approach, saddle point approximation method... where the integral contour C contains a straight line parallel to the imaginary axis Reβ = β'= const and a half-round of radius R = ∞ surrounding the left half-plane of the mentioned straight line (when Q(β) has singularity β with Reβ < β ') or surrounding the right half-plane of the straight line (when Q(β) has the singularity β with Reβ > β'). Integral (2.4) is said to be the inverse Laplace integral and function W(E) is considered as inverse Laplace transform of statistical sum Q(β). Q(β) (see(3.1)) has the singularity at β = 0, therefore, half-round of R =∞ surrounds the left half-plane of the straight line Reβ = β', and, due to (2.2), integral on this half-round must vanish. The integration of inverse Laplace transform now is performed only in the arbitrary straight line parallel to imaginary axis W (E) = 1 2pii β′+i∞∫ β′−i∞ Q(β)eβEdβ, (2.5) or we can write as W (E) = 1 2pii β′+i∞∫ β′−i∞ Q(β)eβEdβ = 1 2pi +∞∫ −∞ Q(β ′ + iβ”)e(β ′+iβ”)Edβ” = 1 2pi +∞∫ −∞ exp[f(β)]dβ”, (2.6) where we have made use of the relation β = β ′ + iβ” and considered the function f(β) as equal to f(β) = lnQ(β) + βE = lnQ(β ′ + iβ”) + (β ′ + iβ”)E. We now discuss the calculation of integral (2.6) of harmonic vibrational oscil- lator system by using saddle point approximation method. 2.2. Saddle point approximation method If f(β) = u+ iv, then it can be expected that most of contribution to integral (2.6) is given by segment of contour in which u has great values. The idea of saddle point method is that, one must deform the contour C so that the region of great values of u is pressed to become the most narrow region that can be made. Con- cretely, we deform contour C so that on surface u(β ′, β”) at the saddle point β0 the contour enters valleys and, hence, values of function u decrease very quickly. This is the steepest descent line and the near neighbourhood of saddle point β0 will be the region with great values of u, which has been pressed. The main contribution to integral (2.6) is given by this region. We now expand function f(β) at the neighbourhood of β ′ in the imaginary direction to get f(β) = lnQ(β′) + β′E + [ E + ( ∂lnQ ∂β ) β=β′ ] (iβ”) + ∞∑ n=2 1 n! ( ∂nlnQ ∂βn ) β=β′ (iβ”)n . (2.7) 83 Tran Vinh Quy, Nguyen Van, To Ba Ha and Nguyen Dinh Do If the function has the minimum [( ∂lnQ ∂β ) β=β′ > 0 ] at the point β∗ = (β∗, 0) given by the condition f ′(β∗) = 0 or ( ∂lnQ ∂β ) β=β∗ = −E, (2.8) then this must be a saddle point in the complex plane. Since f”(β∗) is real, it fol- lows that the second term, the biggest term in expanse (2.7) with iβ” = seiφ can be written as + 1 2 ( ∂2lnQ ∂β2 ) β=β∗ (iβ”)2 = + 1 2 f”(β∗)s2ei2φ = + 1 2 f”(β∗)s2 (cos2φ+ isin2φ) . (2.9) At that time in order that the function f(β) from the saddle point β∗ = (β∗, 0) decreases at fastest speed entering valley, i.e. become steepest descent, one must have the condition cos2φ = −1 or 2φ = ±pi and φ = ±pi 2 . From the mentioned arguments we come to the conclusion that the steepest descent of function f(β) must be in the direction of imaginary axis. In this direction the imaginary part of f(β) almost never varies (because sin2φ = sinpi = 0), therefore the multiplier eiv in the integrand of (2.6) will not create harmful oscillation. Assigning to β ′ the value of saddle point β∗ given by equation (2.8) we obtain f(β) = lnQ(β∗) + β∗E + 1 2 ( ∂2lnQ ∂β2 ) β=β∗ (−β”2) + ∞∑ n=3 1 n! ( ∂nlnQ ∂βn ) β=β∗ (iβ”)n .(2.10) Putting equation (2.10) into (2.6) we get W (E) = Q(β ∗) 2pi eβ ∗E +∞∫ −∞ exp [ −1 2 ( ∂2lnQ ∂β2 ) β=β∗ (iβ”2) + ∑∞ n=3 1 n! ( ∂nlnQ ∂βn ) β=β∗ (iβ”)n ] dβ”. (2.11) Introducing the notation bn = 1 n! ( ∂nlnQ ∂βn ) β=β∗ , (2.12) we can write Eq(11) in the form W (E) = Q(β∗) 2pi eβ ∗E +∞∫ −∞ dβ”e−b2β” 2 ∞∑ n=0 Bnβ” n, (2.13) 84 Laplace transformation approach, saddle point approximation method... where Bn is the coefficient at β” n in the following expression ∞∏ k=3 exp [ bk (iβ”) k ] = ∞∑ n=0 Bnβ” n. (2.14) By the way of direct multiplying out each other all factors (in the Taylor expansion form) in the left hand side of equation (2.14) we can easily define values of coefficient Bn. Concretely, after eliminating all coefficients corresponding imaginary terms, we can obtain some coefficients in the sum of Eq (2.14) as B0 = 1, B2 = 0, B4 = b4, B6 = − ( b6 + 1 2 b23 ) , ..., B2k = ∑ in,jn (−1)k [ bj1i1 .b j2 i2 ... (j1!).(j2!)... ] , (2.15) with i, j satisfying the condition i1j1 + i2j2 + i3j3 + ... = 2k. (2.16) The first term of Eq (2.13) corresponding to the lowest order of approximation, after performating integration, give us the result W (E) = Q(β∗) 2pi eβ ∗E +∞∫ −∞ e−b2β” 2 dβ” = Q(β∗) 2pi eβ ∗E √ pi b2 = Q(β∗)eβ ∗E[ 2pi ( ∂2lnQ ∂β2 ) β=β∗ ] 1 2 .(2.17) All values of integral (2.13) corresponding to odd values of n (all coefficients at imaginary terms (iβ”)n) will vanish due to the oddness of the integrand of (2.13). Therefore W (E) will be expessed as a real function. Performing all intergrations in (2.13) and changing the notation for the sum, we get to W (E) = Q(β∗) 2pi eβ ∗E [ ∞∑ k=1 B2k 1.3.5...(2k − 1) 2k √ pi b2k+12 + √ pi b2 ] , (2.18) where we have used the relation ∞∫ −∞ dβ”β”2ke−b2β” 2 = 1.3.5...(2k − 1) 2k √ pi b2k+12 . (2.19) Equation (2.18) can be rewritten as W (E) =W1(E) [ 1 + ∞∑ k=1 B2k (2k)! 2k.k! ( ∂2lnQ ∂β2 )−k β=β∗ ] . (2.20) 85 Tran Vinh Quy, Nguyen Van, To Ba Ha and Nguyen Dinh Do or can be equal to W (E) = W1(E) [ 1 +B2 +B4 4! 22.2!.22.b22 +B6 6! 23.3!.23.b32 + B8 8! 24.4!.24.b42 + B10 10! 25.5!.25.b52 + ... ] . After substituting values of Bn given by (2.15) we have W (E) = W1(E) [ 1 +B2 + B4 4! 22.2!.22.b22 +B6 6! 23.3!.23.b32 +B8 8! 24.4!.24.b42 + B10 10! 25.5!.25.b52 + ... ] = W1(E) [ 1 + ( 3 4 b4 b22 − 15 16 b23 b32 + 10395 384 b34 b62 − 945 64 b23.b4 b52 ) + 105 32 ( b24 b42 + 2 b3b5 b42 + 99 2 b3.b4.b5 b62 − 9 2 b25 b52 ) + 945 32 ( 11 4 b26 b62 ) − ( 4 63 b6 b32 − b4.b6 b52 + 11 2 b5.b7 b62 + 11 4 b23.b6 b62 − b3.b7 b52 ) + 945 64 ( 11 b3.b9 b62 + 11 b4.b8 b62 + 4 9 b8 b42 − 2 b10 b52 ) + ... ] , (2.21) Using the notation ∂nβ lnQ = ∂ nlnQ/∂βn for expressions (2.12) of bnwe obtain W (E) = W1(E) [ 1 + ( ∂4βlnQ 8 ( ∂2βlnQ )2 − 5 ( ∂3βlnQ )2 24 ( ∂2βlnQ )3 + 385 ( ∂3βlnQ )4 1152 ( ∂2βlnQ )6 + 385 ( ∂4βlnQ )3 3072 ( ∂2βlnQ )6 − 35 ( ∂3βlnQ )2 ( ∂4βlnQ ) 32 ( ∂2βlnQ )6 ) + 105 32 ( ( ∂4βlnQ )2 36 ( ∂2βlnQ )4 + 2 ( ∂3βlnQ ) ( ∂5βlnQ ) 45 ( ∂2βlnQ )4 + 11 ( ∂3βlnQ ) ( ∂4βlnQ ) ( ∂5βlnQ ) 60 ( ∂2βlnQ )6 − ( ∂5βlnQ )2 100 ( ∂2βlnQ )5 ) + 945 32 ( ( 11∂6βlnQ )2 32400 ( ∂2βlnQ )6 − 2∂6βlnQ 2835 ( ∂2βlnQ )3 − ( ∂4βlnQ ) ( ∂6βlnQ ) 45 ( ∂2βlnQ )5 + 11 ( ∂5βlnQ ) ( ∂7βlnQ ) 18900 ( ∂2βlnQ )6 + 11 ( ∂3βlnQ )2 ( ∂6βlnQ ) 270 ( ∂2βlnQ )6 − ( ∂3βlnQ ) ( ∂7βlnQ ) 945 ( ∂2βlnQ )5 ) + 945 32 ( ( 11∂6βlnQ )2 32400 ( ∂2βlnQ )6 + ( 11∂3βlnQ ) ( ∂9βlnQ ) 34020 ( ∂2βlnQ )6 + ( 11∂4βlnQ ) ( ∂8βlnQ ) 15120 ( ∂2βlnQ )6 − ∂10β lnQ 56700 ( ∂2βlnQ )5 ) + ... ] . (2.22) 86 Laplace transformation approach, saddle point approximation method... This is the convergent series and is the equation applied for calculating the density of states W (E) with accuracies up to the fifth, sixth orders of approxima- tions. It will be reduced to accuracy of second order of approximation, when terms of third, fourth, fifth orders are neglected. 3. Application We now make use of the equation (2.22) for calculating the density of states W (E) of methylisocyanide system. Vibrational statistical sum of this system is cal- culated by the aid of formula Q(β) = s∏ i=1 ( exp ( −1 2 hνiβ ) (1− exp (−hνiβ)) ) , (3.1) where νi is vibrational frequency of i-th oscillator. Differentiating lnQ(β) with respect to β at different orders, we obtain ∂lnQ(β) ∂β = s∑ i=1 [ − 1 2 hνi − hνi (ehνiβ − 1) ] ; ∂2lnQ(β) ∂β2 = s∑ i=1 [ (hνi) 2 e hνiβ (ehνiβ − 1)2 ] ; ... (3.2) ∂10lnQ(β) ∂β10 = s∑ i=1 [ (hνi) 10 e hνiβ (ehνiβ − 1)2 − 510. (hνi) 10 e 2hνiβ (ehνiβ − 1)3 + 18150. (hνi) 10 e 3hνiβ (ehνiβ − 1)4 − 186480. (hνi) 10 e 4hνiβ (ehνiβ − 1)5 + 834120. (hνi) 10 e 5hνiβ (ehνiβ − 1)6 − 1905120. (hνi) 10 e 6hνiβ (ehνiβ − 1)7 + 2081520. (hνi) 10 e 7hνiβ (ehνiβ − 1)8 − 1814400. (hνi) 10 e 8hνiβ (ehνiβ − 1)9 + 1088640. (hνi) 10 e 9hνiβ (ehνiβ − 1)10 ] . (3.3) Vibrational frenquencies of (ethylisocyanide C2H5NC) given by spectroskopic experiments are ν(cm−1) = 2953(5), 2161, 1441(4), 1289(2), 1080(3), 945, 783, 531, 270, 224(2). For numerical calculation of density and numbers of quantum states at differ- ent energies Ev we make use of a program written in the language Pascal and we use special algorithm of "Diminishing iteration step" to define value of β [5;41-105]. Results obtained are listed in Table 1. 87 Tran Vinh Quy, Nguyen Van, To Ba Ha and Nguyen Dinh Do These results of the calculation will be applied to evaluate the rate constant of unimolecular reaction which we will present in next paper. In that paper we will have a Figure describing a curve of dependence between ln kuni k∞ and lnP that it may agree with experiments. Thus, in this paper we have reached a good accurary due to the fact that we have used the saddle point method for high orders of approximation and have ovecomed the methodologic unconsistency. Table 1. The results of the numerical calculation of density and numbers of quantum states at different energies Ev applied Turbo Pascal Ev N(E) W(E) Ev N(E) W(E) (kcal/mol) (kcal/mol) 46.00 3.768961E+06 1.994565E+06 57.00 4.559324E+06 7.905007E+06 47.00 1.430794E+06 1.044164E+06 58.00 5.586950E+06 1.010582E+07 48.00 1.168275E+06 1.028400E+06 59.00 6.856920E+06 1.290442E+07 49.00 1.180464E+06 1.188587E+06 60.00 8.424897E+06 1.645740E+07 50.00 1.298132E+06 1.452943E+06 61.00 1.035954E+07 2.096169E+07 51.00 1.489786E+06 1.820571E+06 62.00 1.274549E+07 2.666476E+07 52.00 1.751978E+06 2.308784E+06 63.00 1.568707E+07 3.387750E+07 53.00 2.091530E+06 2.945342E+06 64.00 1.931294E+07 4.299051E+07 54.00 2.521333E+06 3.767749E+06 65.00 2.378181E+07 5.449437E+07 55.00 3.059379E+06 4.824524E+06 66.00 2.928948E+07 6.900530E+07 56.00 3.728945E+06 6.177510E+06 67.00 3.607767E+07 8.729721E+07 REFERENCES [1] H. Eyring, S.H. Lin, S.M. Lin, 1980. Basic Chemical Kinetics. John Whiley and Sons, Inc. [2] P.J. Robinson, K.A. Hoolbruk, 1975.Unimolecular Reactions.Wiley Interscience- London-New York-Sydney-Toronto. [3] R. Kubo, 1967. Statistical Mechanics an Advanced course with problems and solution. North-Holland Publishing House, Amsterdam. [4] J. Mathews, R.L. Walker, 1971. Mathematics for the physical Sciences. Science and Technics Publishing House, Hanoi. [5] Tran Vinh Quy, 2006. Textbook forming a course of the Chemistry Infor- mation. Hanoi National University of Education Publishing House. 88