Coupled-cluster theory and the discovery of astronomical molecules: A pragmatic perspective

Identification of molecules in space

Radio astronomical detection Very sensitive technique that can allow for unambiguous detection of molecular species Searches are usually based on rotational Hamiltonians based on laboratory data Sensitivity goes with square of the molecular dipole moment First radio detection: OH in 1963; NH 3 , H 2 O, H 2 CO before end of decade (first polyatomics); HCO + in 1970 (first positive ion); C 6 H in 2006 (first negative ion) Strength: Observation of several lines at expected frequencies provides excellent evidence for the existence of particular species in space Weakness: Inherent bias towards polar molecules. Questions regarding the existence of fundamental molecules like benzene, C 60 , etc. will never be answered by radio astronomy.

Infrared detection ISO satellite (European Space Agency) Operating range: 2.5 240 µm (42-4000 cm -1 ) Allows for the observation of non-polar polyatomic molecules and condensed phases (ice, for example) Can (in principle) detect near-IR electronic spectra of radicals and metastables with low-lying excited states Some molecules have extremely intense vibrational transitions A few molecules detected this way: H 2 CO, HF (recent), ice, OH, benzene (?), C 3 , C 5 First detection: CO (1979, Kuiper Airborne Observatory Lockheed C141 flying at ca. 12 km) ; first new molecule: C 2 H 4 CH 4 (1983), then SiH (1991), others.

4 (1984), C 2 H 2 and CO 2 (1989), Strength: Observation of several lines at expected frequencies provides excellent evidence for the existence of particular species in space Weakness: Atmosphere largely opaque to long-wavelength radiation, satellites required for good data. Ultimately a less precise (more uncertain) means of molecular identification than radio astronomy. Identification often based on one (1) line!

The

evidence for benzene in the interstellar medium

Optical detection ISO satellite (European Space Agency) Operating range: 2.5 240 µm (42-4000 cm -1 ) Historically the first technique used: CH discovered in space in 1937. CN and CH+ discovered during WWII. CO and H optically in 1970, N 2 in 2004.

2 seen Can (in principle) detect any molecule; selection rules are not as restrictive; e.g. H 2 and N 2 are IR and MW-inactive species.

Some molecules have extremely intense electronic transitions Strength: Less difficult to obtain data than for IR. Ground-based measurements with CCD cameras now yield now yield high quality spectra.

Weakness: Laboratory data for complex molecules not easy to acquire, some molecules do not have intense electronic transitions; dissociative excited states lead to broad features

The general process of identifying

new

molecules

Obtain precise line positions and molecular constants in the laboratory .

Determine or predict potential candidate line positions Quality (robustness) of identification depends on: PRECISION of observed line positions AGREEMENT with laboratory lines NUMBER of lines used in identification Conduct astronomical search for lines, usually choosing the strongest lines, or those that are in an accessible region of the spectrum.

What can a quantum chemist do to help?

A different sort of synthetic chemistry

Obtain precise line positions and molecular constants in the laboratory .

Determine or predict potential candidate line positions Conduct astronomical search for lines, usually choosing the strongest lines, or those that are in an accessible region of the spectrum.

C 2 H 2 , C 4 H 2 Benzene,Toluene SiH 4 , H 2 S, S vapor N 2

C 6 H 6 + 2000 volts and (unquestionably) many, many, many more: both known and unknown species….

The challenge to quantum chemists:

Guide the laboratory searches for particular compounds The laboratory data is then used in turn to guide the astronomical search.

Theoretical predictions Laboratory identification Astronomical search Astronomy and astrochemistry: Models of interstellar chemistry Abundance, isotopic studies Chemical physics: Molecular properties (geometries) Studies of dynamics (S 4 ), etc .

What

ab initio

theory can provide to laboratory searches:

:

Rotational Spectroscopy:

Molecular geometries (r e ) and associated moments of inertia.

Centrifugal distortion constants, dipole moments.

In principle

: Observable rotational constants, but not really important in general for laboratory astrophysics.

``Adequate’’ levels of theory (i.e. CCSD(T)/cc-pVTZ basis set))

