Molecular Vibrations

C. David Sherrill
School of Chemistry and Biochemistry
Georgia Institute of Technology
Why Estimate Molecular Vibrations?

• Simulation of vibrational spectrum (identification of

molecules)

• Vibrational corrections to enthalpy

• (Small) vibrational corrections to polarizability and other

properties

• Understanding of vibrational motion could assist dynamics

experiments and “mode-selective” chemistry
Small Vibrations in Classical Mechanics
The classic reference is Wilson, Decius, and Cross, Molecular
Vibrations (Dover, New York, 1980). Cheap book, makes a
good reference.
Let us focus on purely classical systems at first; all the results
carry over to quantum mechanics.
For small vibrations, the motion of atom α away from its
equilibrium value may be described by ∆xα , ∆yα , ∆zα , with
kinetic energy
Ã !2 !2 !2 
N
à Ã
1X d∆xα d∆yα d∆zα 
T = Mα  + +
2 α=1 dt dt dt
If we switch to mass-weighted coordinates, such as q1 =
√ √ √ √
M1 ∆x1 , q2 = M1 ∆y1 , q3 = M1 ∆z1 , q4 = M2 ∆x2 , etc.,
then the kinetic energy operator becomes simpler since the
mass factors are now absorbed
3N
1X
T = q̇i2
2 i=1
3N 3N 2
à ! à !
X ∂V 1X ∂ V
V = V0 + qi + qi qj + · · ·
i=1 ∂qi 0
2 i=1 ∂qi ∂qj 0
(1)

Remember that at equilibrium, (∂V /∂qi )0 = 0; we can also set

V0 = 0. Also abbreviate (∂ 2 V /∂qi ∂qj )0 as just fij .
Newton’s Equations of Motion
We can rewrite Newton’s equations of motion as
d ∂T ∂V
+ =0 j = 1, 2, · · · , 3N
dt ∂ q̇j ∂qj
or
3N
X
q̈j + fij qi = 0
i=1

A possible solution to this equation is

³√ ´
qi = ai cos λt + φ
√ q
where the angular frequency is λ; this is just k/m in
harmonic oscillator — the m has been absorbed by the mass-
weighted coordinate system used here!
Substitute the last expression into the differential equations to
get
3N
X
(fij − δij λ) ai = 0 j = 1, 2, · · · , 3N
i=1

or in matrix notation, just F a = λ a. This is an eigenvalue

equation! We have a solution to this system of 3N linear
equations only if λ has special values obtainable from the
secular determinant
¯ ¯
¯
¯ f11 − λ f12 f13 · · · f1,3N ¯
¯
¯ ¯
f21 f22 − λ f23 f2,3N
¯ ¯
¯ ··· ¯
¯=0
¯ ¯
¯ .. .. .. .. ..
¯
¯ . . . . . ¯
¯
¯ ¯
f3N,1 f3N,2 f3N,3 · · · f3N,3N −λ ¯
¯ ¯
¯
Normal Modes of Vibration
The matrix eigenvalue equation is equivalent to matrix diago-
nalization which is equivalent to solving the secular determinant
for each λ (N of them). Once we have the eigenvalues λk we
can get the corresponding eigenvectors ak , giving the motion
of each atom for the given eigenvalue λk :
µq ¶
qik = aik cos λk t + φ k .

The eigenvectors ak are the normal modes of vibration. For

each normal mode, all the atoms move with the same frequency
and phase, but with different amplitudes.
Normal Coordinates
We can define a new set of coordinates using the normal modes.
This gives us the “normal coordinates”
3N
X
Qk = aik qi k = 1, 2, · · · , 3N
i=1

Since the eigenvectors of a real, symmetric matrix (F) are

orthogonal, T and V become diagonal (no cross terms):
3N
1X
T = Q̇2k
2 k=1
1
V = λk Q2k
2
The Hamiltonian is separable — in normal mode coordinates,
we can solve for the motion in each coordinate separately!
Polyatomic Molecules
What happens for quantum mechanics, and for polyatomic
molecules? Use Harmonic Oscillator model.

• 3N -6 frequencies (3N -5 for linear molecules); the rest are

translations and rotations with zero frequency

• In normal mode coordinates, Hamiltonian is separable:

wavefunction is a product and energy is a sum. Total
vibrational energy is i ωi h̄(vi + 1/2)
P

