3
What is CanTherm? How is it used? The world’s most compact overview of the theory behind rate theory packages (with emphasis on kinetics) Running CanTherm Complex Pdep Example Calculation, I/O components Outline of this RMG Study Group Objective of this RMG Study Group Provide basic information and conduct a brief overview of topics necessary for computing pressure dependent rates using CanTherm

4
What is CanTherm? CanTherm is an open source python package of utilities for the computation of the following: 1.Thermodynamic properties of stable molecules (H 298, S, C p (T) ) (see Shamel’s study group presentation #5 for more) 2.High pressure limit rate coefficients, k  3.Pressure dependent rate coefficients, k(T,P), for arbitrarily large multiple-well reaction networks using either Modified Strong Collision, Reservoir State or Chemically Significant Eigenvalue (CSE) approximations Notes: CanTherm does not have a GUI There are numerous other similar codes out there, but CanTherm has the nice feature that many molecular properties can be automatically read in from outputs of quantum chemistry jobs If you forked over a copy of RMG-Py from Github, you have CanTherm

6
Electronic Structure and Rates: varying levels of theory Zador et al 2010 Prog. Energy. Combust. Sci. Best practices: always make an attempt to validate or verify the accuracy of your methods, either through comparison with experiments or benchmark calculations

8
Sub-orbital space view: differences between HF, post-HF, and DFT Hartree Fock (HF) theory is a way to variationally estimate the energy of a system of electrons and nuclei, but neglects electron correlation (mean field apprx). post HF methods are advancements of HF that add electron correlation as opposed to simply averaging it out Density Functional Theory (DFT): Computationally faster, scales better with size Focus is on electron density rather than wavefunction Molecular energy is a function of electron density is a function of spacial coordinates (position), hence the name DFT Many DFT methods are semi-empirical (i.e., trained against a experimentally derived dataset) Hybrid or Composite methods: model chemistries involving both HF and DFT components, designed to yield accurate energies at reduced computational costs (e.g., CBS-QB3) Electronic structure calculations only provide geometries, relative energies, force constants, and sometimes, correct point groups necessary for calculation of rates and thermo properties

9
Symmetry Numbers, Point Groups: Important for A-factors and thermo Things to know: 1.Symmetry operations 2.How to identify point groups 3.The rotational symmetry corresponding to various point groups Tips: Flowcharts help. If you can perform basic symmetry operations, you can use a flowchart. Many online resources/tutorials Rotational symmetry reduces a molecule's entropy by a factor of Rln(  ), where  is the rotational symmetry number and R the gas constant. Example: a C 60 Buckminsterfullerene belongs to the I h point group and has a rotational symmetry of 60. Neglecting the rotational contribution to entropy results in an error of over 8 cal/mol-K in an estimation of its standard state entropy. Question. How does the rotational symmetry of cyclohexane change with temperature? Point Group  C1C1 1 CiCi 1 CsCs 1 C 2v 2 CvCv 1 DhDh 2 SmSm m/2 CnCn n C nv n DnDn 2n2n D nh 2n2n D nd 2n2n T12 TdTd OhOh 24 IhIh 60 m = 2, 4, 6,... n = 2, 3, 4,...

10
An effort in futility: statistical mechanics in one slide The canonical partition function (e.g., macroscopic), Q, is summed over all energy levels of a ‘system’ We typically assume that molecular degrees of freedom may be uncoupled: Under the ideal gas assumption, we can rewrite the canonical Partition function as a function of the molecular partition function We use these relations to derive standard thermodynamic properties:

11
Transition state theory gives only the high-pressure limit rate, for most reactions Conventional TST fails for some systems: Barrierless reactions. Must use variational or other methods Systems with many possible transition states

12
RRKM theory is used in the context of the master equation for energy transfer to compute pressure dependence CanTherm counts the density of states using the method of steepest decents, which has been shown to be accurate and faster than direct counting. RRKM rate:

14
rate loss of A i due to reaction rate of collisional production of A at energy level j collisional rate loss of A at energy level i - - Collision rate and frequency:Microcanonical rate constant: The master equation The master equation in chemical kinetics describes the time evolution of a reaction network Consider a reactant, A, with 3N degrees of freedom, depending on the surrounding T and bath gas A is more accurately envisioned as A(E i )

15
The master equation (2) The probability of energy transfer is related to the energy transfer upon collision with bath gas The average downward energy transferred is bath gas (and reactant) dependent and typically a function of temperature Sources: Empirically derived Computed Tuned

16
Collision Frequency, Lennard Jones Parameters The reduced collision integral captures the non-ideality of real colliding molecules by incorporating aspects of the interaction potential between two species. Gas-Kinetic theory is used to compute the collision frequency. Species’ 6-12 Lennard- Jones parameters are needed to compute the reduced collision integral.

17
Online RMG resources make life easier The Joback method is one of corresponding states that relates the critical temperature and pressure of molecules to their LJ-parameters

21
Cantherm input file components – transition states Label ID and location of TS files. Note: no collisional information needed. Reaction cards are needed for each reaction you want to compute the kinetics (one for each TS in your system):

22
kinetics(‘reaction label’): Indicates to CanTherm that you want to compute k  for each of these reactions, which are identified according to labels in the corresponding reaction cards For Pdep reactions, this section is necessary and defines the multiple well reaction network. Include all relevant isomers/wells. The reactant[s] and bath must be included.

23
network label Energy domain discretization rate parametrization: PLOG or Chebyshev Master equation solution method Cantherm input file components – pdep Include External 1D rotor as an active degree of freedom. Specific to assuming that the molecule is a symmetric top with I a  I b  I c By treating it as active, it exchanges energy with other molecular degrees of freedom, convoluted into density of states

25
If all your input parameters are correct, and if CanTherm can accept the level of theory you computed your system at: Look at output files: pdf of reaction network anyFileName.out chem.inp pdfs of 1D rotor potentials and.txts of dihedral angle vs potential energy Run Cantherm. For example, at linux command line: python ~edames/RMG-Py/cantherm.py anyFileName.py

26
Cantherm generates a pdf of your network, which can serve as a good sanity check Make sure your network looks good: No unreasonably large absolute energy values (default units are kJ/mol) All wells are connected as you expect and compare well with your independently created potential energy surface All barriers and relative energies look reasonable compared to your independently performed calculations

34
Hindered rotors Typically can be identified by a vibrational frequencies less than 150 cm -1 Know there are many ways to account for 1-D internal rotors. Cantherm projects out the degree of freedom corresponding to the rotor from the force constant matrix – a good compromise between accuracy and speed. 1-D potential scans typically performed in Gaussian or QChem Care must be taken when preparing cantherm input files If V(  =0  )  0, fourier fit will be inaccurate,  user may ‘shift’ potential to fix this, rather than recompute scan from different starting geometry

35
Hindered rotors not performing thermo calcs so this section is not relevant external rotational symmetry molecular total electronic spin multiplicity (see Shamel’s talk) molecular optical isomers (see Shamel’s talk) Location of Gaussian/QChem output file, and model chemistry used. pivots: two atoms defining axis of rotation top: atoms containing in one of two portions of rotating moiety symmetry: 3 (  CH3), 2 (  CH2), 1 (potato) fit: typically, use ‘best’ Note: atom indices should correspond to those in the geometry file read in by cantherm In this case, I point cantherm to a.txt file for the potential (ScanLog as opposed to GaussianLog or QchemLog)

36
1.Define the reaction network and explore pathways – this can be done using RMG (e.g., via generate reactions); perform a literature search 2.Know what you want to calculate (i.e., relevant T, P) and what you are doing. 3.Conduct quantum chemistry calculations (Gaussian, Qchem, Molpro for CC) at a desired/appropriate level of theory 4.Confirm that your geometries have been optimized properly - look at each structure and ask yourself if the energy is at a minimum - does each saddle point (TS) have one and only one imaginary frequency? - visual inspection via a molecule editor (there are many: GaussView, Avagadro, etc. See http://en.wikipedia.org/wiki/Molecule_editor) Note: avagadro is nice because it can perform isomer searches for you.http://en.wikipedia.org/wiki/Molecule_editor 5.[Very carefully] prepare your CanTherm input files, triple check everything 6.Run CanTherm. 7.Inspect output pdfs: network, 1D HRs 8.Before you use the parametrized rate coefficients in kinetic mechanisms, make sure the fitting errors are acceptable to you, or else consider other options (increase nTemps, use raw output, other fitting methods) Recipe for Reliable Rate Theory Calculations Questions?