A A o e , B , B o e , C , C o e , CD constants Dipole moment ( Dipole moment (   o e ) ) , CD constants Geometry optimization, harmonic force field calculation, perhaps larger basis set energy calculations.

 Evaluation of anharmonic (cubic or quartic) molecular force field, followed by application of vibrational perturbation theory.

  

Increased treatment of correlation Basis set improvement Increases bond lengths I e B e Decreases bond lengths I e B e

FCI CCSD(T) MP2 B e < B e exact CCSD B e  B e exact B e =B e exact SCF B e > B e exact STO-3G 6-31G* cc-pVDZ cc-pVTZ cc-pVQZ Exact Basis set

In general, CCSD(T)/cc-pVTZ overestimates bond lengths a bit:

B

e calc

< B

e expt However, experiments do not measure B e , but rather B v

Usually positive

Molecular vibration effects tend to extend the structure (r g > r e ) B 0 < B e

Exact

In general, CCSD(T)/cc-pVTZ overestimates bond lengths a bit:

B

e calc

< B

e expt However, experiments do not measure B e , but rather B v

Usually positive

Molecular vibration effects tend to extend the structure (r g > r e ) B 0 < B e

Exact CCSD(T)/cc-pVTZ

In general, CCSD(T)/cc-pVTZ overestimates bond lengths a bit:

B

e calc

< B

e expt However, experiments do not measure B e , but rather B v

Usually positive

Molecular vibration effects tend to extend the structure (r g > r e ) B 0 < B e

Exact CCSD(T)/cc-pVTZ

An empirical observation : CCSD(T)/cc-pVTZ B e constants usually within 1% of experimental B 0 constants, and are in better agreement than the theoretical B 0 values Serve as a general-purpose guide To laboratory searches for new molecules A more or less trivial application of quantum chemistry

An example: Isomers of C 5 H 2 Theory

:

An example: Isomers of C 5 H 2 Theory

: ``

eiffelene

’’

However … heroic calculations can act to support important new claims… McCarthy, Gottlieb, Gupta and Thaddeus Ap. J. 652, L141 (2006).

Laboratory spectra of C 4 H and C 8 H recently recorded (Gupta

et al

., preprint)

First molecular anion found in space (carrier of U1377 in IRC+10216)

Vibrational Spectroscopy:

With high quality atomic natural orbital basis sets, CCSD(T) and a good treatment of anharmonicity, it is now possible to predict fundamental and two-quantum level positions and intensities relatively accurately for rigid molecules (frequencies within 10 cm-1, intensities to a few percent). Laboratory spectra can be assigned by this means, a field that is well-ploughed and certainly not limited to laboratory astrophysics.

An example: Overtone spectrum of nitric acid from 4000-6000 cm -1

These calculations require evaluation of

: the harmonic, cubic and quartic force fields the first three derivatives of the dipole moment components vibrational second-order perturbation theory for levels and intensities

A practical point and an advertisement:

Combined with CCSD(T) and VPT2, the ANO basis sets of Taylor and Almlöf work extremely well for vibrational level positions.

But, alas…

While nice, it remains a fact that infrared spectroscopy is not used as often as microwave spectroscopy as an analytical technique to detect new molecules in the astrochemistry community. Exceptions, however, exist: Maier, McMahon, others. Lack of precision in laboratory data (relative to MW) makes it less attractive.

Electronic Spectroscopy:

Precise prediction of band positions is extraordinarily difficult (impossible) here. Accuracy of 0.1 eV (800 cm-1) is in fact rare. Using calculated band positions to guide laboratory searches is a dubious proposition. Moreover, the number of observable vibronic features much less than that typical of a microwave spectrum.

And with regard to the relationship between astrophysics and quantum chemistry, things are largely inverted:

Laboratory lines Astronomical lines Molecular identity

MW Elec

Electronic Spectroscopy:

Precise prediction of band positions is extraordinarily difficult (impossible) here. Accuracy of 0.1 eV (800 cm-1) is in fact rare. Using calculated band positions to guide laboratory searches is a dubious proposition. Moreover, the number of observable vibronic features much less than that typical of a microwave spectrum.