• Minimum energy (due to uncertainty principle) is “zero

point vibrational energy” (ZPVE or ZPE), where vi = 0
1 P
for all i. ZPVE = 2 h̄ i ωi
How Would We Get Harmonic Frequencies for a
Molecule?

• Easy — just diagonalize the second derivative matrix F,

called the Hessian. The frequencies ωi are the square roots
√
of the eigenvalues, λi .

• Where do we get F?

• Recall fij = (∂ 2 V /∂qi ∂qj )

• Potential energy V is just Ee (B.O. approximation!): Need

∂Ee2 /∂qi ∂qj .

• Compute second derivative of Ee in terms of of Carte-

sian displacements (xα , yα , zα , call them q̃i ) and it’s
easy to transform to mass-weighted coordinates, using
F = M−1/2 F̃ M−1/2 .

• How do we get ∂ 2 Ee /∂xα ∂yβ , etc? Need second derivative

of electronic energy vs nuclear coordinates. Can compute
analytically (using formula) or numerically from finite
differences of energies or gradients:
!¯ !¯
2
"Ã Ã #
∂ Ee ∂Ee ¯¯ ∂Ee ¯¯
≈ − /∆x
∂xα ∂yβ ∂yβ ¯ ∂yβ ¯
¯ ¯
xα =xα0 +∆xα xα =xα0 −∆xα
Analytic Hessian Better than Numerical

• Analytic Hessian might cost ∼ 10-30x cost of energy;

analytic gradient costs maybe ∼ 1.5-2x.

• Can need gradients from many displaced geometries —

up to 6N (+ and - for each of 3N coordinates) — unless
reduced by point group symmetry

• Numerical Hessian contains numerical errors (divide small

number by small number)

• (Can land on wrong solution if displacement drops symme-

try)
Availability of Analytic Derivatives
Method Gradient Hessian
HF, DFT Y Y
CI Y N
CCSD, CCSD(T) Y S
MP2 Y S
CASSCF Y S
CIS Y Y
EOM-CCSD Y Y
TD-DFT N N

S = available in some packages; Y = widely available

Approximate Average Errors in Harmonic Frequencies
(Using polarized double and triple zeta basis sets)
Method Error
HF 11%
CISD 4-6%
CCSD 1-4%
CCSD(T) 1-3%

Anharmonicity accounts for another ∼ 2-3% difference from

experimental fundamental frequencies. Many workers employ
scaling factors for each level of theory to better predict
fundamental frequencies.
Spectra of SiC2H2 isomers: DZP CCSD(T)
100
Silapropadienylidene
Silacyclopropyne
Experiment

80

60
Intensity

40

20

0
0 500 1000 1500 2000 2500
Frequency (1/cm)
Scaling ZPVE’s
In an enlightening paper, Grev, Janssen, and Schaefer [J. Chem.
Phys. 95, 5128 (1991)] showed that using scaled fundamental
frequencies to estimate the ZPVE is not necessarily better than
using unscaled frequencies. The reason is anharmonicity.
If ZPVE’s use scaling, they should have a different scaling
factor than the individual frequencies.
1 1 1
X µ ¶ X µ ¶µ ¶
G(v) = ω r vr + + χrs vr + vs + + ···,
r 2 r≥s 2 2
f und f und 3X 1X
∆ = G(0) − ZP V E =− χrr − χrs ,
4 r 4 r>s
1X 1X
∆harm = G(0) − ZP V E harm
= χrr + χrs .
4 r 4 r>s
Characterization of Stationary Points

• A stationary point is a geometry q̃ for which the gradient

∂Ee (q̃)/∂ q̃i for all coordinates q̃i : can be a (global or local)
PES minimum, transition state, or higher order saddle
point

• The Hessian Index is the number of negative force constants

(corresponding to imaginary vibrational frequencies, often
printed as negative frequencies)

• For a minimum, verify that there are no imaginary

frequencies

• For a transition state, verify there is exactly one unique

imaginary frequency
Printed by David Sherr

