week 5 updates

Author: bjmorgan

CH40208/course_contents/week_5.md

@@ -2,7 +2,7 @@
 
 **Geometry Optimisation**
 
-- [Geometry Optimisation](../geometry_optimisation/geometry_optimisation_new.md)
+- [Geometry Optimisation](../geometry_optimisation/geometry_optimisation.md)
     - [Analytical Solution for a Harmonic Potential](../geometry_optimisation/analytical_solution_for_a_harmonic_potential.md)
     - [The Grid Search Method](../geometry_optimisation/grid_search_method.ipynb)
     - [The Gradient Descent Method](../geometry_optimisation/gradient_descent_method.ipynb)

CH40208/geometry_optimisation/newton_raphson_method.ipynb

@@ -76,10 +76,20 @@
     "\n",
     "1. Write a function to calculate the second derivative of a harmonic potential:\n",
     "\n",
-    "$$U^{\\prime\\prime} = 2k$$\n",
+    "$$U^{\\prime\\prime} = k$$\n",
     "\n",
-    "2. Show that the Newton-Raphson method finds the optimal bond length for your harmonic potential in one step, irrespective of your starting point."
+    "2. Show that the Newton-Raphson method finds the optimal bond length for your harmonic potential in one step, irrespective of your starting point.\n",
+    "\n",
+    "*Hint: Consider why the Newton-Raphson method works perfectly for quadratic functions.*"
    ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "id": "55ed4305-f9e0-47c0-b311-c2313f8dea29",