And with regard to the relationship between astrophysics and quantum chemistry, things are largely inverted:

Laboratory lines Astronomical lines Molecular identity

MW Elec

The diffuse interstellar bands (DIBS)

First published report of (two) DIBs in 1921 Found along lines of sight associated with reddened stars Catalog of Herbig (1975): 39 bands Catalog of Jenniskens and Desert (1994): >200 bands

Source(s) of DIBS remain(s) unknown in late 2006!

Some ideas from the past:

Molecular carriers first suspected (Struve, Russell, Saha,Swings) ca. 1937-1944 … but population of polyatomics in typical interstellar clouds thought to be too small. Led to movement towards… Solid-state carriers lots of ideas: impurities in dust grains, solid oxygen, etc.

… but lack of characteristic solid-state absorption features (notably inhomogeneous broadening) led to movement towards Molecules again!

An excellent review: G.H. Herbig Ann. Rev. Astrophys. 33, 19 (1995)

Bis-pyridyl magnesium tetrabenzoporphyrin

(in part)

Postulated to account for ALL diffuse interstellar bands when just 39 were known (1972) Crashed and burned in the 1970’s: Occam’s razor, spectra not reproducible.

Some molecular ideas:

Molecular hydrogen, albeit not simply

} } }

Gaps between intermediate state and higher levels fit many DIBs positions

BUT

Model required photon fluxes of ca. 10 10 s -1 Flux in typical ISM cloud: ca. 1 year -1 … insufficient by 17 orders of magnitude!

Polycyclic Aromatic Hydrocarbons (PAHs)

Polycyclic Aromatic Hydrocarbons (PAHs)

The Unidentified Infrared Bands (UIRs)

NGC 7027     1 .

3 % Orion Bar

Wavelength ( µ)

ISO spectra of the Planet. Neb. NGC 7027 and the PDR region at the Orion Bar

Peeters et al., 2003

The unidentified infrared (UIR) features

4435 4395

Carrier of



4429?

Corannulene, C

20

H

10 

exp B

= 2.07(2)D = 509.8MHz

Lovas et al., 2005

However, there has not been a single unambiguous detection of any PAH in space to date, despite considerable publicity.

Nonpolarity is a problem, although two smallest polar PAH species have not been found

Apart from benzene, there is no evidence of anything other than three-membered ring compounds in the ISM

No hard evidence supports PAH carriers of DIBs

Something quite remarkable transpired in 1998… Tulej, Kirkwood, Pachkov and Maier , Ap.J.Lett. 506, 69 (1998)

Something quite remarkable transpired in 1998… C 7 !?

Tulej, Kirkwood, Pachkov and Maier , Ap.J.Lett. 506, 69 (1998)

Since this time, sobered experimentalists have not yet advanced another candidate An interesting puzzle, though…

Molecule produced in a benzene discharge, also with toluene No REMPI signal (suggests high IP - cation?)

Evidence from isotopic substitution studies

: Similar geometry in lower and upper states Five hydrogen atoms, two pairs of symmetry equivalent hydrogens

What molecule is this?

An interesting twist: No longer are we trying to identify the carrier of a line in the ISM, but the carrier of a line measured in a basement laboratory in Massachusetts

Has not proven to be easy Quantum chemistry not the best approach (we’ll discuss this shortly), only suggestion for carrier thus far seems a rather unlikely molecule.

Still not solved definitively - moreover, carrier of  4429 in the laboratory probably is not carrier of astronomical line

Is my favorite molecule (NO

3

) a DIBs carrier?

6616 6235 6278 X B Absorption spectrum of NO 3 (a very strong transition) From DIBS catalog

:

Lines at 6234.3, 6278.3 and 6613.8

v. 6235 6278 6616

What can quantum chemistry do to help?

Calculate line positions?

Most techniques have accuracies ca. 0.1-0.5 eV

TDDFT CASPT, EOM-CCSD MR-CISD, EOM-CCSDT MR-AQCC

Applicability Accuracy

