CCatomI_Student.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Lab 7 Computational Chemistry - Atoms I"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "In this lab we will explore a variety of techniques used to obtain approximate solutions to the Schrödinger equation for multi-electron systems."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### <font color=green>I. Background Information</font>"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The eigenfunction solution to the Schrodinger equation for a multielectron system depends on the coordinates of all electrons. The orbital approximation says that we can represent a many-electron eigenfunction in terms of individual electron orbitals, each of which depends only on the coordinates of a sinlge electron. A basis set in this context is a set of basis functions used to approximate these orbitals. <br><br>\n",
    "In a strictly mathematical sense, many different kinds of basis set functions could be used. While a variety of different functions have been suggested, two types are commonly used. American physicist J.C. Slater developed equation known as the Slater Type Orbital or (STO) <br><br>\n",
    "<center>$\\text{STO}=\\left(\\frac{\\alpha^3}{\\pi}\\right)^{1/2}e^{-\\alpha r}\\qquad$ (1) </center><br>\n",
    "$\\text{where }\\alpha$ is the “orbital exponent” optimized for a particular atom type.   "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### 1) Why is the STO function a good starting trial function for a basis function? "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Although the STO equation is a good approximation for orbitals, their mathematical evaluation takes a significant amount of computational resources.  Instead, S. F. Boys developed a method of using a combination of Gaussian Type Orbitals in order to express the STO equation. The following is the Gaussian Type Orbital (GTO) equation. <br><br>\n",
    "<center>$\\text{GTO}=\\left(\\frac{2\\alpha}{\\pi}\\right)^{3/4}e^{-\\alpha r^2}\\qquad$ (2)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### 2) Identify the primary functional difference between the STO and GTO functions."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Since the product of the Gaussian \"primitives\" (original Gaussian equations) is another Gaussian, this equation is much easier to use in calculations. Typically multiple Gaussian functions are used to approximate an STO function according to the following summation, \n",
    "<br><br>\n",
    "<center>$\\text{STO-nG}=\\sum_1^n c_n\\text{GTO}_n\\qquad$ (3)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "where the coefficients, $c_n$, specify the weighted contribution of each Gaussian function. The number \"$n$\" before the \"G\" is the number of Gaussian primitives that are used to approximate the STO equation.  So if $n=3$, then the “STO-3G” function will be the sum of 3 Gaussian functions each with their appropriate orbital constant and coefficient.  "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### 3) write out the sum for STO-2G"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "import psi4\n",
    "import numpy as np\n",
    "import matplotlib.pyplot as plt"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### <font color=green>II. Slater and Gaussian type trial wavefunctions for the 1s orbital of hydrogen</font>"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Atomic units assumed throughout. The atomic unit for length is the Bohr radius $a_0$ = 0.529 Å. Thus r = 1 au = 0.529 Å."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "r = np.array([i for i in range(50)])*.1  "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### 3) Plot and compare the STO and STO-1G functions using orbital exponents 1.24 and 0.4166 respectively. How well is an STO function approximated by a single Gaussian function? "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Enter functions"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Plot functions\n",
    "pylab.legend(['STO', 'STO-1G'], loc='upper right')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### 4) Plot the Radial Distribution Function for both the STO and STO-1G functions.  Discuss what type of error (related to the position of the electron) you would expect for the STO-1G based on the form of the wavefunction. "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Enter functions"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Plot functions\n",
    "pylab.legend(['STO', 'STO-1G'], loc='upper right')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Next we will use several Gaussians to model the hydrogen 1s orbital. The numerical constants for a number of different basis sets are stored in a governmental on-line database.\n",
    "\n",
    "Go to the Basis Set Exchange Website: https://www.basissetexchange.org. Select hydrogen from the periodic table and scroll down to the STO-2G basis set and select “Get Basis Set”. The result gives both the journal reference and the orbital constants for a hydrogen 1s orbital.  The first column gives the orbital exponent, $\\alpha$, for each Gaussian function, and the second column gives the corresponding coefficient, $c_n$, used as weighting factors in the sum of the Gaussian functions. "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Enter STO-2G parameters as a list \n",
    "STO2G_a = []   # reference coeff as STO2G_a[0] or STO2G_a[1]\n",
    "STO2G_c = []"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### 5) Create a plot of the following four functions: STO, GTO (1), GTO (2), STO-2G.  Compare and contrast these four functions.  How well is the STO function approximated by each individual Gaussian function?   How well is the STO function approximated by the linear combination of two Gaussian functions?   "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Enter functions"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Plot functions\n",
    "pylab.legend(['STO', 'GTO1', 'GTO2', 'STO-2G'], loc='upper right')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### 6) Obtain the orbital exponents and weighting coefficients for the STO-3G basis set.  Create a plot of the following four functions: STO, STO-1G, STO-2G and STO-3G.  Compare and contrast these four functions.  In comparison to the STO function, which Gaussian approximation is the most accurate?   "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Enter STO-3G parameters as a list \n",
    "STO3G_a = []   \n",
    "STO3G_c = []"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Enter functions"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Plot functions\n",
    "pylab.legend(['STO', 'STO-1G', 'STO-2G', 'STO-3G'], loc='upper right')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### <font color=green>III. Multielectron atom energies for several different basis sets calculated using psi4.</font>"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The following two cells will compute the Hartree Fock Energy (based on the variational method) for the basis set identified. For each basis set specified, be sure to record the energy and number of basis functions used in the calculation."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "basis_set = 'STO-3G' # Change for different Basis Set"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "B = psi4.geometry(\"\"\"\n",
    "0 2\n",
    "B\n",
    "\"\"\")\n",
    "\n",
    "psi4.set_options({'reference': 'uhf','basis': basis_set})\n",
    "e, wfn = psi4.energy('hf',return_wfn=True)\n",
    "\n",
    "print(\"energy =\",e,\" and number of basis functions =\", wfn.basisset().nbf())\n",
    "psi4.core.clean()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "While basis sets of the form STO-nG were popular in the 1980’s, they are not widely used today. For the same reason that multiple Gaussian functions were found to better approximate a Slater type orbital, multiple STO-nG functions are found to better approximate atomic orbitals. For a split-valence basis set, the inner shell (core) electrons are described by a single STO-nG function and only valence electrons are expressed in terms of multiple STO-nG type functions. <br><br>\n",
    "For example, consider the 3-21G basis set. The number before the dash, “3”, represents the 3 Gaussian primitives (i.e. a STO-“3”G) used to represent the inner shell electrons. The next two numbers represent the valence shell split into two sets of STO-nG functions – one with “2” GTOs and one with “1” GTO. "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### 7) Enter the 6-31g basis set for boron. Describe the functions used to represent each orbital and rationalize the number of basis functions in comparison to the STO-3G basis set. (Hint: Consider the number of core vs valence orbitals and the number of STO-nG functions used to respresent each.) "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### 8) Enter the 6-311g basis set for boron. Describe the functions used to represent each orbital and rationalize the total number of basis functions."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "To make an even better approximation of our trial function, we may need to take into account polarization effects. (This accounts for the unequal distribution of electrons when two atoms approach each other.) We can include polarization effects by adding STOs of higher orbital angular momentum, _i.e.,_ d-type functions are added to describe valence electrons in 2p orbitals. The presence of poalrization functions is denoted by an asterisk. Two asterisks indicate polarization functions are added to hydrogens. (This would be important if you wanted to account for hydrogen bonding.) An asterisk in parentheses, ``(*)`` indicates that polarization functions are added only to second-row elements.\n",
    "\n",
    "The presence of diffuse functions is symbolized by the addition of a plus sign to the basis set designator: 6-31+G. Diffuse functions are also useful for systems in an excited stat, for systems with low ionization potential, and systems with some significant negative charge attached. A second +, such as 6-31++G, implies diffuse functions added to hydrogens. The use of double-diffuse basis set is especially useful if you are working with hydrides (H$^-$ ions)."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### 9) Repeat the boron atom energy calculation for each basis set listed below, recording the total HF energy and the number of basis functions used. When approriate, identify the type and number of STO-nG functions used for the core and valence electrons and whether polarized or diffuse functions are included. For the Boron atom, which approximation (choice of basis set) is the most accurate? How does the accuracy relate to the number of basis functions used?"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "| basis set | HF energy | Number of basis functions |\n",
    "|---|---|---|\n",
    "| STO-3G | | |\n",
    "| 6-31g | | |\n",
    "| 6-311g | | |\n",
    "|6-31+G | | |\n",
    "| 6-31G* | | |\n",
    "| 6-31+G* | | |\n",
    "|6-311G* | | |\n",
    "|6-311+G** | | |\n",
    "|cc-pVDZ | | |\n",
    "|cc-pVTZ | | |"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "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.6.10"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}