Evaluation of vibrational partition functions for polyatomic systems:quantum
versus
classical methods for and Ar
ÆÆÆ
CN
¤
H
2
O
Antonio Riganelli, Frederico V. Prudente and J. C. Varandas
*
Anto
nio
Departamento de Universidade de Coimbra
,
P

3049 Coimbra Codex
,
Portugal Qu•
mica
,
Recei
v
ed 3rd March 2000
,
Accepted 14th April 2000Published on the Web 18th May 2000
The vibrational partition function of and Ar
ÉÉÉ
CN systems is calculated within the framework of H
2
Oquantum and classical statistical mechanics. The phase space integral arising in the classical picture is solvedadopting an efficient Monte Carlo technique. The temperature dependence of the partition function for the twomolecules is exploited with a view to study the range of applicability of classical statistical mechanics. For thecase of Ar
ÉÉÉ
CN van der Waals complex the role played by freezing the CN bond is also analyzed.
1 Introduction
The partition function (
Q
) is a fundamental concept in statistical mechanics and accurate values of
Q
are an essentialingredient of any statistical treatment of chemical processes.Values of partition functions for polyatomic species are alsodemanded by astrophysics for modeling cool stellar atmospheres,
1
and are the starting point to derive several thermodynamic quantities. In principle one can obtain thepartition function once the energy levels of the system areknown. With the amazing growth of computational resources,much e†ort has been made
2
in recent years to obtain
Q
theoretically using Ðrstprinciple methods to calculate the energylevels (
i
.
e
., the rotational
È
vibrational spectrum) of small molecules with an accuracy close to high resolution spectroscopicexperiments. In spite of this, the calculations involved whendealing with real systems are formidable and the full description of the rotational motion is still out of reach for mostpolyatomic molecules.For this reason, historically, several strategies have beendeveloped to explore alternative routes for the explicit calculation of the vibrational
È
rotational energy levels. Theseinclude classical methods with quantum corrections,
3
semiclassical methods,
4
path integral methods,
5
and harmonic pluscorrection approximations.
6
In a previous paper,
7
hereafterreferred to as paper I, we developed an alternative approachbased on the solution of the phase space integral which arisesin classical statistical mechanics. The method is based on theMonte Carlo technique and has been shown to yield accurateresults for van der Waals systems over the whole range of temperatures studied when compared with the quantum calculations. It also turns out to be very simple to implement.However, it is well known that classical partition functionscannot give accurate results for strongly bound molecules atthe low temperature regimes. Messina
et al
.
5
have shown howto obtain accurate values of the partition function for stablemolecules such as and by calculating classically theH
2
O SO
2
partition function with the true potential function replaced byan e†ective one. It should be pointed out though that thepotential function employed in their work is a force Ðeld, andhence cannot describe correctly the dissociation limits of themolecule.
¤ Presented at the Workshop on Photodynamics from Isolated Molecules to Condensed Phases, Havana, Cuba, February 13
È
19, 2000.
For molecules described by realistic nonseparable potentials, the solution of the phase space integral is far morecomplex. The traditional way to tackle the problem by transforming the phase space integral
8
h
10
into a conÐgurationalintegral (this approach is implicitly used by Messina
et al
.
5
) isnot in this case a viable route. Indeed, for polyatomic molecules described by realistic potentials, it is not possible tocarry out the integration over momentum space separatelyfrom the coordinates.
11
Note that this route is usuallyadopted in conventional transition state theory
12
and in conformational studies of macromolecules,
13
with a harmonictype potential where
q
is a generalized[lim
q
?
=
V
(
q
)
\
O
,coordinate vector] being often adopted. Vieth
et al
.
14
proposed a method to estimate the partition function and equilibrium constant based on Monte Carlo simulations but herealso the potential function diverges at inÐnity as
q
increases.The aim of this work is to explore the limitations of theapproach suggested in paper I, and further develop it to attaingenerality. Thus, the Ðnal goal will be to obtain a general procedure to calculate the partition function of systems describedby realistic potential functions but avoiding the cumbersomecalculation of the energy levels. We limit our calculations tothe vibrational partition function. This allows us to shed somelight on the classical description limit when a quantity such asthe partition function is needed. Although work is in progressto extend the method to include overall rotation of the system,this is not believed to alter the major conclusions of thepresent work. To test the methodology we consider the prototype systems and Ar
ÉÉÉ
CN (where dots emphasize theH
2
Oweak nature of the van der Waals bond) which have beenmuch studied both experimentally and theoretically. The Ðrstis representative of molecules with a deep well, while thesecond is a Ñoppy molecule of the van der Waals type.The paper is organized as follows. In Section 2 we reviewthe theory, and in Section 3 the adopted Monte Carlo procedure. Section 4 presents the applications to the andH
2
OAr
ÉÉÉ
CN molecules. The conclusions are in Section 5.
2 Theoretical background
2.1 Quantum statistical mechanics
The partition function is usually expressed as
Q
\
Q
trans
Q
int
(1)DOI: 10.1039/b001746i
Phys
.
Chem
.
Chem
.
Phys
., 2000,
2
, 4121
È
4129
4121This journal is The Owner Societies 2000
(
where is the external translational contribution com
Q
trans
monly obtained assuming perfect gas formalism and is the
Q
int
internal (rotational
È
vibrational plus electronic) contributiongiven by
Q
int \
Q
vr
Q
e
(2)where and are the rotational
È
vibrational and the elec
Q
vr
Q
e
tronic partition functions, respectively. In this work, weassume that no electronic excited states are involved,
i
.
e
.,In turn, is obtained in quantum statistical mecha
Q
e \
1.
Q
vr
nics as an explicit summation of Boltzmann factors for all theenergy levels of the system, namely
Q
vrQ \
;
i
g
i
e
~
e
i
@
kT
(3)where is the degeneracy factor associated with level
i
. This
g
i
is given by (2
J
]
1) times the nuclear spin factor where
Jg
n
,is the total angular momentum quantum number. We performcalculations only for
J
\
0, and hence
g
i
\
1.A good example of the intricacies involved in the calculation of the rotational
È
vibrational partition function usingeqn. (3) for bound molecules is given in ref. 15, where Tenny