To be useful as a “screener’’ for DIBS candidates requires accuracy of ca. 100 cm -1

What can quantum chemistry do to help?

Calculate line positions?

Most techniques have accuracies ca. 0.1-0.5 eV

TDDFT CASPT, EOM-CCSD MR-CISD, EOM-CCSDT MR-AQCC

Applicability Accuracy

To be useful as a “screener’’ for DIBS candidates requires accuracy of ca. 100 cm -1

Calculate absorption profiles?

 Experimental spectrum Calculated spectrum for candidate molecule 

Calculate absorption profiles?

 Experimental spectrum Calculated spectrum for candidate molecule Assignment 

Calculate absorption profiles?

Experimental spectrum  Calculated spectrum for candidate molecule Assignment 

This sort of assignment can often provide clear picture of the nature of vibronic states

What is required in the simulation of electronic spectra?

The simplest case: The Franck-Condon approximation

Applicable to most transitions from one “isolated” electronic state to another, particularly in low-energy part of the spectrum.

Assumptions

: 

(r,R) =



(r;R)



(R) Vibronic Electronic Wavefunction Wavefunction Vibrational Wavefunction [calculated from V eff (R)] <



’’(r,R) |



|



’(r,R) > = <



’’(r,R 0 ) |



|



’(r,R 0 ) > “M ij ” ’ Dipole operator Reference geometry

Vibronic level positions Electronic energy difference (adiabatic)

Vibronic level positions Electronic energy difference - ZPE of ground state

Vibronic level positions Electronic energy difference - ZPE of ground state + vibrational energy in final state

Vibronic level positions Electronic energy difference + ZPE of ground state + vibrational energy in final state

Origin

Vibronic level positions Electronic energy difference + ZPE of ground state + vibrational energy in final state

Origin

Intensities Dipole approximation: Intensity

  |

<



’’(r,R) |



M ij 2 F k |



’(r,R)>| 2 Transition Franck-Condon moment factor (FCF) |<



’’(R)|



’(R)>| 2

Vibronic level positions Electronic energy difference + ZPE of ground state + vibrational energy in final state …

Origin

Intensities Dipole approximation: Intensity

  |

<



’’(r,R) |



M ij 2 F k |



’(r,R)>| 2 Transition Franck-Condon moment factor (FCF) |<



’’(R)|



’(R)>| 2

Origin

With several modes …

Origin 1 2 3 0 1 2 3 4 5 6 Origin

…

To do a Franck-Condon simulation, we need:

Quantum chemical calculation of ground state geometry and force field “EASY” - CAN BE DONE WITH ANY METHOD Quantum chemical calculation of excited state geometry and force field MORE CARE REQUIRED WITH RESPECT TO CHOICE OF METHOD (CIS, RPA, CIS(D), MCSCF, CASPT2, EOM-CC, MRCI) (balance is important here) Quantum chemical calculation of transition dipole moment NOT SO HARD - ACCURACY USUALLY NOT VERY IMPORTANT

Franck-Condon simulations account for:

Progressions in totally symmetric vibrations (provides a measure of the geometry change due to excitation) Even-quantum transitions in nonsymmetric vibrations (shows up only if there is an appreciable force constant change)

… but do not account for:

Final states that are not totally symmetric (due to “vibronic” coupling) Spectroscopic manifestations of non-adiabaticity (BO breakdown) (effects of conical intersections, avoided crossings etc.)

Low-lying states usually heavily affected by vibronic coupling!

Vibronic effects on potential energy surfaces

(Model due to Köppel, Domcke and Cederbaum) A model potential for a two coupled-mode, two state system



1



2 = 0.4 eV (323 cm -1 ) [symmetric] = 0.1 eV (807 cm -1 ) [non-symmetric ] Parameters:

  

A



B Vertical energy gap between the states - Linear coupling constant between the two states - Slope of diabatic potential of state A along q 1 - Slope of diabatic potential of state A along q 1 at q 1 =q 2 =0 at q 1 =q 2 =0 Hamiltonian corresponding to model (in diabatic basis) T 1 0