Dec 15, 00 14:54 h2o.freq.out Page 1/5 Dec 15, 00 14:54 h2o.freq.out Page 2/5
There are 5 alpha and 5 beta electrons
*************************************************************** Requested basis set is STO−3G
Running Q−Chem There are 4 shells and 7 basis functions
cicero.chemistry.gatech.edu on Fri Dec 15 14:51:29 EST 2000. A cutoff of 1.0D−10 yielded 10 shell pairs
Version: /usr/local/qchem.6.3.00/bin/qchem There are 34 function pairs
/usr/local/qchem.6.3.00/exe/progman.exe #####################################################
*************************************************************** # Entering gesman.exe on Fri Dec 15 14:51:30 2000 #
Welcome to Q−Chem #####################################################
A Quantum Leap Into The Future Of Chemistry Smallest overlap matrix eigenvalue = 3.63E−01
J. Kong, C. A. White, A. I. Krylov, C. D. Sherrill, Multipole matrices computed through 2nd order
R. D. Adamson, T. R. Furlani, M. S. Lee, A. M. Lee, Guess from superposition of atomic densities
S. R. Gwaltney, T. R. Adams, H. Dachsel, W. M. Zhang, Warning: Energy on first SCF cycle will be non−variational
P. P. Korambath, C. Ochsenfeld, A. T. B. Gilbert, #####################################################
G. S. Kedziora, D. R. Maurice, N. Nair, Y. Shao, # Entering scfman.exe on Fri Dec 15 14:51:30 2000 #
N. A. Besley, P. E. Maslen, J. P. Dombroski, J. Baker, #####################################################
E. F. C. Byrd, T. Van Voorhis, M. Oumi, S. Hirata, A restricted Hartree−Fock SCF calculation will be
C. P. Hsu, N. Ishikawa, J. Florian, A. Warshel, performed using Pulay DIIS extrapolation
B. G. Johnson, P. M. W. Gill, M. Head−Gordon, J. A. Pople, SCF converges when DIIS error is below 1.0E−08
Q−Chem, Version 2.0, Q−Chem, Inc., Export, PA (2000). −−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−
Intel x86 Linux Version Cycle Energy DIIS Error
====================================================================== −−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−
User input: 1 −74.6061773729 4.01E−01 000000
====================================================================== 2 −74.9233301871 6.78E−02 000000
$comment 3 −74.9653425178 8.77E−03 000000
Water 4 −74.9658786194 1.57E−03 000000
$end 5 −74.9659011332 2.76E−05 000000
$molecule 6 −74.9659011457 1.08E−05 000000
0 1 7 −74.9659011480 7.89E−08 000000
O 8 −74.9659011480 2.72E−08 000000
H1 O OH 9 −74.9659011480 1.03E−08 000000
H2 O OH H1 HOH 10 −74.9659011480 3.60E−10 000000 Convergence criterion met
OH = 0.989276 −−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−
HOH = 100.0198 SCF time: CPU 0.08 s wall 0.00 s
$end Final Alpha MO Eigenvalues
$rem 1 2 3 4 5 6
JOBTYPE FREQ 1 −20.2515674 −1.2576205 −0.5939092 −0.4597630 −0.3926280 0.5819312
EXCHANGE HF 7
CORRELATION NONE 1 0.6928062
BASIS STO−3G Final Alpha MO Coefficients
SCF_FINAL_PRINT 2 Energies plus MOS 1 2 3 4 5 6
$end 1 0.9942161 −0.2337592 0.0000000 −0.1040490 0.0000000 0.1258359
====================================================================== 2 0.0258496 0.8443968 −0.0000000 0.5382402 −0.0000000 −0.8203442
Processing $rem in the input. 3 0.0000000 0.0000000 0.6127098 0.0000000 −0.0000000 0.0000000
##################################################### 4 −0.0000000 0.0000000 0.0000000 0.0000000 1.0000000 −0.0000000
# Entering fldman.exe on Fri Dec 15 14:51:30 2000 # 5 −0.0041658 −0.1228911 −0.0000000 0.7558818 −0.0000000 0.7635987
##################################################### 6 −0.0055852 0.1556051 −0.4492190 −0.2950597 0.0000000 0.7692443
−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−− 7 −0.0055852 0.1556051 0.4492190 −0.2950597 0.0000000 0.7692443
Standard Nuclear Orientation (Angstroms) 7
I Atom X Y Z 1 −0.0000000
−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−− 2 0.0000000
1 O 0.000000 0.000000 0.127153 3 0.9598513
2 H −0.757939 0.000000 −0.508611 4 −0.0000000
3 H 0.757939 0.000000 −0.508611 5 −0.0000000
−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−− 6 0.8147663
Molecular Point Group C2v NOp = 4 7 −0.8147663
Largest Abelian Subgroup C2v NOp = 4 Final Alpha density matrix.
Nuclear Repulsion Energy = 8.9077080996 hartrees 1 2 3 4 5 6
1 1.0539352 −0.2276888 0.0000000 −0.0000000 −0.0540635 −0.0112263
2 −0.2276888 1.0033767 −0.0000000 0.0000000 0.3029695 −0.0275650
3 0.0000000 −0.0000000 0.3754133 −0.0000000 −0.0000000 −0.2752409
4 −0.0000000 0.0000000 −0.0000000 1.0000000 −0.0000000 −0.0000000
5 −0.0540635 0.3029695 −0.0000000 −0.0000000 0.5864768 −0.2421295
6 −0.0112263 −0.0275650 −0.2752409 −0.0000000 −0.2421295 0.3131021
7 −0.0112263 −0.0275650 0.2752409 −0.0000000 −0.2421295 −0.0904933
7

