From c75c1293dbbf864e0f3762af275794798cf08b00 Mon Sep 17 00:00:00 2001
From: Andrew Cumming <andrew.cumming@mcgill.ca>
Date: Thu, 14 Sep 2023 11:44:16 -0400
Subject: [PATCH] Updates integration.md

---
 integration.md                | 62 +++++++++++++++++++++++++++++++++++
 interpolation_solutions.ipynb |  2 +-
 2 files changed, 63 insertions(+), 1 deletion(-)

diff --git a/integration.md b/integration.md
index 9ea70c6..edecc66 100644
--- a/integration.md
+++ b/integration.md
@@ -1,2 +1,64 @@
 # Integration
 
+## Newton-Cotes methods
+
+A related problem to interpolation is when we need to calculate the integral of a function 
+
+$$I = \int_{x_1}^{x_N} f(x) dx$$
+
+given values of the function at a discrete set of points $x_1\dots x_N$, which we'll write as $f_i\equiv f(x_i)$. For simplicity here we'll assume that the spacing between points on the grid is constant $x_{i+1}-x_i = \Delta x$, but it's straighforward to generalize to non-uniform sampling if you need to.
+
+To estimate the value of the integral, we need to model how the function behaves in each interval $x_i<x<x_{i+1}$. As with interpolation, we can take $f(x)$ in each interval to be a polynomial, including additional terms in the polynomial to increase the accuracy of our approximation. 
+
+**Rectangle rule**
+
+The first approximation is to assume that $f(x)$ is a constant $f(x)=f_i$ in the interval $x_i\leq x\leq x_{i+1}$. The total integral is then
+
+$$I\approx \sum_{i=1}^{N-1} f_i \Delta x.$$
+
+
+**Trapezoidal rule**
+
+If we make a linear approximation, we can write 
+
+$$f(x) \approx f_i + {f_{i+1}-f_i\over \Delta x} \left(x-x_i\right)\hspace{1cm} x_i\leq x\leq x_{i+1}.$$
+
+Then 
+
+$$\int_{x_i}^{x_{i+1}} f(x) dx =  \int_0^{\Delta x} f(x_i+y) dy \approx f_i \Delta x + (f_{i+1}-f_i){\Delta x\over 2}  = {f_i+f_{i+1}\over 2}\Delta x.$$
+
+This is known as the *trapezoidal rule* because it corresponds to the area of the trapezoid formed by connecting the points $(x_i,f_i)$ and $(x_{i+1}, f_{i+1})$ by a straight line.
+
+When we sum over all intervals, each point gets counted twice except for the left and right boundaries, so
+
+$$I \approx \Delta x \left[{f_1\over 2} + f_2 + f_3 \dots + f_{N-1} + {f_N\over 2}\right].$$
+
+
+**Simpson's rule**
+
+For the next order, we consider the double interval $x_{i-1}<x<x_{i+1}$, and write
+ 
+$$f(x) \approx f_i + (x-x_i) \left.{df\over dx}\right|_{x_i}+ {(x-x_i)^2\over 2}\left.{d^2f\over dx^2}\right|_{x_i}.$$
+
+The reason for considering the double interval from $x_{i-1}$ to $x_{i+1}$ is that the linear term is odd in this interval (antisymmetric about $x_i$) and so does not contribute to the integral, i.e.
+
+$$\int_{x_{i-1}}^{x_{i+1}} f(x) dx \approx 2 f_i \Delta x + {1\over 2} \left.{d^2f\over dx^2}\right|_{x_i} {2(\Delta x)^3\over 3}.$$
+
+Taking 
+
+$$\left.{d^2f\over dx^2}\right|_{x_i} \approx {f_{i+1}-2f_i + f_{i-1}\over (\Delta x)^2}$$
+
+and simplifying then gives
+
+$$\int_{x_{i-1}}^{x_{i+1}} f(x) dx \approx {\Delta x\over 3} \left[f_{i-1} + 4f_i + f_{i+1} \right].$$
+
+Adding up the contributions from each double interval, the total integral is
+
+$$I \approx {\Delta x\over 3}\left[ f_1 + 4f_2 + 2f_3 + 4f_4 + 2f_5 \dots 4f_{N-1} + f_N \right],$$
+
+where we need the total number of points $N$ to be an odd number (so we can divide the domain into a set of double intervals).
+
+
+```{admonition} Exercise: Newton-Cotes
+Implement these three methods to integrate a function $f(x)$ of your choice. How does the error in each method scale with the number of points $N$? What kind of functions are integrated exactly (to machine precision) for each method?
+```
diff --git a/interpolation_solutions.ipynb b/interpolation_solutions.ipynb
index 5b55807..9cec74f 100644
--- a/interpolation_solutions.ipynb
+++ b/interpolation_solutions.ipynb
@@ -157,7 +157,7 @@
     "Things to note:\n",
     "\n",
     "- For $N$ sample points, the error in linear interpolation decreases as $N^2$, cubic and spline decrease as $N^4$.  This is because the largest correction term (next term in the Taylor expansion) is $\\propto (\\Delta x)^2$ for linear and $\\propto (\\Delta x)^4$ for cubic. \n",
-    "- The full $N-1$-order polynomial fit (Lagrange polynomial) does not behave well, often showing high frequency components, particularly near the boundaries of the domain.\n",
+    "- The full $N-1$-order polynomial fit (Lagrange polynomial) does not behave well, often showing high frequency components, particularly near the boundaries of the domain ([Runge's phenomenon](https://en.wikipedia.org/wiki/Runge's_phenomenon)).\n",
     "- The cubic fits show discontinities in the derivatives at the sampled points; the spline (by definition) has continuous first and second derivatives.\n",
     "- The cubic and spline fits develop oscillations near discontinous changes in the function of derivatives.\n",
     "- The cubic and spline interpolations are able to fit polynomials up to 3rd order exactly. The Lagrange fit will fit up to order $N-1$ exactly. "