+   "metadata": {},
+   "outputs": [],
+   "source": []
   }
  ],
  "metadata": {

CH40208/geometry_optimisation/scipy_optimize_minimize.ipynb

@@ -149,7 +149,7 @@
    "id": "515697c6",
    "metadata": {},
    "source": [
-    "Here, we are optimisation a one-dimensional function (the energy varies with bond-length), and we get back a one-dimensional solution array. If we were minimising a function over more dimensions (e.g., the energy of a linear triatomic molecule, which then depends on two bond lengths) then the solution array would contain the minimisation solution for each dimension."
+    "Here, we are optimising a one-dimensional function (the energy varies with bond-length), and we get back a one-dimensional solution array. If we were minimising a function over more dimensions (e.g., the energy of a linear triatomic molecule, which then depends on two bond lengths) then the solution array would contain the minimisation solution for each dimension."
    ]
   },
   {
@@ -159,27 +159,103 @@
    "source": [
     "### Exercises\n",
     "\n",
-    "#### 1.\n",
+    "#### Exercise 1: Lennard-Jones Potential with scipy\n",
     "\n",
-    "Use `scipy.optimize.minimize()` to perform a geometry optimisation for your Lennard-Jones potential, again with starting values of $r$ = 3.2 Å, $r$ = 4.4 Å, and $r$ = 6.0 Å.\n",
+    "Use `scipy.optimize.minimize()` to perform a geometry optimisation for the Lennard-Jones potential:\n",
     "\n",
-    "#### 2.\n",
+    "$$U_\\mathrm{LJ} = \\frac{A}{r^{12}} - \\frac{B}{r^6}$$\n",
     "\n",
-    "The potential energy associated with changes in the Cl-C-C-Cl dihedral angle, $\\theta$ in 1,2-dichloroethane can be modelled as\n",
+    "with $A$ = 1 × 10<sup>5</sup> eV Å<sup>12</sup> and $B$ = 40 eV Å<sup>6</sup>.\n",
     "\n",
-    "$$U = \\frac{1}{2}\\left\\{A_1\\left[1+\\cos\\theta\\right] + A_2\\left[1+\\cos2\\theta\\right]+ A_3\\left[1+\\cos3\\theta\\right]\\right\\}$$\n",
+    "**Part A:**\n",
     "\n",
-    "with $A_1$ = 55.229 kJ mol<sup>&minus;1</sup>, $A_2$ = 3.3472 kJ mol<sup>&minus;1</sup>, and $A_3$ = -58.576 kJ mol<sup>&minus;1</sup>.\n",
+    "Write a function `lennard_jones(r, A, B)` to calculate this potential.\n",
     "\n",
-    "- Plot this function for $-\\pi \\leq \\theta \\leq \\pi$.\n",
+    "**Part B:**\n",
     "\n",
-    "- Using `scipy.optimize.minimize()` test the optimisation behaviour for different initial guesses of $\\theta$. In what ways does the minimisation result depend on the choice of initial guess?"
+    "Run `minimize()` with three different starting values: $r$ = 3.2 Å, $r$ = 4.4 Å, and $r$ = 6.0 Å.\n",
+    "\n",
+    "**Part C:** \n",
+    "\n",
+    "For each starting point, print:\n",
+    "- The starting position\n",
+    "- The final optimised position\n",
+    "- The number of iterations required\n",
+    "\n",
+    "**Part D:** \n",
+    "\n",
+    "Do all starting points converge to the same minimum? Why or why not?"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "id": "712334cc-0c94-4df1-af33-7f5506f2a433",
+   "metadata": {},
+   "source": [
+    "#### Exercise 2. Multiple Minima in 1,2-Dichloroethane\n",
+    "\n",
+    "The potential energy associated with changes in the Cl-C-C-Cl dihedral angle, $\\theta$, in 1,2-dichloroethane can be modelled as\n",
+    "\n",
+    "$$U(\\theta) = \\frac{1}{2}\\left\\{A_1\\left[1+\\cos\\theta\\right] + A_2\\left[1+\\cos2\\theta\\right]+ A_3\\left[1+\\cos3\\theta\\right]\\right\\}$$\n",
+    "\n",
+    "with $A_1$ = 55.229 kJ mol$^{-1}$, $A_2$ = 3.3472 kJ mol$^{-1}$, and $A_3$ = -58.576 kJ mol$^{-1}$.\n",
+    "\n",
+    "This exercise explores how the choice of starting geometry affects which minimum the optimiser finds.\n",
+    "\n",
+    "**Part A: Define and Visualise the Potential Energy Surface**\n",
+    "\n",
+    "1. Write a function `dihedral_potential(theta, A1, A2, A3)` that calculates the potential energy for a given dihedral angle $\\theta$ (in radians).\n",
+    "\n",
+    "2. Create an array of $\\theta$ values from $-\\pi$ to $\\pi$ using `np.linspace()` with at least 200 points.\n",
+    "\n",
+    "3. Calculate the potential energy at each point and create a plot showing $U(\\theta)$ vs $\\theta$.\n",
+    "\n",
+    "4. By examining your plot, identify approximately how many minima exist and roughly where they are located.\n",
+    "\n",
+    "**Part B: Test Individual Starting Points**\n",
+    "\n",
+    "5. Use `scipy.optimize.minimize()` to find the minimum, starting with $\\theta_0 = 0$ radians as your initial guess. What minimum does the optimiser find?\n",
+    "\n",
+    "6. Repeat the optimisation starting from $\\theta_0 = \\pi$ radians. Does it find the same minimum?\n",
+    "\n",
+    "7. Try one more starting point of your choice. Based on these three tests, what pattern do you observe?\n",
+    "\n",
+    "**Part C: Systematic Exploration of Starting Points**\n",
+    "\n",
+    "Now we will systematically test many different starting points to understand the overall behaviour.\n",
+    "\n",
+    "8. Create an array of 20 equally spaced starting angles between $-\\pi$ and $\\pi$ using `np.linspace()`. Store this in a variable called `start_angles`.\n",
+    "\n",
+    "9. Create two empty lists called `start_positions` and `final_positions`.\n",
+    "\n",
+    "10. Write a `for` loop that:\n",
+    "    - Iterates through each angle in `start_angles`\n",
+    "    - Performs an optimisation starting from that angle\n",
+    "    - Appends the starting angle to `start_positions`\n",
+    "    - Appends the optimised angle (from `result.x[0]`) to `final_positions`\n",
+    "\n",
+    "**Part D: Visualise the Basins of Attraction**\n",
+    "\n",
+    "11. Create a scatter plot with:\n",
+    "    - Starting angles on the x-axis (`start_positions`)\n",
+    "    - Final optimised angles on the y-axis (`final_positions`)\n",
+    "    - Appropriate axis labels and title\n",
+    "    \n",
+    "    This plot shows you which starting positions lead to which minima. Each horizontal band of points represents a different minimum (basin of attraction).\n",
+    "\n",
+    "**Part E: Analysis**\n",
+    "\n",
+    "12. How many distinct minima does your scatter plot reveal? (Look for horizontal bands of points.)\n",
+    "\n",
+    "13. For each minimum, estimate the range of starting angles that lead to that minimum. These ranges are called the \"basin of attraction\" for each minimum.\n",
+    "\n",
+    "14. Based on your results, if you wanted to find the global minimum (the lowest energy structure), what strategy would you need to use? (Hint: look back at your potential energy plot from Part A to see which minimum has the lowest energy.)\n"
    ]
   },
   {
    "cell_type": "code",
    "execution_count": null,
-   "id": "840a56c1-c19f-49c4-a083-e8576324a098",
+   "id": "b3926a52-9f65-487b-984e-8e8ee670f0ae",
    "metadata": {},
    "outputs": [],
    "source": []