A q 1 + 1/2 [



1 q 1 2 +



2 q 2 2 ]



q 2 H = ( ) + ( ) 0 T 2



q 2



+



B q 1 + 1/2 [



1 q 1 2 +



2 q 2 2 ]

T + V

Vibronic effects on potential energy surfaces

Diagonalization of potential energy ( V ) gives adiabatic potential energy surfaces *



= 0 eV q 2

* Note that the associated diagonal basis does not necessarily block-diagonalize the Hamiltonian

q 2



= 0.05 eV

q 2



= 0.10 eV

q 2



= 0.15 eV

q 2



= 0.20 eV

Vibronic effects on potential energy surfaces

Example:



= 0.5 eV;



1 = 0.2 eV;



2 = -0.2 eV;



= variable Diabatic surfaces

A state B state

Lowest adiabatic surface (what we calculate in quantum chemistry)

Equivalent minima 

= 0.2 eV

Pseudorotation transition state Conical Intersection

Note that complicated (but realistic) adiabatic surfaces arise from an extremely simple model potential

More visualization…

TOP VIEW PERSPECTIVE VIEW

Lowest adiabatic surface with different coupling strengths

kk 

=0.05 eV



=0.10 eV



=0.15 eV



=0.20 eV



=0.25 eV



=0.30eV

Lowest adiabatic surface with different coupling strengths

kk 

=0.05 eV



=0.10 eV



=0.15 eV



=0.20 eV



=0.25 eV



=0.30eV

In this case, both states are minima on the potential energy surface. But note that two different electronic states lie on the same potential energy surface!!!

Far too infrequently thought about in quantum chemistry. Real world example: 2 A 1 states of NO 2 and 2 B 2 “A state” “B state”

An astrophysical example: Propadienylidene

1 A 1 1 B 1 1 A 2

}

1 A 1 Coupled by modes of b 2 symmetry

Weak vibronic coupling between dark state and bright state 2 B 1 state (minimum) 2 A 2 state (minimum) Conical intersection Lowest adiabatic potential sheet

Vibronic effects on “vibrational” energy levels

Diagonalization of complete Hamiltonian ( T+V ) gives vibronic energy levels Model potential again:



= 0.20 eV;



= 0.5 eV;



1 = 0.2 eV;



2 = -0.2 eV Harmonic A state frequencies from diabatic potential : 806 cm -1 (



s ) ; 323 cm -1 (



a ) Harmonic A state frequencies from adiabatic potential 806 cm -1 (



a ) ; 253 cm -1 (



a ) Exact vibronic levels below 1050 cm-1 250 cm -1 (



n ) , 500 cm -1 ( 2



n ) , 752 cm -1 ( 3



n ) , 802 cm -1 (



s ) ,1004 cm -1 ( 4



n ) , 1043 cm -1 (



s +



n )

Vibronic effects on “vibrational” energy levels

Diagonalization of complete Hamiltonian ( T+V ) gives vibronic energy levels Model potential again:



= 0.20 eV;



= 0.5 eV;



1 = 0.2 eV;



2 = -0.2 eV Harmonic A state frequencies from diabatic potential : 806 cm -1 (



s ) ; 323 cm -1 (



a ) Harmonic A state frequencies from adiabatic potential 806 cm -1 (



a ) ; 107i cm -1 (



a ) Exact vibronic levels below 1050 cm-1 0, 121, 327, 590, 720, 866, 888, 1046 (symmetric levels) 19, 213, 453, 723, 736, 931, 1045 (nonsymmetric levels)

Vibronic effects on potential energy surfaces

Diagonalization of potential energy ( V ) gives adiabatic potential energy surfaces *



n 2



n



s 3



n

Add some slides here with wavefunctions, discussion of nodes, etc.

Wavefunctions and stationary state energies

Eigenstates of system obtained by diagonalizing Hamiltonian Given by Vibrational basis functions



=



A



c i



i +



B



c i



i Diabatic electronic states

Wavefunctions and stationary state energies

Eigenstates of system obtained by diagonalizing Hamiltonian Given by Vibrational basis functions



=



