{
 "cells": [
  {
   "cell_type": "markdown",
   "id": "signal-contrary",
   "metadata": {},
   "source": [
    "# Centered Difference Approximation of the Derivative\n",
    "\n",
    "**Last Revised on Monday March 29**"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "significant-shopper",
   "metadata": {},
   "source": [
    "## Programming Exercise\n",
    "\n",
    "Write a (Python) function to approximate derivatives using the centered difference formula,\n",
    "\n",
    "$$\n",
    "Df(x) \\approx \\delta_hf(x) := \\frac{f(x+h) - f(x-h)}{2h}\n",
    "$$\n",
    "\n",
    "The input variables should include the node spacing $h$; think about and discuss what all the input and output variables should be.\n",
    "(Also, remember that such a function should not do any interactive input, or any output to the screen or files.)\n",
    "\n",
    "Then — as usual — add code that tests this on some examples such as $f(x) = e^x$, $f'(1)$ as suggested in the section on [Test Cases for Differentiation](test-cases-differentiation.ipynb).\n",
    "\n",
    "**Update:** For test cases like the above where the exact result can be calculated, use this to compute and display the error."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "lightweight-treaty",
   "metadata": {},
   "source": [
    "**Note:** For the above function and every function below, create a Jupyter notebook to run test cases: this will handle any output to the screen."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "damaged-phoenix",
   "metadata": {},
   "source": [
    "### Warm-up exercise: forward difference approximation\n",
    "\n",
    "Let us do this for the forward difference approximation\n",
    "\n",
    "$$\n",
    "Df(x) \\approx D_h(x) := \\frac{f(x+h) - f(x)}{h}\n",
    "$$\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "id": "federal-formation",
   "metadata": {},
   "outputs": [],
   "source": [
    "import numpy as np"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "id": "instructional-salad",
   "metadata": {},
   "outputs": [],
   "source": [
    "def D_h(f, x, h):\n",
    "    D_h_of_f_at_x = (f(x+h) - f(x))/h\n",
    "    return D_h_of_f_at_x"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "id": "respiratory-winter",
   "metadata": {},
   "outputs": [],
   "source": [
    "def f(x): return np.exp(x)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "id": "comic-blood",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "\n",
      "For x=0.0:\n",
      "\n",
      "For h=0.1:\n",
      "D_hf=1.0517091807564771\n",
      "error=0.051709180756477124\n",
      "\n",
      "For h=0.01:\n",
      "D_hf=1.005016708416795\n",
      "error=0.005016708416794913\n",
      "\n",
      "For h=0.001:\n",
      "D_hf=1.0005001667083846\n",
      "error=0.0005001667083845973\n",
      "\n",
      "For x=1.0:\n",
      "\n",
      "For h=0.1:\n",
      "D_hf=2.858841954873883\n",
      "error=0.14056012641483795\n",
      "\n",
      "For h=0.01:\n",
      "D_hf=2.7319186557871245\n",
      "error=0.01363682732807936\n",
      "\n",
      "For h=0.001:\n",
      "D_hf=2.7196414225332255\n",
      "error=0.0013595940741804036\n",
      "\n",
      "For x=3.0:\n",
      "\n",
      "For h=0.1:\n",
      "D_hf=21.12414358253968\n",
      "error=1.0386066593520127\n",
      "\n",
      "For h=0.01:\n",
      "D_hf=20.186300205325836\n",
      "error=0.10076328213816765\n",
      "\n",
      "For h=0.001:\n",
      "D_hf=20.095583040074416\n",
      "error=0.01004611688674828\n"
     ]
    }
   ],
   "source": [
    "for x in [0., 1., 3.]:\n",
    "    print()\n",
    "    print(f\"For {x=}:\")\n",
    "    for h in [0.1, 0.01, 0.001]:\n",
    "        print()\n",
    "        print(f\"For {h=}:\")\n",
    "        D_hf = D_h(f, x, h)\n",
    "        print(f\"{D_hf=}\")\n",
    "        Df_exact = f(x)\n",
    "        error = D_hf - Df_exact\n",
    "        print(f\"{error=}\")"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "absolute-shaft",
   "metadata": {},
   "source": [
    "## Mathematical Exercise\n",
    "\n",
    "Verify the order of accuracy result\n",
    "\n",
    "$$\n",
    "Df(x) - \\delta_hf(x) = O(h^2)\n",
    "$$\n",
    "\n",
    "and the (equivalent) fact that this method has degree of precision 2: it is exact fo all quadratics.\n",
    "\n",
    "One could in fact go further and derive the error formula\n",
    "\n",
    "$$\n",
    "Df(x) - \\delta_hf(x) = -\\frac{D^3f(\\xi)}{6} h^2,\\quad\\text{for some } \\xi \\in (x-h, x+h)\n",
    "$$"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "starting-plate",
   "metadata": {},
   "source": [
    "# Improving on the Centered Difference Approximation with Richardson Extrapolation\n",
    "\n",
    "**Slightly revised on Thursday March 25**"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "engaging-jewel",
   "metadata": {},
   "source": [
    "Write a function to apply Richardson extrapolation using the above centered differences function, and a file to test it and handle any interactive input and screen or file output, as noted above.\n",
    "\n",
    "As with almost any practical implementation of a numerical approximation algorithm, the input should include an **error tolerance**, and the output should include an **error estimate**."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "funky-narrative",
   "metadata": {},
   "source": [
    "## Further details and suggestions\n",
    "\n",
    "It could help to break this into to \"sub-tasks\":\n",
    "\n",
    "A) Write a Python function (maybe `CD_richardson_1(f, x, h)`?) that does a single step of Richardson extrapolation forth centered difference approximation of the first derivative.\n",
    "\n",
    "Then as usual add a code cell with testing of this on some examples such as $f(x) = e^x$, $f'(1)$ as suggested in the section on [Test Cases for Differentiation](test-cases-differentiation.ipynb).\n",
    "\n",
    "B) Write a second function that uses a loop to seek a result that meets the error tolerance, by halving $h$ at each iteration.\n",
    "Use either a while loop — or my preferred strategy of using a `for` loop to impose an iteration limit and within that a `break` to stop if and when sufficient accuracy is achieved.\n",
    "\n",
    "Also, see the updated notes on Richardson Extrapolation, where I have added a discussion of how to get error estimates."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "small-calvin",
   "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.8.5"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}