Create Presentation
Download Presentation

Download Presentation

Ch121a Atomic Level Simulations of Materials and Molecules

Download Presentation
## Ch121a Atomic Level Simulations of Materials and Molecules

- - - - - - - - - - - - - - - - - - - - - - - - - - - E N D - - - - - - - - - - - - - - - - - - - - - - - - - - -

**Room BI 115**Lecture: Monday, Wednesday Friday 2-3pm Ch121a Atomic Level Simulations of Materials and Molecules Lecture 6 and 7, April 16 and 21, 2013 MD3: vibrations Lecture 6 Presented by Jason Crowley William A. Goddard III, wag@wag.caltech.edu Charles and Mary Ferkel Professor of Chemistry, Materials Science, and Applied Physics, California Institute of Technology TA’s Jason Crowley and Jialiu Wang**Homework and Research Project**First 5 weeks: The homework each week uses generally available computer software implementing the basic methods on applications aimed at exposing the students to understanding how to use atomistic simulations to solve problems. Each calculation requires making decisions on the specific approaches and parameters relevant and how to analyze the results. Midterm: each student submits proposal for a project using the methods of Ch121a to solve a research problem that can be completed in the final 5 weeks. The homework for the last 5weeks is to turn in a one page report on progress with the project The final is a research report describing the calculations and conclusions**Outline of today’s lecture**• Vibration of molecules • Classical and quantum harmonic oscillators • Internal vibrations and normal modes • Rotations and selection rules • Experimentally probing the vibrations • Dipoles and polarizabilities • IR and Raman spectra • Selection rules • Thermodynamics of molecules • Definition of functions • Relationship to normal modes • Deviations from ideal classical behavior**Simple vibrations**Starting with an atom inside a molecule at equilibrium, we can expand its potential energy as a power series. The second order term gives the local spring constant We conceptualize molecular vibrations as coupled quantum mechanical harmonic oscillators (which have constant differences between energy levels) Including Anharmonicity in the interactions, the energy levels become closer with higher energy Some (but not all) of the vibrational modes of molecules interact with or emit photons This provides a spectroscopic fingerprint to characterize the molecule**Vibration in one dimension – Harmonic Oscillator**No friction Consider a one dimensional spring with equilibrium length xe which is fixed at one end with a mass M at the other. If we extend the spring to some new distance x and let go, it will oscillate with some frequency, w, which is related to the M and spring constant k. To determine the relation we solve Newton’s equation M (d2x/dt2) = F = -k (x-xe) Assume x-x0=d = A cos(wt) then –Mw2 Acos(wt) = -k A cos(wt) Hence –Mw2 = -k or w = Sqrt(k/M). Stiffer force constant k higher w and higher M lower w E= ½ k d2**Reduced Mass**M2 M1 Put M1 at R1 andM2 at R2 CM = Center of mass Fix Rcm = (M1R1 + M2R2)/(M1+ M2) = 0 Relative coordinate R=(R2-R1) Then Pcm = (M1+ M2)*Vcm = 0 And P2 = - P1 Thus KE = ½ P12/M1 + ½ P22/M2 = ½ P12/m Where 1/m = (1/M1 + 1/M2) or m = M1M2/(M1+ M2) Is the reduced mass. Get w = sqrt(k/m). Thus we can treat the diatomic molecule as a simple mass on a spring but with a reduced mass, m**For molecules the energy is harmonic near equilibrium but**for large distortions the bond can break. (n + ½)2 Successive vibrational levels are closer by (n + ½)3 (n + ½)2 Exact solution Real potentials are more complex; in general: (Philip Morse a professor at MIT, did not manufacture cigarettes) The simplest case is the Morse Potential:**Vibration for a molecule with N particles**There are 3N degrees of freedom (dof) which we collect together into the 3N vector, Rk where k=1,2..3N The interactions then lead to 3N net forces, Fk = -(∂E(Rnew)/∂Rk) all of which are zero at equilibrium, R0 Now consider that every particle is moved a small amount leading to a 3N distortion vector, (dR)m = Rnew –R0 Expanding the force in a Taylor’s series leads to Fk = -(∂E(Rnew)/∂Rk) = -(∂E/∂Rk)0 - Sm (∂2E/∂Rk∂Rm) (dR)m Where we have neglected terms of order d2. Writing the 2nd derivatives as a matrix (the Hessian) Hkm= (∂2E/∂Rk∂Rm) and setting (∂E/∂Rk)0 = 0, we get Fk = - SmHkm (dR)m = Mk (∂2Rk/∂t2) To find the normal modes we write (dR)m = Amcoswt leading to Mk(∂2Rk/∂t2) = Mkw2 (Akcoswt) = SmHkm (Amcoswt) Here the coefficient of coswt must be {Mkw2 Ak- SmHkm Am}=0 Newton’s equation**Solving for the Vibrational modes**The normal modes satisfy {Mkw2 Ak- SmHkm Am}=0 To solve this we mass weight the coordinates as Bk = sqrt(Mk)Akleading to Sqrt(Mk) w2 Bk- SmHkm [1/sqrt(Mm)]Bm}=0 leading to SmGkmBm = wk2 Bk where Gkm = Hkm/sqrt(MkMm) G is referred to as the reduced Hessian For M degrees of freedom this has M eigenstates SmGkm Bmp = dkpBk (w2)p where the M eigenvalues (w2)p are the squares of the vibrational energies. If the Hessian includes the 6 translation and rotation modes then there will be 6 zero frequency modes, wp= 0**Saddle points**If the point of interest were a saddle point rather than a minimum, G would have one negative eigenvalue, (w2)p = - A2 where A is a positive number This leads to an imaginary frequency, wp= iA , Saddle point**For practical simulations**We can obtain reasonably accurate vibrational modes from just the classical harmonic oscillators, usually within a few % N atoms => 3N degrees of freedom However, there are 3 degrees for translation, n = 0 3 degrees for rotation for non-linear molecules, n = 0 2 degrees if linear The rest are vibrational modes**Normal Modes of Vibration H2O**H2O D2O Sym. stretch 3657 cm-1 2671 cm-1 Ratio: 0.730 1595 cm-1 1178 cm-1 Bend Ratio: 0.735 Antisym. stretch 3756 cm-1 2788 cm-1 Ratio: 0.742 Isotope effect: n ~ sqrt(k/M): Simple nD/nH ~ 1/sqrt(2) = 0.707: More accurately, reduced masses mOH = MHMO/(MH+MO) mOD = MDMO/(MD+MO) Ratio = sqrt[MD(MH+MO)/MH(MD+MO)] ~ sqrt(2*17/1*18) = 0.728 Most accurately MH=1.007825 MD=2.0141 MO=15.99492 Ratio = 0.728**The Infrared (IR) Spectrum**Characteristic vibrational modes • EM energy absorbed by interatomic bonds in organic compounds • frequencies between 4000 and 400 cm-1 (wavenumbers) • Useful for resolving molecular vibrations 13 http://webbook.nist.gov/chemistry/ http://wwwchem.csustan.edu/Tutorials/INFRARED.HTM**Normal Modes of Vibration CH4**4 independent CH bonds 4 CH stretch modes, by symmetry one is triply degenerate 6 possible angle terms HCH 5 HCH modes, one doubly degenate, on triply deg. Reason only 5 linearly independent HCH 1 3 2 3 Sym. stretch Anti. stretch Sym. bend Sym. bend A1 T2 E T2 3019 cm-1 CH4 2917 cm-1 1534 cm-1 1306 cm-1 CD4 2259 cm-1 1092 cm-1 996 cm-1 1178 cm-1**Fitting force fields to Vibrational frequencies and force**constants • Hessian-Biased Force Fields from Combining Theory and Experiment; S. Dasgupta and W. A. Goddard III; J. Chem. Phys. 90, 7207 (1989) MC: Morse bond stretch and cosine angle bend MCX: include 1 center cross terms CH2scis CO stretch CH symstr H2CO CH2 wag CH2 rock CH asymstr 4 atoms 12-6=6 vibrations**The QM Harmonic Oscillator**The Schrödinger equation H = e for harmonic oscillator H energy wavefunctions Gaussian reference http://hyperphysics.phy-astr.gsu.edu/hbase/shm.html#c1**Raman and InfraRed spectroscopy**• IR • Vibrations at same frequency as radiation • To be observable, there must be a finite dipole derivative • Thus homonuclear diatomic molecule (O2 , N2 ,etc.) does not lead to IR absorption or emission. • Raman spectroscopy is complementary to IR spectroscopy. • radiation at some frequency, n, is scattered by the molecule to frequency, n’, shifted observed frequency shifts are related to vibrational modes in the molecule • IR and Raman have symmetry based selection rules that specify active or inactive modes**IR and Raman selection rules for vibrations**The intensity is proportional to dm/dR averaged over the vibrational state The polarizability is responsible for Raman where e is the external electric field at frequency n • For both, we consider transition matrix elements of the form The electrical dipole moment is responsible for IR**IR selection rules, continued**• We see that the transition elements are • The dipole changes during the vibration • Can show that n can only change 1 level at a time For IR, we expand dipole moment**Raman selection rules**substitute the dipole expression for the induced dipole = Same rules except now it’s the polarizability that has to change For both Raman and IR, our expansion of the dipole and alpha shows higher order effects possible For Raman, we expand polarizability**Translation and Rotation Modes**• center of mass translation • Dxa= DxDya=0 Dza=0 • Dxa=0 Dya=DyDza=0 • Dxa=0 Dya=0 Dza=Dz • center of mass rotation (nonlinear molecules) • Dxa=0 Dya=-caDqxDza=baDqx • Dxa= caDqyDya=0 Dza=-aaDqy • Dxa= -baDqzDya=aaDqxDza=0 • linear molecules have only 2 rotational degrees of freedom • The translational and rotational degrees of freedom can be removed beforehand by using internal coordinates or by transforming to a new coordinate system in which these 6 modes are separated out E is a constantdE/dx = 0 d2E/dx2 = 0 Thusthe eigenmode l=0 E is a constantdE/dx = 0 d2E/dx2 = 0 Thusthe eigenmode l=0 21**Classical Rotations**xk(q) is the perpendicular distance to the axis q Can also define a moment of inertia tensor where (just replace the mass density with point masses and the integral with a summation. Diagonalization of this matrix gives the principle moments of inertia! the rotational energy has the form The moment of inertia about an axis q is defined as**Quantum Rotations**For symmetric rotors, two of the moments of inertia are equivalent, combine: Eigenfunctions are spherical harmonic functions YJ,K or Zlm with eigenvalues The rotational Hamiltonian has no associated potential energy**Transition rules for rotations**• For rotations • Wavefunctions are spherical harmonics • Project the dipole and polarizability due to rotation • It can be shown that for IR • Delta J changes by +/- 1 • Delta MJ changes by 0 or +/-1 • Delta K does not change • For Raman • Delta J could be 1 or 2 • Delta K = 0 • But for K=0, delta J cannot be +/- 1**Raman scattering**Phonons are the normal modes of lattice vibrations (thermal + zero point energy) When a photon absorbs/emits a single phonon, momentum and energy conservation the photon gains/loses the energy and the crystal momentum of the phonon. q ~ q` => K = 0 The process is called anti-Stokes for absorption and Stokes for emission. Alternatively, one could look at the process as a Doppler shift in the incident photon caused by a first order Bragg reflection off the phonon with group velocity v = (ω/ k)*k**Raman selection rules**substitute the dipole expression for the induced dipole = Same rules except now it’s the polarizability that has to change For both Raman and IR, our expansion of the dipole and alpha shows higher order effects possible For Raman, we expand polarizability**Another simple way of looking at Raman**Take our earlier expression for the time dependent dipole and expose it to an ideal monochromatic light (electric field) We get the Stokes lines when we add the frequency and the anti-Stokes when we substract The peak of the incident light is called the Rayleigh line**Skip The Sorption lineshape - 1**• The external EM field is monochromatic • Dipole moment of the system • Interaction between the field and the molecules • Probability for a transition from the state i to the state f (the Golden Rule) • Rate of energy loss from the radiation to the system • The flux of the incident radiation c: speed of light n: index of refraction of the medium 28**Skip The Sorption lineshape - II**• Absorption cross section a(w) • Define absorption linshape I(w) as • It is more convenient to express I(w) in the time domain Beer-Lambert law Log(P/P0)=abc I(w) is just the Fourier transform of the autocorrelation function of the dipole moment ensemble average 29**Non idealities and surprising behavior**• Anharmonicity – bonds do eventually dissociate • Coriolis forces • Interaction between vibration and rotation • Inversion doubling • Identical atoms on rotation – need to obey the Pauli Principle • Total wavefunction symmetric for Boson and antisymmetric for Fermion**Electromagnetic Spectrum**How does a Molecule response to an oscillating external electric field (of frequency w)? Absorption of radiation via exciting to a higher energy state ħw ~ (Ef- Ei) Figure taken from Streitwiser & Heathcock, Introduction to Organic Chemistry, Chapter 14, 1976 31**Using the vibrational modes: thermodynamics**In QM and MM the Energy at minima = motionless state at 0K BUT, experiments are made at finite T, hence corrections are required to allow for rotational, translational and vibrational motion. The internal energy of the system: U(T)=Urot(T)+Utran(T)+Uvib(T)+Uvib(0) • From equipartition theorem: Urot(T) = (3/2)KBT , Utran(T) = (3/2)KBT per molecule (except Urot(T)=KBT for linear molecules) • BUT, vibrational energy levels are often only partially excited at room T, thus Uvib(T) requires knowledge of vibrational frequencies • Uvib(T) = vibrational enthalpy @ T - vibrational enthalpy @ 0K Vibrational frequencies can be used to calculate entropies and free energies, or to compare with results of spectroscopic experiments The vibrational frequencies (νi) of the normal modes (Nmod) calculated from the eigenvalues (λi) of the force-constant equivalent of Hessian matrix of second derivatives 32-AJB**Thermodynamics**Describe a system in terms of Hamiltonian H(p,q) where p is generalized momentum and q is generalized coordinate For a system in equilibrium, probability of a state with energy H(p,q) is P(p,q) = exp[-H(p,q)/kBT]/Q which is referred to as a Boltzmann distribution, Here Q, the Partition function, is a normalization constant Q = S exp[-H(p,q)/kBT] summed over all states of the system**The partition function for translation**Assume a cubic periodic box of side L The QM Hamiltonian is The QM eigenfunctions are just periodic functions for x, y, and z directions, sin(nxxp/L) etc Leading to Thus the partition function for translation becomes**Thermodynamic functions for translation**Q = equipartition = (3/2) kT = k = -kT Ideal gas = = (3/2) k**The partition function for rotation**This is 2 the Laplacian This leads to energy levels of I = moment of inertia Thus the partition function becomes**Thermodynamic functions for rotation (non linear)**Q = equipartition = (3/2) kT = k = -kT = 0 equipartition = (3/2) k**The partition function for vibrations**An isolated harmonic oscillator with vibrational frequency ω Has a spectrum of energies Substituting into the Boltzmann expression leads to q = Summing over all normal modes leads to**Thermodynamic functions for vibration (harmonic oscillator)**Q = = (3/2) kT = k = -kT = 0 = k**Thermodynamic functions for electronic stateswe will assume**qelect=1 Q = Assuming the reference state has free atoms = - kT De = k = -De - kT = 0 = 0**Thermodynamic Properties for a Crystal**Write partition function of the system As a continuous superposition of oscillators Harmonic oscillatorPartition function where Thermodynamic properties Weighting functions Reference energy Zero point energy**Where do we get the vibrational density of states DoS(n)?**Experimentally from Inelastic neutron scattering Can use to calculate thermodynamic properties Compare to phonon dispersion curves. Peak is for phonons with little dispersion “Phonon Densities of States and Related Thermodynamic Properties of High Temperature Ceramics” C.-K Loong, J.European Ceramic Society, 1998**Validation of 2PT method for entropy and free energy for**amino acid side chains Tod A Pascal Ravi Abrol William A Goddard III MSC 2013 Research Conference**General Outline**• Overview of the 2PT method for calculating accurate entropies and free energies • Application to the solvation thermodynamics of Amino Acid sidechains**Entropy**“Any method involving the notion of entropy, the very existence of which depends on the second law of thermodynamics, will doubtless seem to many far-fetched, and may repel beginners as obscure and difficult of comprehension.” Rudolf Clausius - originator of the concept of "entropy". --Willard Gibbs, Graphical Methods in the Thermodynamics of Fluids (1873)**Entropy: How to calculate it?**• Test Particle Methods (insertion or deletion) • Good for low density systems • Perturbation Methods • Thermodynamic integration, Thermodynamic perturbation • Applicable to most problems • Require long simulations to maintain “reversibility” • Nonequilibrium Methods (Jarzinski’s equality) • Obtaining differential equilibrium properties from irreversible processes • Require multiple samplings to ensure good statistics • Normal mode Methods • Good for gas and solids • Fast • Not applicable for liquids References: Frenkel, D.; Smit, B. Understanding Molecular Simulation from Algorithms to Applications. Academic press: Ed., New York, 2002. McQuarrie, A. A. Statistical Mechanics. Harper & Row: Ed., New York, 1976. Jarzynski, C. Nonequilibrium Equality for Free Energy Differences. Phys. Rev. Lett. 1997, 78, 2690.**Remaining issue: experimental energies are free energies,**need to calculate entropy General approach to predict Entropy, S, and Free Energy Free Energy, F = U – TS = − kBTlnQ(N,V,T) J. G. Kirkwood. Statistical mechanics of fluid mixtures, J. Chem. Phys., 3:300-313,1935 Kirkwood thermodynamic integration Enormous computational cost required for complete sampling of the thermally relevant configurations of the system often makes this impractical for realistic systems. Additional complexities, choice of the appropriate approximation formalism or somewhat ad-hoc parameterization of the “reaction coordinate”**Solvation free energies amino acid side chains**Pande and Shirts (JCP 122 134508 (2005) Thermodynamic integration leads to accurate differential free energies