Calculating Reaction Energetics for Molecular Systems
Tutorial Created with Software Release: 2024-4
Topics: Catalysis & Reactivity , Energy Capture & Storage , Organic Electronics , Thin Film Processing
Methodology: Molecular Quantum Mechanics
Products Used: Jaguar , MS Maestro , Maestro
|
1.0 MB |
This tutorial is written for use with a 3-button mouse with a scroll wheel.
Words found in the Glossary of Terms are shown like this: Workspacethe 3D display area in the center of the main window, where molecular structures are displayed
Abstract:
In this tutorial, we will demonstrate performing a batch of Jaguar quantum mechanical geometry optimizations in order to predict the thermodynamics of a Diels-Alder reaction. The energetics of the reactions will then be visualized using the Reaction Profile Viewer panel.
Tutorial Content
1. Introduction to Reaction Energetics
Quantum mechanical calculations can enable prediction of reaction energetics. Specifically, quantum mechanical calculations can be used to determine various thermodynamics quantities of interest, including reaction enthalpy (∆H), Gibbs free energy (∆G), and entropy (∆S).
In this tutorial, we will study a prototypical Diels-Alder reaction. We will perform geometry optimization calculations and frequency calculations on each species to predict the thermodynamics of endo versus exo addition using the QM Multistage Workflow panel. See the For Further Reading section for more information about how these quantities are calculated. We will then use the Reaction Profile Viewer panel to easily analyze whether this reaction is thermodynamically or kinetically controlled.
Figure 1: Diels Alder reaction we will explore in this tutorial.
For a basic review on geometry optimizations, see the Introduction to Geometry Optimizations, Functionals and Basis Sets tutorial. While the workflow demonstrated herein is relevant for calculating thermodynamic quantities, it is often of interest to also explore the kinetics of a reaction, which requires transition state searching. Transition state searching for determining activation energy is demonstrated in the Locating Transition States: Part 1 tutorial, but is also summarized here.
Finally, when exploring thermodynamics (or kinetics) with molecular modeling, it is important to note that absolute energies that are compared must be computed with the same functional and basis set combination. If the functional and the basis set are not consistent between two energy calculations, the comparisons will produce meaningless results.
Various tools are available for automating the calculation of thermodynamic quantities for reactions. For more automated procedures, see the Reaction Energetics documentation and the Reaction Network Profiler panel.
In addition to manually defined reactions, automated tools for predicting reactivity are also available for common processes: Bond and Ligand Dissociation, Beta Elimination and pKa (micro-pKa and macro-pKa).
2. Creating Projects and Importing Structures
At the start of the session, change the file path to your chosen Working Directorythe location that files are saved in MS Maestro to make file navigation easier. Each session in MS Maestro begins with a default Scratch Projecta temporary project in which work is not saved, closing a scratch project removes all current work and begins a new scratch project, which is not saved. A MS Maestro project stores all your data and has a .prj extension. A project may contain numerous entries corresponding to imported structures, as well as the output of modeling-related tasks. Once a project is saved, the project is automatically saved each time a change is made.
Structures can be built in MS Maestro or can be imported using File > Import Structures (or drag-and-dropped), and are added to the Entry Lista simplified view of the Project Table that allows you to perform basic operations such as selection and inclusion and Project Tabledisplays the contents of a project and is also an interface for performing operations on selected entries, viewing properties, and organizing structures and data . The Entry Lista simplified view of the Project Table that allows you to perform basic operations such as selection and inclusion is located to the left of the Workspacethe 3D display area in the center of the main window, where molecular structures are displayed. The Project Tabledisplays the contents of a project and is also an interface for performing operations on selected entries, viewing properties, and organizing structures and data can be accessed by Ctrl+T (Cmd+T) or Window > Project Table if you would like to see an expanded view of your project data.
- Double-click the Maestro or Materials Science icon
- (No icon? See Starting Maestro)
Note: This Jaguar workflow can be performed in Maestro or Materials Science Maestro. Use whichever interface you are comfortable with or typically use for your projects.
- Go to File > Change Working Directory
- Find your directory, and click Choose
- Pre-generated input and results files are included for running jobs or examining output. Download the zip file here: schrodinger.com/sites/default/files/s3/release/current/Tutorials/zip/reaction_energy.zip
- After downloading the zip file, unzip the contents in your Working Directorythe location that files are saved for ease of access throughout the tutorial
- Go to File > Save Project As
- Change the File name to reaction_energy_tutorial, click Save
- The project is now named
reaction_profile_tutorial.prj
- The project is now named
For this example, structures of the reactants and products have been built in advance and are available for importing in the next step. These structures were also used in the Locating Transition States: Part 2 tutorial. For more information about how to build and optimize these structures, see the Introduction to Geometry Optimizations, Functionals, and Basis Sets, Introduction to Multistage Quantum Mechanical Workflows, and Locating Transition States: Part 1 tutorials.
If you are interested in building the reactants and products yourself, feel free to do so. For a refresher on using the 2D Sketcher or importing structures, see the Introduction to Maestro for Materials Science tutorial.
- Go to File > Import Structures
- Navigate to where you downloaded the tutorial files, and choose
starting_structures.mae. Click Open- A new entry group titled unoptimized_structures (4) appears in the entry list, containing four entries associated with the reactants and products of our reaction of interest
3. Performing the Geometry Optimization
We will perform the geometry optimization on the two reactants and two possible products with a batch quantum mechanical calculation using the QM Multistage Workflow panel. If you have not used the panel before, an Introduction to Multistage Quantum Mechanical Workflows tutorial is available.
- Select the entire unoptimized_structures (4) entry group from the entry list
- Go to Tasks > Materials > Quantum Mechanics > Molecular Quantum Mechanics > QM Multistage Workflow
- The QM Multistage Workflow panel opens
- For Use structures from, ensure that Project Table (4 selected entries) is chosen
Note: All of the stages and parameters defined in the subsequent steps will apply to all of the input structures, i.e. the panel enables not only multiple QM jobs per molecule but also multiple molecules. This is ideal for energy comparisons, because you can be sure that each molecule was optimized with the same functional and basis set.
Note: The procedure outlined here can also be performed with the Jaguar - Optimization panel.
- For Stage type, choose Optimization
- For Theory (functional), retain B3LYP-D3
- Choose LACVP** for the Basis Set
- These settings match those used to calculate the transition state properties we will use later on in our analysis
- Go to the Properties tab
- Select Vibrational frequencies
- A vibrational frequency calculation is required in order to compute ∆H, ∆G, and ∆S of a reaction
- Go to the Solvation tab
Note: Within the Vibrational frequencies selection, a Thermochemistry section is available below where pressure and temperature can be specified (as well as a ramp of pressure or temperature if you would like to predict thermochemistry at several non-standard states).
In solution the standard state is 1M. To reflect 1 atm >>> 1 M change at 298.15K, one can either adjust the final free energy by 0.00301 Hartree, or calculate free energy under 24.47 atm. One can also use the Thermochemistry Viewer panel for that.
In this tutorial, we will retain standard conditions (1.0 atm and 298.15 K).
- For Solvent model, choose PBF
- Further discussion of solvent models is provided below
- For Solvent, choose water
- Check Optimize in gas-phase
- The structures will also be optimized in the gas-phase. This can be useful on its own (we will be able to analyze the Gas Phase Energy of the reaction) and it also enables the calculation of solvation energy
- Change the Job name to optimization_DA
- Adjust the job settings (
) as needed
- This job requires a CPU host. The job can be completed in about 30 minutes on a 4 CPU host
- If you would like to run the job yourself, click Run. Otherwise, import the pregenerated
Section_03 > optimization_DA > optimization_DA-out.maegzfile - Close the QM Multistage Workflow panel
With respect to the solvent model, we recommend performing geometry optimizations of large, flexible molecules, or molecules containing multiple fragments with the PCM (polarizable continuum model). The PBF (Poisson Boltzmann Finite-element) model occasionally has convergence issues, although it should behave well for smaller, less flexible molecules. In this tutorial, small, inflexible molecules are studied, so the PBF model is a good first choice.
Solvation energies calculated by PCM methods cannot be directly compared to those computed with PBF. This is because the implementation of PCM only covers electrostatic effects and does not include first shell terms that correct for the formation of H-bonds, and other explicit interactions between the solvent and the solute.
For large and flexible molecules, we recommend optimizing geometries with PCM, and if either solvation energy or solution phase energy is important, then perform a single point energy calculation with PBF, acting on the PCM-optimized geometry.
For more information about these models, see the For further reading section.
For a complete description of solvation choices, see the help documentation.
Figure 3-5. Entry list after the job incorporates (or after importing). The optimized diene is shown in the workspace.
When the job is complete (or after importing) a new entry group is incorporated titled optimization_DA-out (4) which contains the four geometry optimized structures as well as the associated vibrations from the frequency calculation.
It is a best practice to visualize the outputs to ensure that the geometries have converged to reasonable structures. It is also best practice to check the vibrational frequencies and ensure that there are no imaginary (listed as negative) frequencies associated with any of the molecules.
4. Analyzing Energetics Using the Reaction Profile Viewer
The geometry optimization and frequency calculation that we performed generate several useful values for quantifying reaction thermodynamics. Here we will look at two ways to parse the data: the Project Tabledisplays the contents of a project and is also an interface for performing operations on selected entries, viewing properties, and organizing structures and data and the Reaction Profile Viewer panels. With the latter approach, we will plot the reaction diagram to visualize the thermodynamics.
A Jaguar calculation generates many energy quantities. We here focus on four values that are most practically useful for understanding reaction thermodynamics:
Gas phase energy: the electronic energy (in Hartree units) of the molecule in the gas phase
Solution phase energy: the electronic energy (in Hartrees) of the molecule surrounded by the solvent continuum. If no solvation is requested in the input, this value will not be returned.
Total enthalpy: the total enthalpy (in Hartrees) of the molecule - includes the electronic energy (gas phase or solution phase), the zero point energy (ZPE), the internal vibrational energy, the internal translational energy, the internal rotational energy, the internal electronic energy and the pV term. Computing these vibrational contributions requires a frequency calculation.
Total Gibbs free energy: the total Gibbs free energy (in Hartrees) of the molecule, which is the sum of enthalpy and of temperature multiplied by entropy. As previously mentioned, enthalpy includes the electronic energy (gas phase or solution phase) and vibration-dependent terms such as ZPE. Computing the entropy and enthalpy thus requires a frequency calculation (which is used to build partition functions and hence the thermodynamic state variables).
Note that quantum chemistry programs, including Jaguar, typically return energies in Hartrees, also known as atomic units (au). Be sure to convert to kcal/mol, kJ/mol or eV as needed. The Reaction Profile Viewer (which we will use momentarily) can help with this.
We will proceed to demonstrate how to locate these values in MS Maestro in the Project Tabledisplays the contents of a project and is also an interface for performing operations on selected entries, viewing properties, and organizing structures and data .
- Open the Project Tabledisplays the contents of a project and is also an interface for performing operations on selected entries, viewing properties, and organizing structures and data
- By default, the Gas Phase Energy and Solution Phase Energy are displayed
- Click on the Property Tree icon (
)
- The Property Tree opens on the side of the Project Table, enabling additional properties to be shown or hidden
As described above, in addition to the gas and solution phase energies, we are also typically interested in the total enthalpy and total Gibbs free energy.
- In the Property Tree, under Jaguar > Secondary, choose Total Enthalpy (au) 298.15 K 1 atm and Total Free Energy (au) 298.15 K 1 atm
- The energy values (again in Hartrees) are added to the Project Tabledisplays the contents of a project and is also an interface for performing operations on selected entries, viewing properties, and organizing structures and data
Optional: Feel free to add other properties to the Project Table if interested. For example, adding the Lowest Frequency is a convenient way to check that you have a minimum-energy structure (the lowest frequency should be > 0 cm-1). Many other properties are available.
In practice, one could now parse the Project Tabledisplays the contents of a project and is also an interface for performing operations on selected entries, viewing properties, and organizing structures and data to compute the relative energies, ∆H and/or ∆G for the reaction. The Property calculator (
) can be used to facilitate energy conversions. However, we now introduce the Reaction Profile Viewer as another way to look at the relative energies and directly provide us with the reaction thermodynamics.
Finally, we will demonstrate using the Reaction Profile Viewer panel to calculate and visualize the reaction thermodynamics:
- Close the Project Tabledisplays the contents of a project and is also an interface for performing operations on selected entries, viewing properties, and organizing structures and data
- Go to File > Import Structures and select
Section_04>modified_TS_DA_endo_ts_search > modified_TS_DA_endo_ts_search.01.maeto import the transition state for the endo-product -
Repeat the previous step for the transition state for the exo-product
- These transition states were calculated in the Locating Transition States: Part 2 tutorial
- Double click on the modified_TS_DA_preopt entry to rename it to modified_TS_DA_exo
- Go to Tasks > Materials > Quantum Mechanics > Molecular Quantum Mechanics > Reaction Profile Viewer
- The Reaction Profile Viewer panel opens
Let’s first look at the reaction to give the endo (kinetic) product.
- With the Reaction Profile Viewer panel open, select(1) the atoms are chosen in the Workspace. These atoms are referred to as "the selection" or "the atom selection". Workspace operations are performed on the selected atoms. (2) The entry is chosen in the Entry List (and Project Table) and the row for the entry is highlighted. Project operations are performed on all selected entries the optimized cyclopentadiene and furan entries from the entry lista simplified view of the Project Table that allows you to perform basic operations such as selection and inclusion
- Next to Reactant, under Import, click Selected Entries
- The diene and dienophile are loaded into the panel. Their Gas Phase Energies are shown, as well as the sum of the energies
- Because the stoichiometry of the reaction is 1:1:1 (reactant:reactant:product), we can retain the Multiplier
- Click Add New Step
- A new step called Product is added to the panel
- Rename this row to Transition State by clicking in the title cell
- Select(1) the atoms are chosen in the Workspace. These atoms are referred to as "the selection" or "the atom selection". Workspace operations are performed on the selected atoms. (2) The entry is chosen in the Entry List (and Project Table) and the row for the entry is highlighted. Project operations are performed on all selected entries the modified_TS_DA_endo entry from the entry lista simplified view of the Project Table that allows you to perform basic operations such as selection and inclusion
- Keep the Reaction Profile Viewer panel open
- Click Selected Entries in the Transition State row in the Import column
- The modified_TS_DA_endo structure is loaded into the panel
- Click Add New Step
- We will leave the title as Product this time
- Select(1) the atoms are chosen in the Workspace. These atoms are referred to as "the selection" or "the atom selection". Workspace operations are performed on the selected atoms. (2) The entry is chosen in the Entry List (and Project Table) and the row for the entry is highlighted. Project operations are performed on all selected entries the optimized DA_kinetic_prod structure from the entry lista simplified view of the Project Table that allows you to perform basic operations such as selection and inclusion
- Next to Product, under Import, click Selected Entries
- The optimized structure for the endo-product is loaded into the panel. The Gas Phase Energy is shown.
- Because the stoichiometry of the reaction is 1:1:1 (reactant:reactant:product), we can retain the Multiplier
We can now plot the reaction profile.
- Next to Jaguar, retain Gas Phase Energy
- In a moment we will look at other thermodynamic quantities
- Next to Reaction Energy units, retain kcal/mol
- eV, kJ/mol and au (Hartrees) are also available
- The Reaction Profile Viewer panel handles the unit conversion
- Click Plot
- The Reaction Profile is displayed in the plot window
Optional: To save an image of the plot, use the Save (
) or Copy (
) icons
- Next to Jaguar, choose Solution Phase Energy
- The table updates to include the Solution Phase Energies
- Note that the transition state structures do not have solution phase energies, so only the reactants and products will be plotted
- Click Plot
- The Reaction Profile is displayed in the plot window
- Next to Jaguar, choose Total Enthalpy 298.15 K 1.00E+00atm
- The table updates to include the Total Enthalpy values
- Click Plot
- The Reaction Profile is displayed in the plot window
- Next to Jaguar, choose Total Free Energy 298.15 K 1.00E+00atm
- The table updates to include the Total Free Energy values
- Next to Reaction Energy units, retain kcal/mol
- Click Plot
- The Reaction Profile is displayed in the plot window
To analyze the reaction thermodynamics for the exo product, repeat the previous steps, but use the modified_TS_DA_exo structure for the transition state and the DA_exo_thermo_prod as the product structure. Alternatively, click Reset (
) in the bottom-left corner of the panel and simply start again importing the second reaction.
Here is a summary of the same four plots for the exo reaction and the previous plots for the endo reaction:
Analysis of the reaction thermodynamics indicates that, as we expect, the exo product is the (slightly) thermodynamically preferred product. It is important to note that the addition of entropy to our calculations had a strong impact on the conclusions drawn: compared to the total enthalpy plots, the free energy plots (which are usually more similar to experimental results) shows an increase in the energetic barrier by 12 kcal/mol and it shows that the reactions are endothermic.
Note: If reporting ∆H or ∆G, be sure to always specify the functional, basis set, solvation model, solvent, temperature and pressure. The thermodynamics quantities depend on all of these factors.
As alluded to above, the Reaction Profile Viewer can handle multi-step reactions as well. Be sure to account for stoichiometry when building out more complex reaction profiles.
While in this tutorial we performed the quantum mechanical calculations in vacuum and then corrected the values to study the reactions at 298.15 K and 1 atm, these parameters can be adjusted in the QM Multistage Workflow panel. Note that reaction thermodynamics at different temperatures and pressures can be obtained without actually running an additional QM calculation. Specifically, use the Thermochemistry Viewer panel to calculate enthalpy and free energy corrections at other temperatures or pressures.
5. Conclusion and References
In this tutorial, we learned how to perform Jaguar quantum mechanical geometry optimizations in order to predict the thermodynamics of a Diels-Alder reaction. Comparison of the reaction profiles confirmed that the exo product is the thermodynamically preferred product at the temperature and pressure of interest.
For further learning:
For introductory content, focused on navigating the Schrödinger Materials Science interface, an Introduction to Materials Science Maestro tutorial is available. Please visit the materials science training website for access to 70+ tutorials. For scientific inquiries or technical troubleshooting, submit a ticket to our Technical Support Scientists at help@schrodinger.com.
For self-paced, asynchronous, online courses in Materials Science modeling, including access to Schrödinger software, please visit the Schrödinger Online Learning portal on our website.
For some related practice, proceed to explore other relevant tutorials:
- Introduction to Geometry Optimizations, Functionals and Basis Sets
- Introduction to Multistage Quantum Mechanical Workflows
- Locating Transition States: Part 1
- Locating Transition States: Part 2
- Optoelectronics
- Band Shape
- Calculating Transition Dipole Moments
- pKa Predictions with Jaguar pKa
- pKa Prediction with Macro-pKa
- Rigid and Relaxed Coordinate Scans
- Computing Atomic Charges
- Optoelectronics Active Learning
- Activation Energies for Reactivity in Solids and on Surfaces
- Design of Asymmetric Catalysts with Reaction Network Enumeration Profiler
For further reading:
- Introduction to Computational Chemistry, 3rd Edition
- Essentials of Computational Chemistry: Theories and Models, 2nd Edition
- Molecular Modelling: Principles and Applications, 2nd Edition
- Best-Practice DFT Protocols for Basic Molecular Computational Chemistry. DOI:10.1002/anie.202205735
- Polarizable Continuum Model. DOI:10.1002/wcms.1086
- The Poisson-Boltzmann model for implicit solvation of electrolyte solutions: Quantum chemical implementation and assessment via Sechenov coefficients. DOI:10.1063/1.5131020
- Explanation for how Jaguar results can be used to calculate free energies. DOI:10.1021/bk-1998-0677.cho22
6. Glossary of Terms
Entry List - a simplified view of the Project Table that allows you to perform basic operations such as selection and inclusion
Included - the entry is represented in the Workspace, the circle in the In column is blue
Project Table - displays the contents of a project and is also an interface for performing operations on selected entries, viewing properties, and organizing structures and data
Recent actions - This is a list of your recent actions, which you can use to reopen a panel, displayed below the Browse row. (Right-click to delete.)
Scratch Project - a temporary project in which work is not saved, closing a scratch project removes all current work and begins a new scratch project
Selected - (1) the atoms are chosen in the Workspace. These atoms are referred to as "the selection" or "the atom selection". Workspace operations are performed on the selected atoms. (2) The entry is chosen in the Entry List (and Project Table) and the row for the entry is highlighted. Project operations are performed on all selected entries
Working Directory - the location that files are saved
Workspace - the 3D display area in the center of the main window, where molecular structures are displayed