{
 "cells": [
  {
   "cell_type": "markdown",
   "id": "b7c92b5b-0a0f-40fa-8039-bdf8898e8691",
   "metadata": {},
   "source": [
    "This notebook is based on section (c) on **computing with differential numbers**\n",
    "running on pages 7 to 10 of *Introduction to Numerical Analysis*\n",
    "by Arnold Neuimaier, Cambridge 2001.\n",
    "The operator overloading encodes the elementary derivative rules\n",
    "for the addition, subtraction, multiplication, and division of\n",
    "any two real numbers.  As we evaluate any expression written in\n",
    "these four elementary operations, we also compute its derivative.\n",
    "The notebook illustrates the forward mode of algorithmic differentiation."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "b46f3ed4-93e3-4773-b3c1-16c44f53905a",
   "metadata": {},
   "source": [
    "# 1. A Class Definition"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "id": "4247427d-d19b-4dd3-9a14-08ba9a849955",
   "metadata": {},
   "outputs": [],
   "source": [
    "class DifferentialNumber(object):\n",
    "    \"\"\"\n",
    "    A differential number is a tuple df = (f, f') of two doubles\n",
    "    where f' is the evaluated derivative of f.\n",
    "    \"\"\"\n",
    "    def __init__(self, xv, xp=1):\n",
    "        \"\"\"\n",
    "        Initializes a differential number\n",
    "        to the tuple (xv, xp).\n",
    "        \"\"\"\n",
    "        self.val = float(xv)\n",
    "        self.der = float(xp)\n",
    "\n",
    "    def __str__(self):\n",
    "        \"\"\"\n",
    "        Returns the tuple representation of the differential number.\n",
    "        \"\"\"\n",
    "        return str((self.val, self.der))\n",
    "\n",
    "    def __div__(self, other):\n",
    "        return self\n",
    "\n",
    "    def __add__(self, other):\n",
    "        \"\"\"\n",
    "        Defines the addition of two differential numbers.\n",
    "        \"\"\"\n",
    "        if isinstance(other, float):\n",
    "            return DifferentialNumber(self.val + other, self.der)\n",
    "        else:\n",
    "            return DifferentialNumber(self.val + other.val, \\\n",
    "                                      self.der + other.der)\n",
    "\n",
    "    def __radd__(self, other):\n",
    "        \"\"\"\n",
    "        Addition when operand is not a differential number\n",
    "        as in 2.7 + x or 3 + x (reflected operand).\n",
    "        \"\"\"\n",
    "        result = self + other\n",
    "        return result\n",
    "\n",
    "    def __sub__(self, other):\n",
    "        \"\"\"\n",
    "        Defines the subtraction of two differential numbers.\n",
    "        \"\"\"\n",
    "        if isinstance(other, float):\n",
    "            return DifferentialNumber(self.val - other, self.der)\n",
    "        else:\n",
    "            return DifferentialNumber(self.val - other.val, \\\n",
    "                                      self.der - other.der)\n",
    "\n",
    "    def __mul__(self, other):\n",
    "        \"\"\"\n",
    "        Defines the product of two differential numbers.\n",
    "        \"\"\"\n",
    "        if isinstance(other, float):\n",
    "            return DifferentialNumber(self.val*other, self.der*other)\n",
    "        else:\n",
    "            return DifferentialNumber(self.val*other.val, \\\n",
    "                                      self.der*other.val + self.val*other.der)\n",
    "\n",
    "    def __truediv__(self, other):\n",
    "        \"\"\"\n",
    "        Defines the division of two differential numbers.\n",
    "        \"\"\"\n",
    "        if isinstance(other, float):\n",
    "            return DifferentialNumber(self.val/other, self.der/other)\n",
    "        else:\n",
    "            val = self.val/other.val\n",
    "            return DifferentialNumber(val, \\\n",
    "                                     (self.der - val*other.der)/other.val)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "df99fad9-8e7b-4f60-b467-e23d36314ebf",
   "metadata": {},
   "source": [
    "# 2. Using the Class"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "4727dc55-7113-4313-af78-6cfa9f162838",
   "metadata": {},
   "source": [
    "Consider the expression\n",
    "$$\n",
    "    \\frac{(x-1) (x+3)}{x+2}.\n",
    "$$\n",
    "We want to evaluate and differentiate this at 3."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "id": "3372cb34-d478-4e56-b7f4-1cd92dcc64cc",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "'(3.0, 1.0)'"
      ]
     },
     "execution_count": 6,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "dfx = DifferentialNumber(3, 1)\n",
    "str(dfx)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "88d088b0-ee07-49d4-98c7-9909d4fb7b02",
   "metadata": {},
   "source": [
    "Let us test the reflected addition."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "id": "698487b3-5950-45b6-9d53-c5acc6eb576b",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "x + 1.0 at  (3.0, 1.0) : (4.0, 1.0)\n",
      "1.0 + x at  (3.0, 1.0) : (4.0, 1.0)\n"
     ]
    }
   ],
   "source": [
    "f1a = lambda x: x + 1.0\n",
    "f1b = lambda x: 1.0 + x\n",
    "print('x + 1.0 at ', dfx, ':', f1a(dfx))\n",
    "print('1.0 + x at ', dfx, ':', f1b(dfx))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "id": "59d81d7e-0797-4cda-ab41-452dd06f8878",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "(x-1)*(x+3) at (3.0, 1.0) : (12.0, 8.0)\n"
     ]
    }
   ],
   "source": [
    "fun1 = lambda x: (x-1.0)*(x+3.0)\n",
    "print('(x-1)*(x+3) at', dfx, ':', fun1(dfx))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "id": "105286d3-4c2a-436e-939e-878ea20e462a",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "((x-1)*(x+3))/(x+2) at (3.0, 1.0) : (2.4, 1.1199999999999999)\n"
     ]
    }
   ],
   "source": [
    "fun2 = lambda x: ((x-1.0)*(x+3.0))/(x+2.0)\n",
    "print('((x-1)*(x+3))/(x+2) at', dfx, ':', fun2(dfx))"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "1b80ad13-d6f2-4d62-a6af-72e4a6fd888c",
   "metadata": {},
   "source": [
    "Let us verify with sympy."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "id": "408ed3ea-cf3c-4720-90cf-903ef4d44020",
   "metadata": {},
   "outputs": [],
   "source": [
    "from sympy import var, diff"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "id": "2d15e893-f709-4378-8d0f-0411a7243eaf",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/latex": [
       "$\\displaystyle \\frac{\\left(x - 1.0\\right) \\left(x + 3.0\\right)}{x + 2.0}$"
      ],
      "text/plain": [
       "(x - 1.0)*(x + 3.0)/(x + 2.0)"
      ]
     },
     "execution_count": 10,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "x = var('x')\n",
    "e = fun2(x)\n",
    "e"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "id": "e2111eab-41b8-4810-8d60-6c457e83f0c8",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/latex": [
       "$\\displaystyle 2.4$"
      ],
      "text/plain": [
       "2.40000000000000"
      ]
     },
     "execution_count": 11,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "e.subs({x: 3})"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "id": "f23328d5-b5a0-49da-8605-e870d20eda61",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/latex": [
       "$\\displaystyle 1.12$"
      ],
      "text/plain": [
       "1.12000000000000"
      ]
     },
     "execution_count": 12,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "d = diff(e, x)\n",
    "d.subs({x: 3})"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "23317cfc-b321-44b9-b3ec-2e69376f6366",
   "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.4"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}