son and coworkers report a study for the molecule. TheH
2
Ocalculation of the energy levels was performed using the PJT2potential
16
up to 30000 cm
~1
, and the surface of Ho
et al
.
17
for the calculations above this energy. It is important to pointout that these calculations are extremely time consuming evenwhen adopting modern parallel computer technologies, andhence approximate strategies are needed to extend the calculations to high
J
values. The summation in eqn. (3) is thereforecommonly carried out only up to a cuto† energy value. Inthis work we assume as cuto† value the energy at which dissociation of the molecule occurs as will be discussed in thenext section. The role of metastable states will not beaddressed here (see ref. 18 and ref. 15 for a detailed discussionon this topic).
2.2 Classical statistical mechanics
In the classical statistical mechanics approach, the summationof energy levels of eqn. (3) is replaced by the multidimensionalphase space integral
19
Q
vrC \
1
h
n
P
ÉÉÉ
P
B
e
~
H
(
q
,
p
)@
kT
d
V
(4)with
H
being the classical Hamiltonian,
n
the number of degrees of freedom, and d
V
the volume element or ““extensionin phase spaceÏÏd
V
\
d
q
n
d
p
n
(5)where
q
stands for generalized coordinates and
p
for theassociated conjugate momenta. Note that from the transformation theory of the Hamiltonian
20,21
the volume element is aninvariant for any canonical transformation of variables whichdoes not explicitly involve the time, and hence d
V
in eqn. (4)always takes the formd
V
\
d
q
n
d
p
n
\
d
q
1
d
q
2
... d
q
n
d
p
1
d
p
2
... d
p
n
(6)This is independent of the coordinates adopted to express theHamiltonian, a result which will be used in the following.The subscript
B
in eqn. (4) implies that the hypervolume of integration is restricted to the phase space region corresponding to a bound state situation: 0
O
H
(
q
, where is
p
)
O
E
d
,
E
d
the dissociation energy of the molecule (we assume as reference energy the minimum of the potential well). In otherwords, we must select the region of phase space for which theinternal energy is smaller than the dissociation energy. Thus,we will assume as reference atom
È
diatom dissociation energythe one associated to the channel of lowest energy. For thetitle systems the diatomic fragments will be CN and OH.A question that arises naturally is which coordinate systemis more suitable to express the Hamiltonian and, as a result, tosolve efficiently the phase space integral. In paper I, we havechosen the set of Jacobi coordinates (
R
,
r
) to describe a vander Waals complex like Ar
ÉÉÉ
CN (BC is a bound diatomicmolecule) where only one dissociation channel is possible. Inour notation,
r
(modulus of
r
) represents the BC distance, and
R
(modulus of
R
) is the distance between Ar and the center of mass of CN. However, in the case of molecules such as H
2
O,the presence of two equivalent OH channels (we have two distinct fragments) makes the Jacobi coordinates unsuitable forour purpose. In such cases, the hyperspherical system of coordinates, where a unique radial coordinate allows a democraticdescription of all rearrangement or dissociation processes,suggests itself as the ideal choice. In the next section, wepresent the Monte Carlo technique used to evaluate the integral in eqn. (4) for the two di†erent coordinate systems: massweighted Jacobi, and hyperspherical coordinates.
3 Monte Carlo method
If only the vibrational motion is considered, one has
n
\
3and the integral in eqn. (4) is sixdimensional. The integrationcan then be performed using a Monte Carlo type approachfollowing the method described in paper I. In that work, itwas shown that the problem of efficiency in the calculation isa crucial one. For this reason, we follow an adaptation
7
of BarkerÏs
22
algorithm to calculate the density of states in thetransition state theory context. We discuss here the basic principles of the method, the algorithm for its implementation,and the new features introduced. The idea is near in spirit tothe well known importance sampling technique: the variablesare not sampled independently of each other but instead somekind of dependence is introduced. This leads to a nonuniformdistribution, and hence appropriate weighting factors need tobe used. The goal is to choose the sampling volume as close aspossible to the integrated one in order to maximize the efficiency of the Monte Carlo simulation. The main idea of thisalgorithm is to deÐne the hypervolume
B
of eqn. (4) by introducing simple bounds on each variable:
B
\
M
x
o
x
i
min
O
x
i
O
x
i
max
#
i
\
1, ..., 2
n
N
(7)where
x
\
(
q
, ..., ..., is a vector with 2
n
p
)
\
(
q
1
,
q
n
,
p
1
,
p
n
)components. The basic ingredients of the algorithm may besummarized in the following steps (for a formal demonstration, see ref. 22).1. Find the equilibrium values of
x
that deÐne the globalminimum of the Hamiltonian function
H
(
x
):
x
eq \
(
q
eq
,
p
eq
)where ..., and
p
eq \
(0, ..., 0). In other words
q
eq \
(
q
1eq
,
q
n
eq
)the molecule is set at its equilibrium geometry.2. Find the integration domain of with all the(
x
1min
,
x
1max
)
x
1
other variables at their equilibrium geometry:
H
(
x
1min
,
x
2eq
, ...,
x
2
n
eq
)
\
H
(
x
1max
,
x
2eq
, ...,
x
2
n
eq
)
\
E
d
(8)3. Sample randomly the variable within this range
x
1
according to:
x
1
S
\
x
1min ]
(
x
1max [
x
1min
)
m
(9)where
m
is a random number in the range [0, 1].4. Find the integration domain for with the(
x
i
min
,
x
i
max
)
x
i
Ðrst (
i
[
1) variables Ðxed at the sampled value and the other(2
n
[
i
) variables at the equilibrium value:
H
(
x
1
S
, ...
x
i
~1
S
,
x
i
min
,
x
i`
1eq
, ...,
x
2
n
eq
)
\
H
(
x
1
S
, ...,
x
i
~1
S
,
x
i
max
,
x
i`
1eq
, ...,
x
2
n
eq
)
\
E
d
(10)5. Sample randomly the variable inside this range accord
x
i
ing to eqn. (9).6. Repeat steps 4 and 5 for all the variables until the valueof the last variable is sampled.(
x
2
nS
)
4122
Phys
.
Chem
.
Chem
.
Phys
., 2000,
2
, 4121
È
4129
7. Calculate the weight for each sampled point
g
,..., according to
x
gS
\
(
x
1
S
,
x
2
nS
),
w
g
\
%
i
/12
n
(
x
i
max[
x
i
min
) (11)which represents the ““phase spaceÏÏ volume associated to
x
gS
.8. Calculate the integral in eqn. (4) as
S
\
;
g
/1
M
f
g
w
g
/
M
(12)where is the function to be integrated
f
g
\
exp[
[
H
(
x
gS
)/
kT
]and
M
is the total number of sampled points.The standard deviation associated with eqn. (12) is givenby
22
p
2\
[
M
(
M
[
1)]
~1
;
g
/1
M
(
f
g
w
g
[
S
)
2
(13)In the following subsections, we present the details of the sampling procedure in massweighted Jacobi coordinates andhyperspherical coordinates.
3.1 Sampling of phase space in Jacobi coordinates
The classical Hamiltonian in massweighted Jacobi coordinates (
R
,
r
,
r
where
r
is the angle between the two Jacobivectors) for
J
\
0 assumes the form (see Appendix for details)
H
\
12
k
C
P
R
2 ]
P
r
2]
A
1
R
2]
1
r
2
B
P
r
2
D
]
V
(
R
,
r
,
r
) (14)where the volume element is given by d
V
\
d
R
d
r
d
r
d
P
R
d
P
r
Note that the interatomic distances andd
P
r
. (
r
12
,
r
13
,
r
23
)internal massweighted coordinates are related by
r
23\
drr
132 \
A
m
2
drm
2]
m
3
B
2]
R
2
d
2 [
2
m
2
m
2]
m
3
rR
cos
r
r
122 \
A
m
3
drm
2]
m
3
B
2]
R
2
d
2 ]
2
m
3
m
2]
m
3
rR
cos
r
(15)where
d
is a mass scaling factor, and is the mass of the
i
th
m
i
atom [see Appendix, eqn. (24)].The approach followed to sample the coordinates (
R
,
r
,
r
)and associated momenta is that indicated in the(
P
R
,
P
r
,
P
r
)algorithm described above. When applying such a procedure,we should keep in mind that, for Ar
ÉÉÉ
CN, the
R
coordinateneeds to be sampled before
r
. Instead, for the orderH
2
O,of the sampling in
R
and
r
is irrelevant and can be donerandomly.
3.2 Sampling of phase space in hyperspherical coordinates
Except for an additive constant in the hyperangle
/
, we use amodiÐed version of JohnsonÏs hyperspherical coordinates
23
suggested by Varandas and Yu.
24
They consist of a hyperradius
o
, and two hyperangles (
h
,
/
) for the internal coordinates. Three external or rotational hyperangles (
a
,
b
,
c
)complete the set. As shown in the Appendix, the classicalHamiltonian for
J
\
0 assumes in terms of these coordinatesthe form
H
\
12
k
C
P
o
2]
1
o
2
A
16
P
h
2]
4
P
Õ
2
sin
2
(
h
/2)
BD
(16)where the volume element is d
V
\
d
o
d
h
d
/
d
P
o
d
P
h
d
P
Õ
,and the coordinate ranges are0
OoO
O
; 0
OhOn
; 0
O
/
O
2
n
(17)In terms of these coordinates, the interatomic distances aregiven by
r
122 \ 12
d
32
o
2
[1
]
sin(
h
/2)sin(
/
]
m
3
)]
r
132 \ 12
d
22
o
2
[1
]
sin(
h
/2)sin(
/
[
m
2
)]
r
232 \ 12
d
12
o
2
[1
]
sin(
h
/2)sin
/
] (18)where
k
\
1, 2, 3, are scaling param
d
k
2\
(
m
k
/
k
)(1
[
m
k
/
M
),eters, with
k
and
M
being the reduced and total mass of thesystem [see Appendix, eqn. (25)], and
m
3\
2tan
~1
(
m
2
/
k
)are the channel angles.
m
2\
2tan
~1
(
m
3
/
k
)To clarify the advantage of using this coordinate system totreat molecules with two or more symmetric dissociationchannels, we show in Fig. 1 the potential energy surface of at Ðxed values of the hyperradius
o
as a function of theH
2
Oother two hyperangles (
h
,
/
). Note that the potentialH
2
Oenergy surface has a deep minimum and two equivalentH
]
OH dissociation channels. As a result, three distinct situations corresponding to equilibrium, intermediate, and large
o
values are visible from Fig. 1. For the equilibrium
o
value (
i
.
e
.,
o
\
2.5813 panel (a) shows that the minimum of the
a
0
),potential energy surface occurs at
h
\
0.6205 rad and
/
\
3.2009 rad. At intermediate
o
values (
i
.
e
.,
o
\
4.0 the
a
0
)two equivalent H
]
OH channels begin to appear. Theasymptotic situation where the two OH channels are completely separated is illustrated for
o
\
5.0 in panel (c). Note
a
0
Fig. 1
Potential energy surface for the molecule as a functionH
2
Oof the hyperangles (
h
,
/
) for di†erent values of
o
: (a)
o
\
2.58 (b)
a
0
;
o
\
4.00 (c)
o
\
5.00
a
0
;
a
0
.
Phys
.
Chem
.
Chem
.
Phys
., 2000,
2
, 4121
È
4129
4123
Fig. 2
Values of the equilibrium hyperangles as a function of
o
forthe potential energy surface: (a)
/
eq
; (b)
h
eq
.H
2
O
that in both situations (intermediate and large
o
values) thetwo H
]
OH channels are described symmetrically.The sampling of phase space in these coordinates presentssome new features. The major points are: the equilibriumhyperangle values (
h
eq
,
/
eq
) are not constant for di†erent
o
values, and hence for molecules such as there will not beH
2
Oa unique pair of equilibrium hyperangle values (
h
eq
,
/
eq
). Thiscan be seen in Figs. 2 and 3. In particular, Fig. 2 shows thevalues of the equilibrium hyperangles for the potentialH
2
Oenergy surface plotted as function of the hyperradius
o
.Clearly, panel (b) of Fig. 2 shows that there is only one
h
eq
foreach
o
value. Conversely, panel (a) shows that there are twoequivalent values of
o
starting at
o
\
3.1 which are associ
a
0
,ated to the two possible OH
]
H dissociation channels. Of course, the plot of
/
eq
for the Ar
ÉÉÉ
CN molecule in Fig. 3shows only one dissociation channel.For the sampling procedure in hyperspherical coordinateswe start for convenience with the variable
o
. For this, the
o
min
and
o
max
values were obtained according to step 2 [eqn. (8)] of our algorithm. The equilibrium values of the two other internal coordinates (
h
,
/
) were then obtained (see Figs. 2 and 3) bya cubic spline interpolation of a uniform
o
grid. To warrant auniform sampling of phase space, for each sampled
o
, the twochannels in were selected randomly.H
2
O
Fig. 3
Values of
/
eq
as a function of
o
for the Ar
ÉÉÉ
CN potentialenergy surface.
4 Results and discussion
The strategies developed in the previous section are tested inthis section for two representative systems, andH
2
OAr
ÉÉÉ
CN.
4.1 H
2
O
The molecule is characterized by a deep well and playsH
2
Othe role of a benchmark system for both boundstate andreactive scattering calculations reaction]. This[O(
1
D)
]
H
2
may be understood from the intrinsic interest of the title molecule due to its importance in atmospheric chemistry, and alsoto the complexity of its electronic structure, which makes themodeling of the potential energy surface (PES) particularlychallenging (see ref. 25). In the present work we use the potential energy surface for the water molecule obtained by one of us
26
using the energy switching (ES) method. Such a PES hasspectroscopic accuracy, being also suitable for reactivedynamics calculations.
27
We calculate the vibrational partition function using(
Q
v
)di†erent methodologies: quantum statistical mechanics(QSM) and classical Monte Carlo (CMC) simulations, withthe Hamiltonian being expressed both in hyperspherical (HCMC) and massweighted Jacobi (JCMC) coordinates. Theenergy levels needed to get the QSM result according to eqn.(3) were calculated for each parity (even and odd) using theDVR3D program of Tennyson
et al
.
28
For these calculations,we have employed Radau coordinates and the followingoptimum parameters for the DVR basis set: 50 DVR pointsfor each radial coordinate, and 64 DVR points for the angularcoordinate. The Ðnal secular problem has a dimension 4000. Atotal of 460 and 373 vibrational states have been obtained forthe even and odd parities (respectively) at energies below dissociation cm
~1
. The calculated zero point
E
d\
42064.9energy is cm
~1
. In the HCMC calculations, we
e
0\
4710.329have sampled 10
7
points to calculate the integral in eqn. (4)using the procedure described in Section 3.2 with
o
min\
1.939and
o
max\
9.0 In turn, for the JCMC calculations, we
a
0
a
0
.have adopted the procedure described in Section 3.1 to samplethe same number of points using
R
min\
1.262
a
0
,
R
max\
8.210
r
min\
1.231 and
r
max\
6.532 (the
a
0
,
a
0
,
a
0
range of
r
is implicitly deÐned to be 0
È
n
). The HCMC andJCMC results, together with their respective errors calculatedusing eqn. (13), are reported in Table 1. In addition, wepresent the widely used results based on the harmonic oscillator (HO) approximation, which is (in conjunction with therigid rotor assumption) routinely used in the calculation of most of the JANAF thermodynamic functions
29
for polyatomic molecules. Note that the vibrational partition function isin this case given by
20,30
Q
vHO\
%
i
/13
exp(
[
hc
l
i
/2
kT
)1
[
exp(
[
hc
l
i
/
kT
) (19)where denotes the three fundamental frequencies of
l
i
H
2
O:
29
cm
~1
, cm
~1
, cm
~1
.
l
1\
3651.1
l
2\
1594.7
l
3\
3755.9Table 1 shows that the classical description in HCMC andJCMC formulations gives values of of the same order of
Q
v
magnitude only for temperatures above 1000 K, whileapproaching the quantum mechanical results as the temperature increases. For lower temperatures (
T
\
1000 K), theclassical values overestimate the exact ones by several ordersof magnitude. It is important to point out that the classicalcurve always lies above the quantum results, as it should. Acomparison between the results of the two di†erent classicalprocedures also shows clearly that HCMC produces betterresults than JCMC. This validates the assertion that thehyperspherical coordinate system is more adequate when twoor more equivalent dissociation channels are present. On theother hand, the HO approximation gives acceptable results
4124
Phys
.
Chem
.
Chem
.
Phys
., 2000,
2
, 4121
È
4129
Table 1
Vibrational partition function for the molecule
a
H
2
OHCMC
c
JCMC
e
T
/K QSM
b
Q
vC
e
d
Q
vC
e
d
HO
f
100 3.69(
[
30)
g
3.78(
[
17) 1.90(
[
17) 4.06(
[
11) 2.87(
[
11) 5.99(
[
35)500 1.31(
[
6) 4.48(
[
5) 1.50(
[
5) 1.70(
[
3) 1.05(
[
3) 1.44(
[
7)1000 1.28(
[
3) 4.61(
[
3) 7.49(
[
4) 1.96(
[
2) 9.71(
[
3) 4.25(
[
4)2000 5.77(
[
2) 9.38(
[
2) 8.59(
[
3) 1.28(
[
1) 3.29(
[
2) 3.29(
[
2)3000 2.91(
[
1) 3.80(
[
1) 2.24(
[
2) 4.16(
[
1) 5.67(
[
2) 1.96(
[
1)4000 8.21(
[
1) 9.70(
[
1) 3.9(
[
2) 1.01 8.33(
[
2) 5.89(
[
1)5000 1.77 1.98 5.75(
[
2) 2.04 1.13(
[
1) 1.316000 3.28 3.55 7.74(
[
2) 3.63 1.45(
[
1) 2.447000 5.47 5.79 9.84(
[
2) 5.93 1.80(
[
1) 4.098000 8.44 8.79 1.21(
[
1) 9.03 2.17(
[
1) 6.339000 12.22 12.61 1.44(
[
1) 12.00 2.558(
[
1) 9.2710000 16.84 17.24 1.69(
[
1) 17.87 2.94(
[
1) 13.0011000 22.26 22.70 1.94(
[
1) 23.61 3.35(
[
1) 17.6012000 28.42 28.86 2.20(
[
1) 30.20 3.76(
[
1) 23.18
a
The calculations are performed with the zero of energy at the minimum of the potential energy surface.
b
Quantum statistical mechanical result,this work.
c
Classical Monte Carlo simulation in hyperspherical coordinates, this work.
d
Error calculated with eqn. (13).
e
Classical Monte Carlosimulation in Jacobi coordinates, this work.
f
Harmonic oscillator approximation, this work.
g
The number in parentheses indicates the power of 10.
for this range of temperatures but worse results than the classical simulation for higher temperatures (
T
[
1000 K), asexpected for molecules characterized by deep potential wells.
4.2 Ar
ÆÆÆ
CN
The Ar
ÉÉÉ
CN van der Waals molecule is a weakly bound andvery Ñoppy system. The potential used to describe this molecule is a pairwise potential where the twobody terms havebeen represented by using the realistic extended Hartree
È
Fockapproximate correlation energy (EHFACE)
31,32
model. Theparameters of the potential for the Ar
ÉÉÉ
C and Ar
ÉÉÉ
N fragments can be found in ref. 33, and ref. 34 for the CN molecule.We calculate the vibrational partition function by adoptingtwo di†erent methodologies: QSM and CMC. Usually, thestudy of van der Waals molecules (such as Ar
ÉÉÉ
CN) is carriedout by considering the chemically stable diatomic as a rigidrotor
35
(in the present case, the CN molecule). To tackle thissystem, we have performed the classical calculation with theHamiltonian (
J
\
0) expressed in three di†erent ways: hyperspherical coordinates (HCMC), full (3D) massweightedJacobi coordinates (JCMC3D), and 2D massweightedJacobi coordinates (JCMC2D) where the coordinate
r
isÐxed at its equilibrium value. To get the QSM results according to eqn. (3) we have used the energy levels calculated inpaper I using the BOUND computer code.
36,37
The numberof energy levels obtained below the dissociation energy (
E
d\
58.0 cm
~1
) is 13, and the calculated zero point energy is
e
0\
cm
~1
. Note that this zero point energy is associated16.893only with the bending and the
R
stretching vibrational modesof the pseudodiatomic Ar
È
(CN). For the HCMC calculation,10
7
points have been sampled to obtain the integral in eqn. (4)using the procedure described in Section 3.2 with
o
min\
8.9and
o
max\
55.0 while, for the JCMC3D calculation, the procedure described in Section 3.1 has been adopted by samplingthe same number of points using
R
min\
8.72
a
0
,
R
max\
100.00
r
min\
1.755 and
r
max\
1.791 Note
a
0
,
a
0
a
0
.that the JCMC2D calculations di†er from the JCMC3Dones only because the coordinate
r
has been kept Ðxed at itsequilibrium value (
r
eq\
1.772 The HCMC, JCMC3D,
a
0
).JCMC2D, and the QSM results are presented in Table 2. Inturn, Fig. 4 shows the logarithm of the di†erent values as a
Q
vC
function of temperature. Clearly, the HCMC and JCMC3Dresults di†er by orders of magnitude from the QSM results,and only the JCMC2D calculations give good agreementwith the QSM ones at all temperatures. However, it is alsoclear from Fig. 4 that the HCMC and JCMC3D resultsdi†er from the JCMC2D ones approximately by a constantfactor. This suggests that the phase space volume associatedto the zero point energy has not been properly treated in theclassical calculations (for a discussion on this topic in relationto classical reaction dynamics, see ref. 38 and referencestherein). In order to rationalize this observation, one may estimate how many states of the Ar
ÉÉÉ
CN complex have been
Table 2
Vibrational partition function for the Ar
ÉÉÉ
CN molecule
a
JCMC2D
d
T
/K HCMC
b
JCMC3D
c
Q
vC
e
e
QSM
f
100 4.48(
[
3)
g
4.90(
[
2) 7.28 2.94(
[
2) 7.07200 5.98(
[
3) 6.74(
[
2) 9.86 3.75(
[
2) 9.55300 6.60(
[
3) 7.51(
[
2) 10.94 4.09(
[
2) 10.57400 6.94(
[
3) 7.92(
[
2) 11.52 4.27(
[
2) 11.13500 7.14(
[
3) 8.13(
[
2) 11.89 4.38(
[
2) 11.48600 7.28(
[
3) 8.36(
[
2) 12.14 4.46(
[
2) 11.72700 7.39(
[
3) 8.49(
[
2) 12.32 4.52(
[
2) 11.89800 7.47(
[
3) 8.59(
[
2) 12.46 4.56(
[
2) 12.03900 7.53(
[
3) 8.67(
[
2) 12.57 4.59(
[
2) 12.131000 7.58(
[
3) 8.73(
[
2) 12.66 4.62(
[
2) 12.21
a
The calculations are performed with the zero of energy at the minimum of the potential energy surface.
b
Classical Monte Carlo simulation inhyperspherical coordinates, this work.
c
Classical Monte Carlo simulation (3D) in Jacobi coordinates, this work.
d
Classical Monte Carlo simulation (2D) in Jacobi coordinates, this work.
e
Error calculated with eqn. (13).
f
Quantum statistical mechanics, this work.
g
The number inparentheses indicates the power of 10.
Phys
.
Chem
.
Chem
.
Phys
., 2000,
2
, 4121
È
4129
4125