WIXLIB

Quantum Mechanics in BiologyTodd J. Martinez

Quantum Biology Is quantum mechanics necessary for biology? Yes, but mostly for “light” particles Electrons Force Fields Bond-Rearrangement Electron Transfer Nuclei Tunneling – Proton Transfer Multiple electronic states – Photobiology

Cytochrome c Oxidase –ET/PT/Bond RearrangementTerminal Enzyme inRespiratory ChainO24e-MembraneMembrane4H 2H2O4H

Proton Transfer

Bacteriorhodopsin – Light-InducedProton PumpH H hνMembraneMembraneH Need QM to describe excited stateAnd bond-rearrangement associated with H pump

Force Fields – The Building Blockof Biomolecular SimulationsGV ( R) i bonds 12ki ,bond ( ri ,bond rieq,bond ) i ji , j atomsqi q jri rj 12 i ji , j atoms i angles(ki ,angle (θi ,angle θieq,angle ) "VLJij ri rj)But where does this come from? In reality,GGGG GGH electronic ( R)ψ i ,electronic (relectronic ; R) Ei ,electronic ( R)ψ i ,electronic (relectronic ; R)Electronic Schrodinger Equation

Electronic HamiltonianZAZ AZB1ˆˆH electronic T (i ) iiA riAi , j rijA, B RABKinetic Energyof onrepulsionNuclear-nuclearrepulsion

Ab Initio Quantum Chemistry The Good Well-defined hierarchy – in principle always know routeto improve results Prescriptions for thermochemistry with kcal/molaccuracy exist (but may not always be practical) Excited electronic states without special treatment The Bad Periodic boundary conditions are difficult Can be computationally costly; even “showcase”calculations on 200 atoms are rare

Quantum Chemical “Canon” Two-pronged HierarchyElectron CorrelationMinimal Basis Set Full CI“Right Answer”Minimal Basis Set/Hartree-FockComplete Basis Set/Hartree-FockBasis set

The Never-Ending ContractionχAO ,CBFk cOne-particlebasis set MOφiNkp 1 Pr im CBFpkχ AO , primitivep( r, R )N CBFMOAO ,CBFcRχ ik ( ) k ( r , R )k 1 MOψ µ ( r , R ) A φi ( r , R ) i Molecular orbitals areorthogonal contractionsof AOsAntisymmetrized productsof MOsMany-particleN APΨ elec ( r , R ) cµCI ( R )ψ µ ( r , R )Basis setµ 1Every atomic orbitalis a fixed contractionof GaussiansTotal electronic wfn iscontraction of APs

Basis Sets (One-Particle) Centered on atoms – this means we need fewer functionsbecause geometry of molecule is embedded in basis set Ideally, exponentially-decaying. This is the form of Hatom solutions and is also the correct decay behaviorfor the density of a molecule. But then integrals areintractable This is the reason for the fixed contractions ofGaussians – try to mimic exponential decay and cuspwith l.c. of GaussiansAdding Basis Functions:Reeves and Harrison, JCP 39 11 (1963)Bardo and Ruedenberg, JCP 59 5956 (1973)Schmidt and Ruedenberg, JCP 71 3951 (1979)

Gaussians vs. Plane WavesAtom-centered Places basis functions in the important regions Gradient of energy with respect to atom coordinateswill be complicated (need derivatives of basisfunctions) Linear dependence could be a problem Localized – Good for reducing scaling Plane Waves Force periodic description (could be good) Gradients are trivial Need many more basis functions Required integrals are easier

Basis Set ClassificationMinimal Basis Set (MBS)One CBF per occupied orbital on an atomE.g., H has one s function, C has 2s and 1pn-zetan CBF per occupied orbital on an atomValence n-zetaMBS for core (1s of C), n-zeta for valencePolarizedAdd higher angular momentum functions thanMBS – e.g., d functions on CDiffuse or augmentedAdd much wider functions to describe weaklybound electrons and/or Rydberg states

Physical Interpretation Could just say more functions more complete, but thisgives no insight n-zeta:csmall clargeAllows orbitals to “breathe,”i.e. to change their radial extent

Physical Interpretation IIPolarization functions:cs cpIt should be clear thatextra valence andpolarization functionswill be most importantwhen bonds are stretchedor atoms are overcoordinatedExample for H atom; generallypolarization functions alloworbitals to “bend”

Alphabet Soup of Basis SetsAfter 30 years, only a handful of basis sets still used: STO-3G – The last MBS standing “Pople-style” – m-n1 nXGX-zetam # prim in core ni # prim in ith valence AO3-21G – Pathologically good geometries for closedshell molecules w/HF (cancellation of errors)6-31G, 6-31G*, 6-31G**, 6-31G , 6-31G * polarization on non-H ** polarization on all diffuse on non-H diffuse on all cc-pvXz, aug-cc-pvXz – X-zeta - “correlation-consistent”best, but tend to be larger than Pople sets

Hartree-Fock Truncating the many-particle basis set at one term givesHartree-FockΨ HF NCBF MO AO ,CBF A cki ( R ) χ k( r , R ) i k 1 Can be shown that this implies a nonlinear effectiveone-particle problem MOMOMO ˆˆˆˆF ( c ) h (i ) 2J j ( c ) K j ( c ) i j occ GGGHˆ R ψ R Eψ R Fˆ {ϕ ( r )} ϕi ( r ) ε iϕi ( r )( ) ( )( )()

Self-Consistent FieldGuess solution (cMO)Build Fock MatrixSolve eigenvalue equation Fc EcIf coefficients are stil changing

“Static” CorrelationConsider HF wavefunction at dissociation for H2:φσ χ left χ rightMOs: φ * χ χleftrightσψ HF A (φσ φσ (αβ βα ) )Infinite separationεσ*εσFinite RH-HorExpand in AOs:?ψ HF χ l χ l χ r χ r χ l χ r χ r χ lNeed more than one determinant!

Restricted vs. UnrestrictedCan solve the previous problem by allowing orbitals tobe singly occupied (unrestricted HF)ψ UHF A (φσ φσ′ (αβ βα ) )Problem: This is not a spin eigenfunctionSˆ 2ψ UHF S ( S 1)ψ UHFWhy didn’t we write:ψ UHF A (φσ′ φσ (αβ βα ) )?In fact, pure spin state is l.c. of the two ψ singlet ψ UHF ψ UHFψ triplet ψ UHF ψ UHF

Describing CorrelationEasiest Way: Moller-Plesset Perturbation Theory (MPn)N elHˆ 0 Hˆ Fˆ ( i )i 1Series diverges for stretched bonds!?!Only first correction (MP2) is worthwhileMore stable: configuration interaction (CI)Solve for CI coefficients variationallyψ CIcreation/annihilationoperators CI CICI†† † c0 aa ai cia aa ab ai a j cijab " ψ SCFi ,aijab truncated atsome excitationlevel (FCI no truncation)may be HF or multi-determinant

Multi-Determinant HF (MCSCF)HF solves only for cMO – Add cCI and solve for both“Active Space” – the set of orbitals where electronicoccupation variese.g. for H2:φσ φσ ; φσ φσ ;φσ φσ***CASSCF – “Complete” active space – allrearrangements of electrons allowed within activespace

Size ConsistencyE(AN) for A infinitely separated should be NE(A) This simple requirement is not met by truncated CI. E should be additive for noninteracting systems ψ should be a productExponential maps products to sums Alternative (Coupled Cluster):ψ CC e CCCC aa† ai cia aa† ab† ai a j cijab " i ,a ijab ψ0When exponential ansatz is expanded, find contributionsfrom excitations up to all orders 1 kcal/mol accuracy possible, but can fail for bond-breakingbecause there are no good multi-reference versions

Density Functional TheoryIs there another way?DFT replaces the wavefunction with charge densityas the fundamental unknownρ ( r1 ) dr2 " drnψ * ( r1 " rn )ψ ( r1 " rn )Wavefunction – 3n coordinatesCharge Density – 3 coordinatesDFT can be better than HF. How can this be?

DFT – FunctionalsDFT expression for the energy:E[ ρ ] T [ ρ ] ρVnuclei ρ (r ) ρ (r ') r r' drdr ' K XC [ ρ ]e- e- repulsionKinetic energye-/nuclei attractionExchange / Correlation[] denotes functional – take function and return a numberFor example, a definite integral is a type of functional bI [ f ] f (r )dra

So How Can This Work?KXC is UNKNOWN!! (And is unlikely to ever beknown in a form which is simpler than solving theelectronic Schrodinger equation)T is also unknown, but can be approximated if thedensity is associated with a wavefunction.Kohn-Sham procedure:ψ KS Aˆ φiρ (φi )2iT [ ρ ] ψ KS Tˆ ψ KS

DFT and HF Need to define KXC Exactly the same ansatz is used as HF – the onlydifference is in the Fockian operator MOMOMO ˆˆˆˆFHF ( c ) h ( i ) 2 J j ( c ) K j ( c ) i j occ MOMOMO ˆˆˆˆFKS ( c ) h ( i ) 2 J j ( c ) aK K j ( c ) Kˆ xc [ ρ , ρ ]i j occ Same SCF procedure as in HF since the equation isnonlinear

Local Density Approximation (LDA)KXC is known numerically for homogeneous gas ofelectronsAssume density is slowly varying:hom ogenous ,exactK[ ρ ]drXC K XC [ ρ ] drAssume constantdensity in each regionProblem: Errors are large (up to 30kcal/mol)

Gradient CorrectionsPiecewise-linear approximation to densityExact results not known; hence there are several“gradient-corrected” functionalsKXC KXC [ρ, ρ]Examples: BLYP, PW91Much improved approximation, but errors can still be aslarge as 10 kcal/mol

Hybrid FunctionalsThe Coulomb interaction we wrote counts theinteraction of electrons with themselvesIn Hartree-Fock, this is exactly canceled by exchangeintegralsTry adding in some Hartree-Fock exchangeB3LYP is most popular functional of this typeErrors go down to 3-5 kcal/mol in most casesCost still roughly same as HF

Behavior of HF and DFT By definition, HF has no electron correlationAs we saw earlier, this implies more serious errorsfor stretched/distorted bonds, i.e. disfavorsovercoordination Pure DFT overestimates correlationPreference for overcoordination Hence success of hybrid functionals which add exchangeto DFT, e.g. B3LYP Hartree-Fock alone is not very useful – barriers are usuallyoverestimated by more than DFT underestimates

Problems with DFT Is DFT a panacea? No! Even the best DFT often yield errors of 5 kcal/mol No hierarchy for improvement Different functionals Different answers Poor for proton transfer and bond rearrangment Tendency to overcoordinate Extreme example: LDA predicts no protontransfer barrier in malonaldehydeinstead of No satisfactory route to excited electronic states

Semiempirical MethodsBasic approximation: φµ1 ( r1 )φν 2 ( r1 )φη 3 ( r2 )φτ 4 ( r2 )r1 r2 δ µν δητAtomic indices for basis functions Hartree-Fock type of SCF using this (and related)integral approximations Problem: Need to parameterize remaining integrals tomodel correlation Many variants (MNDO, AM1, PM3)

Semiempirical MethodsAdvantages Cheaper than DFTOnly truly viable QM-like methods for entireproteins, but even small proteins are barely withinreachCan be reparameterized for each system/processDisadvantages H-bond strengths often wrong by several kcal/molStill expensive

Summary of MethodsRHFUHFCASSCFCICCMP2DFTVar?Multi SizeApprox ErrorRef? Consistent? in 10 kcal/molbarrier heightYYYYNNNNNYYNNNNYNearlyOnly Full-CIYYY/N5-155-153-71-50.1-34-101-5N.B. There are multi-reference perturbation and CC theories, esp.CASPT2 has been successful but sometimes has technical problems

PES TopographyTransition StateConicalIntersectionGlobalMinimumLocal Minima

Important Points Normally, only look for stationary pointsG E ( R )G 0 R RGstationary These geometries may be local minima, global minima,transition states or higher order saddle points How to check? Build and diagonalize the “Hessian” matrix 2 E2 R1 % 2 R N E R1 2E R1 RN 2E RN2 Count negative eigenvalues0 local minimum1 saddle point 1 useless

Hessian MatrixGenerally built in Cartesian coordinates Will have 6 zero eigenvalues corresponding to rotationand translationThese must be identifiedand ignored in the analysisHow to identify? Animatenormal modes, e.g. withMolDenDisadvantage – Expensive(10x Energy Calculation)

Special Warning! When a molecule has symmetry beware of optimizing tosaddle points! If you enforce symmetry, obviously will maintainsymmetry But, just starting from a high symmetry geometry isenough, because symmetry requires that gradient isnonzero only with respect to totally-symmetricmodes Example: Try optimizing the geometry of water startingwith perfectly linear molecule for initial guess Conclusions: Avoid high symmetry starting points Always verify that stationary points are minima, atleast by perturbing geometry (but Hessian is best)