Monday April 23, 2001 h2o.freq.out 1/

Printed by David Sherr

Dec 15, 00 14:54 h2o.freq.out Page 3/5 Dec 15, 00 14:54 h2o.freq.out Page 4/5
1 −0.0112263 # Entering drvman.exe on Fri Dec 15 14:51:30 2000 #
2 −0.0275650 #####################################################
3 0.2752409 Calculating MO derivatives via CPHF
4 −0.0000000 1 10 2 0.009892 0.007094
5 −0.2421295 2 12 0 0.000000 0.000000 Roots Converged
6 −0.0904933 Calculating analytic Hessian of the SCF energy
7 0.3131021 Polarizability Matrix (a.u.)
##################################################### 1 2 3
# Entering anlman.exe on Fri Dec 15 14:51:30 2000 # 1 −5.5054487 −0.0000000 −0.0000000
##################################################### 2 −0.0000000 −0.0400454 −0.0000000
Analysis of SCF Wavefunction 3 −0.0000000 −0.0000000 −2.5654448
−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−− Direct stationary perturbation theory relativistic correction:
Orbital Energies (a.u.) and Symmetries rels = 0.031221304489
−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−− relv = −0.096844180832
Alpha MOs, Restricted rel2e = 0.023322614724
−− Occupied −− E_rel = −0.042300261619
−20.252 −1.258 −0.594 −0.460 −0.393 Hessian of the SCF Energy
1 A1 2 A1 1 B1 3 A1 1 B2 1 2 3 4 5 6
−− Virtual −− 1 0.8044023 0.0000000 −0.0000000 −0.4022012 −0.0000000 −0.3374251
0.582 0.693 2 0.0000000 −0.0001352 0.0000000 −0.0000000 0.0000676 −0.0000000
4 A1 2 B1 3 −0.0000000 0.0000000 0.6352335 −0.2165519 −0.0000000 −0.3176168
Beta MOs, Restricted 4 −0.4022012 −0.0000000 −0.2165519 0.4391751 0.0000000 0.2769885
−− Occupied −− 5 −0.0000000 0.0000676 −0.0000000 0.0000000 −0.0000838 0.0000000
−20.252 −1.258 −0.594 −0.460 −0.393 6 −0.3374251 −0.0000000 −0.3176168 0.2769885 0.0000000 0.3002754
1 A1 2 A1 1 B1 3 A1 1 B2 7 −0.4022012 −0.0000000 0.2165519 −0.0369740 −0.0000000 0.0604366
−− Virtual −− 8 0.0000000 0.0000676 −0.0000000 −0.0000000 0.0000162 0.0000000
0.582 0.693 9 0.3374251 −0.0000000 −0.3176168 −0.0604366 −0.0000000 0.0173414
4 A1 2 B1 7 8 9
−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−− 1 −0.4022012 0.0000000 0.3374251
Mulliken Net Atomic Charges 2 −0.0000000 0.0000676 −0.0000000
Atom Charge (a.u.) 3 0.2165519 −0.0000000 −0.3176168
−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−− 4 −0.0369740 −0.0000000 −0.0604366
1 O −0.330636 5 −0.0000000 0.0000162 −0.0000000
2 H 0.165318 6 0.0604366 0.0000000 0.0173414
3 H 0.165318 7 0.4391751 0.0000000 −0.2769885
−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−− 8 0.0000000 −0.0000838 0.0000000
Sum of atomic charges = 0.000000 9 −0.2769885 0.0000000 0.3002754
−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−− Gradient time: CPU 0.59 s wall 1.00 s
Cartesian Multipole Moments #####################################################
−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−− # Entering vibman.exe on Fri Dec 15 14:51:31 2000 #
Charge (ESU x 10^10) #####################################################
0.0000 **********************************************************************
Dipole Moment (Debye) ** **
X 0.0000 Y 0.0000 Z −1.7094 ** VIBRATIONAL ANALYSIS **
Tot 1.7094 ** −−−−−−−−−−−−−−−−−−−− **
Quadrupole Moments (Debye−Ang) ** **
XX −4.4859 XY 0.0000 YY −6.1256 ** VIBRATIONAL FREQUENCIES (CM**−1) AND NORMAL MODES **
XZ 0.0000 YZ 0.0000 ZZ −5.3331 ** INFRARED INTENSITIES (KM/MOL) **
Octapole Moments (Debye−Ang^2) ** **
XXX 0.0000 XXY 0.0000 XYY 0.0000 **********************************************************************
YYY 0.0000 XXZ −0.5313 XYZ 0.0000 Frequency: 2169.95 4141.60 4392.63
YYZ 0.0191 XZZ 0.0000 YZZ 0.0000 IR Active: YES YES YES
ZZZ −0.1747 IR Intens: 7.245 44.303 29.972
Hexadecapole Moments (Debye−Ang^3) Raman Active: YES YES YES
XXXX −6.7321 XXXY 0.0000 XXYY −1.8081 X Y Z X Y Z X Y Z
XYYY 0.0000 YYYY −3.2652 XXXZ 0.0000 O 0.000 0.000 −0.069 0.000 0.000 0.052 −0.068 0.000 0.000
XXYZ 0.0000 XYYZ 0.0000 YYYZ 0.0000 H −0.448 0.000 0.545 −0.570 0.000 −0.416 0.540 0.000 0.453
XXZZ −1.7388 XYZZ 0.0000 YYZZ −1.4596 H 0.448 0.000 0.545 0.570 0.000 −0.416 0.540 0.000 −0.453
XZZZ 0.0000 YZZZ 0.0000 ZZZZ −5.2191 STANDARD THERMODYNAMIC QUANTITIES AT 298.18 K AND 1.00 ATM
−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−− This Molecule has 0 Imaginary Frequencies
#####################################################