A



c i



i +



B



c i



i true only if



=0 Diabatic electronic states



=



A



c i



I or



=



B



c i



i “Vibrational level of electronic state A” “Vibrational level of electronic state A” Electronic states are coupled by the off-diagonal matrix element “Breakdown of the Born-Oppenheimer Approximation”

Add more slides here with wavefunctions, densities, projection onto diabatic states, etc.

The Calculation of Electronic Spectra including Vibronic Coupling Energies given by Eigenvalues of model Hamiltonian

‘

Relative intensities given by

|

<



’’(r,R) |



|



’(r,R)>| 2 Ground state Final state

Model system: Adiabatic perspective

Corresponding Hamiltonian very complicated: Potential energy matrix is diagonal but not simple (discontinuities); kinetic energy matrix is clearly not diagonal; transition dipole moment very sensitive wrt geometry

Model system: Adiabatic perspective

Corresponding Hamiltonian very complicated: Potential energy matrix is diagonal but not simple (discontinuities); kinetic energy matrix is clearly not diagonal; transition dipole moment very sensitive wrt geometry Green arrow transition to “bright state

”

Model system: Adiabatic perspective

Corresponding Hamiltonian very complicated: Potential energy matrix is diagonal but not simple (discontinuities); kinetic energy matrix is clearly not diagonal; transition dipole moment very sensitive wrt geometry Red arrow transition to “dark state” Green arrow transition to “bright state

”

An aside:

Traditional quantum chemistry assumes: T A

)

0 V A 0 T B

) )

0 0 V B

)

Adiabatic potential energy surfaces Vibrational energy levels calculated from the Schrödinger equations (T A + V A )



(T B + V B )



= E vib



= E vib



and total (vibronic energies) given by: E ev (A) = (V A ) min + E vib E ev (B) = (V B ) min + E vib

Diabatic perspective (KDC Hamiltonian) conceptually (and computationally a much simpler approach

T 1 0



A q 1 + 1/2 [



1 q 1 2 +



2 q 2 2 ]



q 2 H = ( ) + ( ) 0 T 2



q 2



+



B q 1 + 1/2 [



1 q 1 2 +



2 q 2 2 ] Treatment: 1. Assume initial state not coupled to final states (not necessary, but a simple place to start)



’’ =



0



000… 2. Assume transition moments between diabatic states are constant <



0 |



<



0 |



| |

 

A B > = M > = M A B 3. Diagonalize Hamiltonian (Lanczos recursion is best choice)



’ = c A 0



A



000… + c B 0



B



000… +



’[ c A i



A



i + c B i



B



i ]i

T 1 0



A q 1 + 1/2 [



1 q 1 2 +



2 q 2 2 ]



q 2 H = ( ) + ( ) 0 T 2



q 2



+



B q 1 + 1/2 [



1 q 1 2 +



2 q 2 2 ] 4. Stick spectra given by [c A 0 + c B 0 ] 2



(E -



) Basis set and symmetry considerations Direct product basis



A



00,



A



01,



A



02 …



B



00,



B



01,



B



02 … Symmetry of vibronic level



ve

=



v

x



e

T 1 0



A q 1 + 1/2 [



1 q 1 2 +



2 q 2 2 ]



q 2 H = ( ) + ( ) 0 T 2



q 2



+



B q 1 + 1/2 [



1 q 1 2 +



2 q 2 2 ] 4. Stick spectra given by Makes entire contribution to intensity only ONE element of eigenvector matters [c A 0 + c B 0 ] 2



(E -



) Basis set and symmetry considerations Direct product basis



A



00 ,



A



01 ,



A



02



B



00 ,



B



01 ,



B



02 … … Symmetry of vibronic level



ve

=



v

x



e

Appearance of eigenvectors - pictorial view

 S  N  A Franck-Condon Vibronic coupling (weaker)  S  N  B

“vibronically allowed level”

Appearance of eigenvectors - pictorial view

 S  N  A Franck-Condon Vibronic coupling (stronger)  S  N  B

“vibronically allowed level”

Franck-Condon Vibronic coupling