Chapter09_Part1_finish.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "id": "ce2d6e30",
   "metadata": {},
   "source": [
    "# Chapter 9: Analysis of Chemical Reactions"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "295a0f8f",
   "metadata": {},
   "source": [
    "## 9.1. Calculation of Reaction Thermodynamics"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "7a326b8f",
   "metadata": {},
   "source": [
    "In the following section, we will calculate the changes in enthalpy, entropy and Gibbs free energy of Diels-Alder reactions:"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "4eb5538b",
   "metadata": {},
   "source": [
    "### 9.1.1. Calculation of Thermodynamic State Functions of Reactants and Products"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "id": "7dc680d6",
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "# Import modules\n",
    "import numpy as np\n",
    "import matplotlib.pyplot as plt\n",
    "from rdkit import Chem\n",
    "from rdkit.Chem import AllChem\n",
    "from utils import View3DModel\n",
    "import psi4\n",
    "import py3Dmol"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "id": "5cfdb69a",
   "metadata": {
    "tags": []
   },
   "outputs": [
    {
     "data": {
      "application/3dmoljs_load.v0": "<div id=\"3dmolviewer_17174969940651393\"  style=\"position: relative; width: 800px; height: 400px;\">\n        <p id=\"3dmolwarning_17174969940651393\" style=\"background-color:#ffcccc;color:black\">3Dmol.js failed to load for some reason.  Please check your browser console for error messages.<br></p>\n        </div>\n<script>\n\nvar loadScriptAsync = function(uri){\n  return new Promise((resolve, reject) => {\n    //this is to ignore the existence of requirejs amd\n    var savedexports, savedmodule;\n    if (typeof exports !== 'undefined') savedexports = exports;\n    else exports = {}\n    if (typeof module !== 'undefined') savedmodule = module;\n    else module = {}\n\n    var tag = document.createElement('script');\n    tag.src = uri;\n    tag.async = true;\n    tag.onload = () => {\n        exports = savedexports;\n        module = savedmodule;\n        resolve();\n    };\n  var firstScriptTag = document.getElementsByTagName('script')[0];\n  firstScriptTag.parentNode.insertBefore(tag, firstScriptTag);\n});\n};\n\nif(typeof $3Dmolpromise === 'undefined') {\n$3Dmolpromise = null;\n  $3Dmolpromise = loadScriptAsync('https://cdnjs.cloudflare.com/ajax/libs/3Dmol/2.1.0/3Dmol-min.js');\n}\n\nvar viewer_17174969940651393 = null;\nvar warn = document.getElementById(\"3dmolwarning_17174969940651393\");\nif(warn) {\n    warn.parentNode.removeChild(warn);\n}\n$3Dmolpromise.then(function() {\nviewer_17174969940651393 = $3Dmol.createViewer(document.getElementById(\"3dmolviewer_17174969940651393\"),{backgroundColor:\"white\"});\nviewer_17174969940651393.zoomTo();\n\tviewer_17174969940651393.addModel(\"\\n     RDKit          3D\\n\\n 10  9  0  0  0  0  0  0  0  0999 V2000\\n   -0.4493    2.2322    0.7688 C   0  0  0  0  0  0  0  0  0  0  0  0\\n    0.2433    1.4125    1.5616 C   0  0  0  0  0  0  0  0  0  0  0  0\\n    0.3088   -0.0585    1.4114 C   0  0  0  0  0  0  0  0  0  0  0  0\\n   -0.3158   -0.7591    0.4630 C   0  0  0  0  0  0  0  0  0  0  0  0\\n   -1.0335    1.8809   -0.0684 H   0  0  0  0  0  0  0  0  0  0  0  0\\n   -0.4348    3.2974    0.9605 H   0  0  0  0  0  0  0  0  0  0  0  0\\n    0.8023    1.8504    2.3792 H   0  0  0  0  0  0  0  0  0  0  0  0\\n    0.9113   -0.6015    2.1289 H   0  0  0  0  0  0  0  0  0  0  0  0\\n   -0.2059   -1.8355    0.4359 H   0  0  0  0  0  0  0  0  0  0  0  0\\n   -0.9359   -0.2988   -0.2915 H   0  0  0  0  0  0  0  0  0  0  0  0\\n  1  2  2  0\\n  2  3  1  0\\n  3  4  2  0\\n  1  5  1  0\\n  1  6  1  0\\n  2  7  1  0\\n  3  8  1  0\\n  4  9  1  0\\n  4 10  1  0\\nM  END\\n\",\"molecule\",{\"keepH\": true});\n\tviewer_17174969940651393.setBackgroundColor(\"white\");\n\tviewer_17174969940651393.setStyle({\"stick\": {\"scale\": 0.3}, \"sphere\": {\"scale\": 0.3}});\n\tviewer_17174969940651393.zoomTo();\nviewer_17174969940651393.render();\n});\n</script>",
      "text/html": [
       "<div id=\"3dmolviewer_17174969940651393\"  style=\"position: relative; width: 800px; height: 400px;\">\n",
       "        <p id=\"3dmolwarning_17174969940651393\" style=\"background-color:#ffcccc;color:black\">3Dmol.js failed to load for some reason.  Please check your browser console for error messages.<br></p>\n",
       "        </div>\n",
       "<script>\n",
       "\n",
       "var loadScriptAsync = function(uri){\n",
       "  return new Promise((resolve, reject) => {\n",
       "    //this is to ignore the existence of requirejs amd\n",
       "    var savedexports, savedmodule;\n",
       "    if (typeof exports !== 'undefined') savedexports = exports;\n",
       "    else exports = {}\n",
       "    if (typeof module !== 'undefined') savedmodule = module;\n",
       "    else module = {}\n",
       "\n",
       "    var tag = document.createElement('script');\n",
       "    tag.src = uri;\n",
       "    tag.async = true;\n",
       "    tag.onload = () => {\n",
       "        exports = savedexports;\n",
       "        module = savedmodule;\n",
       "        resolve();\n",
       "    };\n",
       "  var firstScriptTag = document.getElementsByTagName('script')[0];\n",
       "  firstScriptTag.parentNode.insertBefore(tag, firstScriptTag);\n",
       "});\n",
       "};\n",
       "\n",
       "if(typeof $3Dmolpromise === 'undefined') {\n",
       "$3Dmolpromise = null;\n",
       "  $3Dmolpromise = loadScriptAsync('https://cdnjs.cloudflare.com/ajax/libs/3Dmol/2.1.0/3Dmol-min.js');\n",
       "}\n",
       "\n",
       "var viewer_17174969940651393 = null;\n",
       "var warn = document.getElementById(\"3dmolwarning_17174969940651393\");\n",
       "if(warn) {\n",
       "    warn.parentNode.removeChild(warn);\n",
       "}\n",
       "$3Dmolpromise.then(function() {\n",
       "viewer_17174969940651393 = $3Dmol.createViewer(document.getElementById(\"3dmolviewer_17174969940651393\"),{backgroundColor:\"white\"});\n",
       "viewer_17174969940651393.zoomTo();\n",
       "\tviewer_17174969940651393.addModel(\"\\n     RDKit          3D\\n\\n 10  9  0  0  0  0  0  0  0  0999 V2000\\n   -0.4493    2.2322    0.7688 C   0  0  0  0  0  0  0  0  0  0  0  0\\n    0.2433    1.4125    1.5616 C   0  0  0  0  0  0  0  0  0  0  0  0\\n    0.3088   -0.0585    1.4114 C   0  0  0  0  0  0  0  0  0  0  0  0\\n   -0.3158   -0.7591    0.4630 C   0  0  0  0  0  0  0  0  0  0  0  0\\n   -1.0335    1.8809   -0.0684 H   0  0  0  0  0  0  0  0  0  0  0  0\\n   -0.4348    3.2974    0.9605 H   0  0  0  0  0  0  0  0  0  0  0  0\\n    0.8023    1.8504    2.3792 H   0  0  0  0  0  0  0  0  0  0  0  0\\n    0.9113   -0.6015    2.1289 H   0  0  0  0  0  0  0  0  0  0  0  0\\n   -0.2059   -1.8355    0.4359 H   0  0  0  0  0  0  0  0  0  0  0  0\\n   -0.9359   -0.2988   -0.2915 H   0  0  0  0  0  0  0  0  0  0  0  0\\n  1  2  2  0\\n  2  3  1  0\\n  3  4  2  0\\n  1  5  1  0\\n  1  6  1  0\\n  2  7  1  0\\n  3  8  1  0\\n  4  9  1  0\\n  4 10  1  0\\nM  END\\n\",\"molecule\",{\"keepH\": true});\n",
       "\tviewer_17174969940651393.setBackgroundColor(\"white\");\n",
       "\tviewer_17174969940651393.setStyle({\"stick\": {\"scale\": 0.3}, \"sphere\": {\"scale\": 0.3}});\n",
       "\tviewer_17174969940651393.zoomTo();\n",
       "viewer_17174969940651393.render();\n",
       "});\n",
       "</script>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "# Create a molecule of a diene (buta-1,3-diene)\n",
    "diene_mol = Chem.MolFromSmiles('C=CC=C')\n",
    "\n",
    "# Prepare the molecule\n",
    "diene_mol = Chem.AddHs(diene_mol)\n",
    "AllChem.EmbedMolecule(diene_mol, useRandomCoords=True)\n",
    "Chem.rdMolTransforms.SetDihedralDeg(diene_mol.GetConformer(), 0, 1, 2, 3, 0) # Set the C-C-C-C dihedral angle to 0 (s-cis conformation)\n",
    "AllChem.UFFOptimizeMolecule(diene_mol, maxIters=200)\n",
    "\n",
    "# View 3D model with py3Dmol\n",
    "View3DModel(diene_mol)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "id": "a766364e",
   "metadata": {
    "tags": []
   },
   "outputs": [
    {
     "data": {
      "application/3dmoljs_load.v0": "<div id=\"3dmolviewer_17174969940722795\"  style=\"position: relative; width: 800px; height: 400px;\">\n        <p id=\"3dmolwarning_17174969940722795\" style=\"background-color:#ffcccc;color:black\">3Dmol.js failed to load for some reason.  Please check your browser console for error messages.<br></p>\n        </div>\n<script>\n\nvar loadScriptAsync = function(uri){\n  return new Promise((resolve, reject) => {\n    //this is to ignore the existence of requirejs amd\n    var savedexports, savedmodule;\n    if (typeof exports !== 'undefined') savedexports = exports;\n    else exports = {}\n    if (typeof module !== 'undefined') savedmodule = module;\n    else module = {}\n\n    var tag = document.createElement('script');\n    tag.src = uri;\n    tag.async = true;\n    tag.onload = () => {\n        exports = savedexports;\n        module = savedmodule;\n        resolve();\n    };\n  var firstScriptTag = document.getElementsByTagName('script')[0];\n  firstScriptTag.parentNode.insertBefore(tag, firstScriptTag);\n});\n};\n\nif(typeof $3Dmolpromise === 'undefined') {\n$3Dmolpromise = null;\n  $3Dmolpromise = loadScriptAsync('https://cdnjs.cloudflare.com/ajax/libs/3Dmol/2.1.0/3Dmol-min.js');\n}\n\nvar viewer_17174969940722795 = null;\nvar warn = document.getElementById(\"3dmolwarning_17174969940722795\");\nif(warn) {\n    warn.parentNode.removeChild(warn);\n}\n$3Dmolpromise.then(function() {\nviewer_17174969940722795 = $3Dmol.createViewer(document.getElementById(\"3dmolviewer_17174969940722795\"),{backgroundColor:\"white\"});\nviewer_17174969940722795.zoomTo();\n\tviewer_17174969940722795.addModel(\"\\n     RDKit          3D\\n\\n  6  5  0  0  0  0  0  0  0  0999 V2000\\n    0.8740   -0.6170   -2.7614 C   0  0  0  0  0  0  0  0  0  0  0  0\\n   -0.0771    0.2900   -2.9725 C   0  0  0  0  0  0  0  0  0  0  0  0\\n    1.6666   -0.4216   -2.0464 H   0  0  0  0  0  0  0  0  0  0  0  0\\n    0.8625   -1.5573   -3.3030 H   0  0  0  0  0  0  0  0  0  0  0  0\\n   -0.0656    1.2302   -2.4309 H   0  0  0  0  0  0  0  0  0  0  0  0\\n   -0.8697    0.0945   -3.6875 H   0  0  0  0  0  0  0  0  0  0  0  0\\n  1  2  2  0\\n  1  3  1  0\\n  1  4  1  0\\n  2  5  1  0\\n  2  6  1  0\\nM  END\\n\",\"molecule\",{\"keepH\": true});\n\tviewer_17174969940722795.setBackgroundColor(\"white\");\n\tviewer_17174969940722795.setStyle({\"stick\": {\"scale\": 0.3}, \"sphere\": {\"scale\": 0.3}});\n\tviewer_17174969940722795.zoomTo();\nviewer_17174969940722795.render();\n});\n</script>",
      "text/html": [
       "<div id=\"3dmolviewer_17174969940722795\"  style=\"position: relative; width: 800px; height: 400px;\">\n",
       "        <p id=\"3dmolwarning_17174969940722795\" style=\"background-color:#ffcccc;color:black\">3Dmol.js failed to load for some reason.  Please check your browser console for error messages.<br></p>\n",
       "        </div>\n",
       "<script>\n",
       "\n",
       "var loadScriptAsync = function(uri){\n",
       "  return new Promise((resolve, reject) => {\n",
       "    //this is to ignore the existence of requirejs amd\n",
       "    var savedexports, savedmodule;\n",
       "    if (typeof exports !== 'undefined') savedexports = exports;\n",
       "    else exports = {}\n",
       "    if (typeof module !== 'undefined') savedmodule = module;\n",
       "    else module = {}\n",
       "\n",
       "    var tag = document.createElement('script');\n",
       "    tag.src = uri;\n",
       "    tag.async = true;\n",
       "    tag.onload = () => {\n",
       "        exports = savedexports;\n",
       "        module = savedmodule;\n",
       "        resolve();\n",
       "    };\n",
       "  var firstScriptTag = document.getElementsByTagName('script')[0];\n",
       "  firstScriptTag.parentNode.insertBefore(tag, firstScriptTag);\n",
       "});\n",
       "};\n",
       "\n",
       "if(typeof $3Dmolpromise === 'undefined') {\n",
       "$3Dmolpromise = null;\n",
       "  $3Dmolpromise = loadScriptAsync('https://cdnjs.cloudflare.com/ajax/libs/3Dmol/2.1.0/3Dmol-min.js');\n",
       "}\n",
       "\n",
       "var viewer_17174969940722795 = null;\n",
       "var warn = document.getElementById(\"3dmolwarning_17174969940722795\");\n",
       "if(warn) {\n",
       "    warn.parentNode.removeChild(warn);\n",
       "}\n",
       "$3Dmolpromise.then(function() {\n",
       "viewer_17174969940722795 = $3Dmol.createViewer(document.getElementById(\"3dmolviewer_17174969940722795\"),{backgroundColor:\"white\"});\n",
       "viewer_17174969940722795.zoomTo();\n",
       "\tviewer_17174969940722795.addModel(\"\\n     RDKit          3D\\n\\n  6  5  0  0  0  0  0  0  0  0999 V2000\\n    0.8740   -0.6170   -2.7614 C   0  0  0  0  0  0  0  0  0  0  0  0\\n   -0.0771    0.2900   -2.9725 C   0  0  0  0  0  0  0  0  0  0  0  0\\n    1.6666   -0.4216   -2.0464 H   0  0  0  0  0  0  0  0  0  0  0  0\\n    0.8625   -1.5573   -3.3030 H   0  0  0  0  0  0  0  0  0  0  0  0\\n   -0.0656    1.2302   -2.4309 H   0  0  0  0  0  0  0  0  0  0  0  0\\n   -0.8697    0.0945   -3.6875 H   0  0  0  0  0  0  0  0  0  0  0  0\\n  1  2  2  0\\n  1  3  1  0\\n  1  4  1  0\\n  2  5  1  0\\n  2  6  1  0\\nM  END\\n\",\"molecule\",{\"keepH\": true});\n",
       "\tviewer_17174969940722795.setBackgroundColor(\"white\");\n",
       "\tviewer_17174969940722795.setStyle({\"stick\": {\"scale\": 0.3}, \"sphere\": {\"scale\": 0.3}});\n",
       "\tviewer_17174969940722795.zoomTo();\n",
       "viewer_17174969940722795.render();\n",
       "});\n",
       "</script>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "# Create a molecule of a dienophile (ethylene)\n",
    "dienophile_mol = Chem.MolFromSmiles('C=C')\n",
    "\n",
    "# Prepare the molecule\n",
    "dienophile_mol = Chem.AddHs(dienophile_mol)\n",
    "AllChem.EmbedMolecule(dienophile_mol, useRandomCoords=True)\n",
    "AllChem.UFFOptimizeMolecule(dienophile_mol, maxIters=200)\n",
    "\n",
    "# View 3D model with py3Dmol\n",
    "View3DModel(dienophile_mol)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "id": "c01e98c2",
   "metadata": {
    "tags": []
   },
   "outputs": [
    {
     "data": {
      "application/3dmoljs_load.v0": "<div id=\"3dmolviewer_17174969940804062\"  style=\"position: relative; width: 800px; height: 400px;\">\n        <p id=\"3dmolwarning_17174969940804062\" style=\"background-color:#ffcccc;color:black\">3Dmol.js failed to load for some reason.  Please check your browser console for error messages.<br></p>\n        </div>\n<script>\n\nvar loadScriptAsync = function(uri){\n  return new Promise((resolve, reject) => {\n    //this is to ignore the existence of requirejs amd\n    var savedexports, savedmodule;\n    if (typeof exports !== 'undefined') savedexports = exports;\n    else exports = {}\n    if (typeof module !== 'undefined') savedmodule = module;\n    else module = {}\n\n    var tag = document.createElement('script');\n    tag.src = uri;\n    tag.async = true;\n    tag.onload = () => {\n        exports = savedexports;\n        module = savedmodule;\n        resolve();\n    };\n  var firstScriptTag = document.getElementsByTagName('script')[0];\n  firstScriptTag.parentNode.insertBefore(tag, firstScriptTag);\n});\n};\n\nif(typeof $3Dmolpromise === 'undefined') {\n$3Dmolpromise = null;\n  $3Dmolpromise = loadScriptAsync('https://cdnjs.cloudflare.com/ajax/libs/3Dmol/2.1.0/3Dmol-min.js');\n}\n\nvar viewer_17174969940804062 = null;\nvar warn = document.getElementById(\"3dmolwarning_17174969940804062\");\nif(warn) {\n    warn.parentNode.removeChild(warn);\n}\n$3Dmolpromise.then(function() {\nviewer_17174969940804062 = $3Dmol.createViewer(document.getElementById(\"3dmolviewer_17174969940804062\"),{backgroundColor:\"white\"});\nviewer_17174969940804062.zoomTo();\n\tviewer_17174969940804062.addModel(\"\\n     RDKit          3D\\n\\n 16 16  0  0  0  0  0  0  0  0999 V2000\\n   -1.7211   -1.0087   -0.5372 C   0  0  0  0  0  0  0  0  0  0  0  0\\n   -0.6818   -0.7873   -1.6719 C   0  0  0  0  0  0  0  0  0  0  0  0\\n    0.7479   -0.5134   -1.1672 C   0  0  0  0  0  0  0  0  0  0  0  0\\n    1.0897   -1.4325   -0.0270 C   0  0  0  0  0  0  0  0  0  0  0  0\\n    0.1951   -1.6230    0.9497 C   0  0  0  0  0  0  0  0  0  0  0  0\\n   -1.1290   -0.9132    0.8820 C   0  0  0  0  0  0  0  0  0  0  0  0\\n   -2.5428   -0.2662   -0.6334 H   0  0  0  0  0  0  0  0  0  0  0  0\\n   -2.1772   -2.0157   -0.6605 H   0  0  0  0  0  0  0  0  0  0  0  0\\n   -0.6588   -1.6922   -2.3183 H   0  0  0  0  0  0  0  0  0  0  0  0\\n   -1.0043    0.0615   -2.3132 H   0  0  0  0  0  0  0  0  0  0  0  0\\n    0.8158    0.5379   -0.8122 H   0  0  0  0  0  0  0  0  0  0  0  0\\n    1.4771   -0.6479   -1.9944 H   0  0  0  0  0  0  0  0  0  0  0  0\\n    2.0297   -1.9755   -0.0222 H   0  0  0  0  0  0  0  0  0  0  0  0\\n    0.4029   -2.3221    1.7539 H   0  0  0  0  0  0  0  0  0  0  0  0\\n   -1.8333   -1.3531    1.6199 H   0  0  0  0  0  0  0  0  0  0  0  0\\n   -0.9755    0.1563    1.1436 H   0  0  0  0  0  0  0  0  0  0  0  0\\n  1  2  1  0\\n  2  3  1  0\\n  3  4  1  0\\n  4  5  2  0\\n  5  6  1  0\\n  6  1  1  0\\n  1  7  1  0\\n  1  8  1  0\\n  2  9  1  0\\n  2 10  1  0\\n  3 11  1  0\\n  3 12  1  0\\n  4 13  1  0\\n  5 14  1  0\\n  6 15  1  0\\n  6 16  1  0\\nM  END\\n\",\"molecule\",{\"keepH\": true});\n\tviewer_17174969940804062.setBackgroundColor(\"white\");\n\tviewer_17174969940804062.setStyle({\"stick\": {\"scale\": 0.3}, \"sphere\": {\"scale\": 0.3}});\n\tviewer_17174969940804062.zoomTo();\nviewer_17174969940804062.render();\n});\n</script>",
      "text/html": [
       "<div id=\"3dmolviewer_17174969940804062\"  style=\"position: relative; width: 800px; height: 400px;\">\n",
       "        <p id=\"3dmolwarning_17174969940804062\" style=\"background-color:#ffcccc;color:black\">3Dmol.js failed to load for some reason.  Please check your browser console for error messages.<br></p>\n",
       "        </div>\n",
       "<script>\n",
       "\n",
       "var loadScriptAsync = function(uri){\n",
       "  return new Promise((resolve, reject) => {\n",
       "    //this is to ignore the existence of requirejs amd\n",
       "    var savedexports, savedmodule;\n",
       "    if (typeof exports !== 'undefined') savedexports = exports;\n",
       "    else exports = {}\n",
       "    if (typeof module !== 'undefined') savedmodule = module;\n",
       "    else module = {}\n",
       "\n",
       "    var tag = document.createElement('script');\n",
       "    tag.src = uri;\n",
       "    tag.async = true;\n",
       "    tag.onload = () => {\n",
       "        exports = savedexports;\n",
       "        module = savedmodule;\n",
       "        resolve();\n",
       "    };\n",
       "  var firstScriptTag = document.getElementsByTagName('script')[0];\n",
       "  firstScriptTag.parentNode.insertBefore(tag, firstScriptTag);\n",
       "});\n",
       "};\n",
       "\n",
       "if(typeof $3Dmolpromise === 'undefined') {\n",
       "$3Dmolpromise = null;\n",
       "  $3Dmolpromise = loadScriptAsync('https://cdnjs.cloudflare.com/ajax/libs/3Dmol/2.1.0/3Dmol-min.js');\n",
       "}\n",
       "\n",
       "var viewer_17174969940804062 = null;\n",
       "var warn = document.getElementById(\"3dmolwarning_17174969940804062\");\n",
       "if(warn) {\n",
       "    warn.parentNode.removeChild(warn);\n",
       "}\n",
       "$3Dmolpromise.then(function() {\n",
       "viewer_17174969940804062 = $3Dmol.createViewer(document.getElementById(\"3dmolviewer_17174969940804062\"),{backgroundColor:\"white\"});\n",
       "viewer_17174969940804062.zoomTo();\n",
       "\tviewer_17174969940804062.addModel(\"\\n     RDKit          3D\\n\\n 16 16  0  0  0  0  0  0  0  0999 V2000\\n   -1.7211   -1.0087   -0.5372 C   0  0  0  0  0  0  0  0  0  0  0  0\\n   -0.6818   -0.7873   -1.6719 C   0  0  0  0  0  0  0  0  0  0  0  0\\n    0.7479   -0.5134   -1.1672 C   0  0  0  0  0  0  0  0  0  0  0  0\\n    1.0897   -1.4325   -0.0270 C   0  0  0  0  0  0  0  0  0  0  0  0\\n    0.1951   -1.6230    0.9497 C   0  0  0  0  0  0  0  0  0  0  0  0\\n   -1.1290   -0.9132    0.8820 C   0  0  0  0  0  0  0  0  0  0  0  0\\n   -2.5428   -0.2662   -0.6334 H   0  0  0  0  0  0  0  0  0  0  0  0\\n   -2.1772   -2.0157   -0.6605 H   0  0  0  0  0  0  0  0  0  0  0  0\\n   -0.6588   -1.6922   -2.3183 H   0  0  0  0  0  0  0  0  0  0  0  0\\n   -1.0043    0.0615   -2.3132 H   0  0  0  0  0  0  0  0  0  0  0  0\\n    0.8158    0.5379   -0.8122 H   0  0  0  0  0  0  0  0  0  0  0  0\\n    1.4771   -0.6479   -1.9944 H   0  0  0  0  0  0  0  0  0  0  0  0\\n    2.0297   -1.9755   -0.0222 H   0  0  0  0  0  0  0  0  0  0  0  0\\n    0.4029   -2.3221    1.7539 H   0  0  0  0  0  0  0  0  0  0  0  0\\n   -1.8333   -1.3531    1.6199 H   0  0  0  0  0  0  0  0  0  0  0  0\\n   -0.9755    0.1563    1.1436 H   0  0  0  0  0  0  0  0  0  0  0  0\\n  1  2  1  0\\n  2  3  1  0\\n  3  4  1  0\\n  4  5  2  0\\n  5  6  1  0\\n  6  1  1  0\\n  1  7  1  0\\n  1  8  1  0\\n  2  9  1  0\\n  2 10  1  0\\n  3 11  1  0\\n  3 12  1  0\\n  4 13  1  0\\n  5 14  1  0\\n  6 15  1  0\\n  6 16  1  0\\nM  END\\n\",\"molecule\",{\"keepH\": true});\n",
       "\tviewer_17174969940804062.setBackgroundColor(\"white\");\n",
       "\tviewer_17174969940804062.setStyle({\"stick\": {\"scale\": 0.3}, \"sphere\": {\"scale\": 0.3}});\n",
       "\tviewer_17174969940804062.zoomTo();\n",
       "viewer_17174969940804062.render();\n",
       "});\n",
       "</script>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "# Create a molecule of the product\n",
    "product_mol = Chem.MolFromSmiles('C1CCC=CC1')\n",
    "\n",
    "# Prepare the molecule\n",
    "product_mol = Chem.AddHs(product_mol)\n",
    "AllChem.EmbedMolecule(product_mol, useRandomCoords=True)\n",
    "AllChem.UFFOptimizeMolecule(product_mol, maxIters=200)\n",
    "\n",
    "# View 3D model with py3Dmol\n",
    "View3DModel(product_mol)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "id": "c00ee585",
   "metadata": {
    "tags": []
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "17179869184"
      ]
     },
     "execution_count": 5,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# Set the number of threads and memory limit\n",
    "psi4.set_num_threads(16)\n",
    "psi4.set_memory(16*1024*1024*1024) # 16 GB"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "id": "f7c89412",
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "# Set calculation options\n",
    "psi4.set_options({\n",
    "    'basis': '6-31G*',\n",
    "    'scf_type': 'DF',\n",
    "    'reference': 'RHF'\n",
    "})"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "id": "89410102",
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "# Define a function that calculate the thermodynamic properties of molecule\n",
    "def calculate_thermo_properties(mol, temperature, pressure):\n",
    "    # Write the geometry to XYZ string\n",
    "    xyz_string = Chem.MolToXYZBlock(mol)\n",
    "\n",
    "    # Get the psi4 geometry\n",
    "    geometry = psi4.geometry(xyz_string)\n",
    "        \n",
    "    # Set the temperature and pressure\n",
    "    psi4.set_options({\n",
    "        'T': temperature,\n",
    "        'P': pressure * 101325\n",
    "    })\n",
    "    \n",
    "    # Perform geometry optimization\n",
    "    psi4.optimize('b3lyp', molecule=geometry)\n",
    "\n",
    "    # Perform frequency analysis\n",
    "    psi4.frequencies('b3lyp', molecule=geometry)\n",
    "    \n",
    "    # Extract thermodynamic parameters\n",
    "    h = psi4.core.variable('ENTHALPY') * psi4.constants.hartree2kcalmol # Equal to E_h + h_corr\n",
    "    g = psi4.core.variable('GIBBS FREE ENERGY') * psi4.constants.hartree2kcalmol # Equal to E_h + g_corr\n",
    "    s = (h - g) / temperature\n",
    "    return h, s, g"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "id": "6ba836fe",
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "# Set the temperature and pressure\n",
    "temperature = 298.15 # K\n",
    "pressure = 1 # atm"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "id": "4e137353",
   "metadata": {
    "tags": []
   },
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "\tPrevious geometry is closer to target in internal coordinates, so using that one.\n",
      "\n",
      "\tBest geometry has RMS(Delta(q)) = 1.00e-07\n",
      "\n",
      "\tUnable to completely converge to displaced geometry.\n",
      "\tRMS(dx):  2.116e-06 \tMax(dx):  7.045e-06 \tRMS(dq):  2.418e-06\n",
      "\tPrevious geometry is closer to target in internal coordinates, so using that one.\n",
      "\n",
      "\tBest geometry has RMS(Delta(q)) = 1.91e-06\n",
      "\n",
      "\tPrevious geometry is closer to target in internal coordinates, so using that one.\n",
      "\n",
      "\tBest geometry has RMS(Delta(q)) = 3.09e-10\n",
      "\n"
     ]
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Optimizer: Optimization complete!\n",
      "Warning: thermodynamics relations excluded imaginary frequencies: ['124.0203i']\n",
      "Warning: used thermodynamics relations inappropriate for low-frequency modes: ['295.5477' '518.3975' '568.2898']\n",
      "[-97826.11646807921, 0.06536914227027177, -97845.6062778471]\n"
     ]
    }
   ],
   "source": [
    "# Calculate thermodynamic properties of the diene\n",
    "diene_h, diene_s, diene_g = calculate_thermo_properties(diene_mol, temperature, pressure)\n",
    "print([diene_h, diene_s, diene_g])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "id": "b8a8c18f",
   "metadata": {
    "scrolled": true,
    "tags": []
   },
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "\tUnable to completely converge to displaced geometry.\n",
      "\tRMS(dx):  5.843e-08 \tMax(dx):  1.188e-07 \tRMS(dq):  1.462e-07\n",
      "\tPrevious geometry is closer to target in internal coordinates, so using that one.\n",
      "\n",
      "\tBest geometry has RMS(Delta(q)) = 1.46e-07\n",
      "\n",
      "\tUnable to completely converge to displaced geometry.\n",
      "\tRMS(dx):  1.753e-07 \tMax(dx):  3.564e-07 \tRMS(dq):  4.383e-07\n"
     ]
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Optimizer: Optimization complete!\n",
      "[-49279.73318766262, 0.0550811453014369, -49296.15563113424]\n"
     ]
    }
   ],
   "source": [
    "# Calculate thermodynamic properties of the dienophile\n",
    "dienophile_h, dienophile_s, dienophile_g = calculate_thermo_properties(dienophile_mol, temperature, pressure)\n",
    "print([dienophile_h, dienophile_s, dienophile_g])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "id": "35688e6b",
   "metadata": {
    "scrolled": true,
    "tags": []
   },
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "\tPrevious geometry is closer to target in internal coordinates, so using that one.\n",
      "\n",
      "\tBest geometry has RMS(Delta(q)) = 1.99e-04\n",
      "\n",
      "\tPrevious geometry is closer to target in internal coordinates, so using that one.\n",
      "\n",
      "\tBest geometry has RMS(Delta(q)) = 2.25e-06\n",
      "\n",
      "\tPrevious geometry is closer to target in internal coordinates, so using that one.\n",
      "\n",
      "\tBest geometry has RMS(Delta(q)) = 1.04e-06\n",
      "\n",
      "\tPrevious geometry is closer to target in internal coordinates, so using that one.\n",
      "\n",
      "\tBest geometry has RMS(Delta(q)) = 1.21e-07\n",
      "\n",
      "\tPrevious geometry is closer to target in internal coordinates, so using that one.\n",
      "\n",
      "\tBest geometry has RMS(Delta(q)) = 7.21e-08\n",
      "\n"
     ]
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Optimizer: Optimization complete!\n",
      "Warning: thermodynamics relations excluded imaginary frequencies: ['27.9239i']\n",
      "Warning: used thermodynamics relations inappropriate for low-frequency modes: ['185.6355' '419.7076' '446.2481' '534.1510']\n",
      "[-147142.65947350374, 0.07199455859312066, -147164.12465114828]\n"
     ]
    }
   ],
   "source": [
    "# Calculate thermodynamic properties of the product\n",
    "product_h, product_s, product_g = calculate_thermo_properties(product_mol, temperature, pressure)\n",
    "print([product_h, product_s, product_g])"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "0a176762",
   "metadata": {},
   "source": [
    "### 9.1.2. Calculation of Thermodynamic State Functions of the Reaction"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "cbaa9c40",
   "metadata": {},
   "source": [
    "Using Hess's law, we can calculate $\\Delta{H}$, $\\Delta{S}$, and $\\Delta{G}$ of the reaction:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "id": "23636bc3",
   "metadata": {
    "tags": []
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "ΔH = -36.8098 kcal/mol\n",
      "ΔS = -0.0485 kcal/mol\n",
      "ΔG = -22.3627 kcal/mol\n"
     ]
    }
   ],
   "source": [
    "# Calculate the changes in enthalpy, entropy and Gibbs free energy for the reaction\n",
    "dh = product_h - (diene_h + dienophile_h)\n",
    "ds = product_s - (diene_s + dienophile_s)\n",
    "dg = product_g - (diene_g + dienophile_g)\n",
    "\n",
    "print(f'ΔH = {dh:.4f} kcal/mol')\n",
    "print(f'ΔS = {ds:.4f} kcal/mol')\n",
    "print(f'ΔG = {dg:.4f} kcal/mol')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "c60d6963",
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3 (ipykernel)",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.11.7"
  },
  "toc": {
   "base_numbering": "",
   "nav_menu": {},
   "number_sections": false,
   "sideBar": true,
   "skip_h1_title": false,
   "title_cell": "Table of Contents",
   "title_sidebar": "Table of Contents",
   "toc_cell": false,
   "toc_position": {
    "height": "calc(100% - 180px)",
    "left": "10px",
    "top": "150px",
    "width": "203px"
   },
   "toc_section_display": true,
   "toc_window_display": true
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}