Monday April 23, 2001 h2o.freq.out 1/

Printed by David Sherr

Dec 15, 00 14:54 h2o.freq.out Page 5/5

Zero point vibrational energy: 15.302 kcal/mol
Atom 1 Element O Has Mass 15.99491
Atom 2 Element H Has Mass 1.00783
Atom 3 Element H Has Mass 1.00783
Molecular Mass: 18.010570 amu
Principal axes and moments of inertia in atomic units:
1 2 3
Eigenvalues −− 2.58381 4.13508 6.71889
X 1.00000 0.00000 0.00000
Y 0.00000 0.00000 1.00000
Z 0.00000 1.00000 0.00000
Rotational Symmetry Number is 2
The Molecule is an Asymmetric Top
Translational Enthalpy: 0.889 kcal/mol
Rotational Enthalpy: 0.889 kcal/mol
Vibrational Enthalpy: 15.303 kcal/mol
gas constant (RT): 0.593 kcal/mol
Translational Entropy: 34.609 cal/mol.K
Rotational Entropy: 10.672 cal/mol.K
Vibrational Entropy: 0.001 cal/mol.K
Total Enthalpy: 17.673 kcal/mol
Total Entropy: 45.281 cal/mol.K
Archival summary:
1\1\cicero.chemistry.gatech.edu\FREQ\HF\STO−3G\O1H2\sherrill\15Dec2000\0\\#,FREQ
,HF,STO−3G,\\Water\\0,1\O\H,1,0.989276\H,1,0.989276,2,100.02\\HF=−74.9659011\\@
*************************************************************
* *
* Thank you very much for using Q−Chem. Have a nice day. *
* *
*************************************************************
Total job wall time: 3.00 s
***************************************************************
Q−Chem finished
cicero.chemistry.gatech.edu on Fri Dec 15 14:51:32 EST 2000.
***************************************************************

Monday April 23, 2001 h2o.freq.out 1/