Assignment 1 - vF.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Assignment 1 - Intro to Modeling\n",
    "## Vinícius Miranda"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 1.a)\n",
    "\n",
    "Since $(a^b)^c = a^{b\\times c}$, we can set up the linear systems of equations for each of the fundamental dimensions as the sum of the multiplied exponents. In other words,\n",
    "\n",
    "$$ [M] = \n",
    "\\begin{cases}\n",
    "\\alpha a_1 &+ \\beta a_2 &+ \\gamma a_3 &= 1 \\\\\n",
    "\\alpha b_1 &+ \\beta b_2 &+ \\gamma b_3 &= 0 \\\\\n",
    "\\alpha c_1 &+ \\beta c_2 &+ \\gamma c_3 &= 0\n",
    "\\end{cases}, [L] = \n",
    "\\begin{cases}\n",
    "\\alpha a_1 &+ \\beta a_2 &+ \\gamma a_3 &= 0 \\\\\n",
    "\\alpha b_1 &+ \\beta b_2 &+ \\gamma b_3 &= 1 \\\\\n",
    "\\alpha c_1 &+ \\beta c_2 &+ \\gamma c_3 &= 0\n",
    "\\end{cases}, [T] = \n",
    "\\begin{cases}\n",
    "\\alpha a_1 &+ \\beta a_2 &+ \\gamma a_3 &= 0 \\\\\n",
    "\\alpha b_1 &+ \\beta b_2 &+ \\gamma b_3 &= 0 \\\\\n",
    "\\alpha c_1 &+ \\beta c_2 &+ \\gamma c_3 &= 1\n",
    "\\end{cases}.\n",
    "$$\n",
    "\n",
    "Once we know the dimensions of $x_1, x_2\\text{ and }x_3$, we can use any methods at our disposal to find the appropriate values of $\\alpha, \\beta\\text{ and }\\gamma$. For example, we could use the substitution method, `solve` in Sage Math, or find the RREF of the augmented matrix.\n",
    "\n",
    "The coefficient matrices are all the same. Only the RHS vectors change. The rows in the coefficient matrix represent the exponent of a fundamental dimension. In our case, the first row represents mass, the second length, and the third time. The columns represent our physical quantities $x_i$. The RHS vector describes the final dimensionality we seek to compose. If we try to model a single fundamental unit with the exponent of 1, the RHS will be the zero vector with a single entry changed to one.\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 1.b)\n",
    "The dimensionality of the physical quantities is as follows: \n",
    "\n",
    "- $w \\rightarrow [M][T]^{-2}[L]^{-2}$\n",
    "- $h \\rightarrow [L]$\n",
    "- $q \\rightarrow [L]^3[T]^{-1}$ \n",
    "\n",
    "We find combinations by plugging the exponents in the linear systems outlined above and solving for $\\alpha, \\beta\\text{ and }\\gamma.$ Thus,\n",
    "\n",
    "$$\n",
    "\\begin{align}\n",
    "[M] &= \\begin{cases}\n",
    "1\\alpha  &+ 0\\beta &+ 0\\gamma &= 1\\\\\n",
    "-2\\alpha  &+ 1\\beta &+ 3\\gamma &= 0\\\\\n",
    "-2\\alpha  &+ 0\\beta &+ (-1)\\gamma &= 0\n",
    "\\end{cases} \\Rightarrow \\begin{cases}\n",
    "\\alpha  & & &= 1 \\\\\n",
    "-2\\alpha  &+ \\beta &+ 3\\gamma &= 0  \\\\\n",
    "-2\\alpha  & &-\\gamma &= 0  \n",
    "\\end{cases} \\Rightarrow (\\alpha, \\beta, \\gamma) = (1,8,-2) \\Rightarrow [M] = wh^8q^{-2} \\\\ \\\\\n",
    "[L] &= \\begin{cases}\n",
    "1\\alpha  &+ 0\\beta &+ 0\\gamma &= 0\\\\\n",
    "-2\\alpha  &+ 1\\beta &+ 3\\gamma &= 1\\\\\n",
    "-2\\alpha  &+ 0\\beta &+ (-1)\\gamma &= 0\n",
    "\\end{cases} \\Rightarrow \\begin{cases}\n",
    "\\alpha  & & &= 0 \\\\\n",
    "-2\\alpha  &+ \\beta &+ 3\\gamma &= 1  \\\\\n",
    "-2\\alpha  & &-\\gamma &= 0  \n",
    "\\end{cases} \\Rightarrow (\\alpha, \\beta, \\gamma) = (0,1,0) \\Rightarrow [L] = h \\\\ \\\\\n",
    "[T] &= \\begin{cases}\n",
    "1\\alpha  &+ 0\\beta &+ 0\\gamma &= 0\\\\\n",
    "-2\\alpha  &+ 1\\beta &+ 3\\gamma &= 0\\\\\n",
    "-2\\alpha  &+ 0\\beta &+ (-1)\\gamma &= 1\n",
    "\\end{cases} \\Rightarrow \\begin{cases}\n",
    "\\alpha  & & &= 0 \\\\\n",
    "-2\\alpha  &+ \\beta &+ 3\\gamma &= 0  \\\\\n",
    "-2\\alpha  & &-\\gamma &= 1  \n",
    "\\end{cases} \\Rightarrow (\\alpha, \\beta, \\gamma) = (0,3,-1) \\Rightarrow [T] = h^3q^{-1}\n",
    "\\end{align}\n",
    "$$\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 1.c)\n",
    "\n",
    "We can condense our requirements into a single linear system of equations outlined below.\n",
    "\n",
    "$$\n",
    "\\begin{cases}\n",
    "\\alpha a_1 &+ \\beta a_2 &+ \\gamma a_3 &= a_4 \\\\\n",
    "\\alpha b_1 &+ \\beta b_2 &+ \\gamma b_3 &= b_4 \\\\\n",
    "\\alpha c_1 &+ \\beta c_2 &+ \\gamma c_3 &= c_4\n",
    "\\end{cases}\n",
    "$$\n",
    "\n",
    "We should be able to use the same tools as described in 1.a). We should also note that once we have a solution for $\\alpha, \\beta\\text{ and }\\gamma$ for each of the fundamental units we require, the solution for $[x_4]$ will be simply a linear combination of these fundamental units. That is the case since what we effectively do in solving the linear systems outlined above is to find the standard basis (i.e., a basis of unit vectors) of a three-dimensionial vector space which contains all possible combinations of the three fundamental units as linear combinations of our physical quantities.\n",
    "\n",
    "Our method would fail if our physical quantities are linearly dependent and hence unable to form the basis of the n-dimensional vector space we require. For example, if instead of $w$ we had $A$ (area), our method would fail since $A$ is a linear combination of $h$. We would only be able to form a two-dimensional vector space since we would only have two linearly independent base vectors."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 1.d)\n",
    "\n",
    "Power has the units of $[M][L]^2[T]^{-3}$. Plugging the appropriate values in the system above, we find\n",
    "\n",
    "$$\n",
    "[P] = \\begin{cases}\n",
    "1\\alpha  &+ 0\\beta &+ 0\\gamma &= 1\\\\\n",
    "-2\\alpha  &+ 1\\beta &+ 3\\gamma &= 2\\\\\n",
    "-2\\alpha  &+ 0\\beta &+ (-1)\\gamma &= -3\n",
    "\\end{cases} \\Rightarrow \\begin{cases}\n",
    "\\alpha  & & &= 1 \\\\\n",
    "-2\\alpha  &+ \\beta &+ 3\\gamma &= 2  \\\\\n",
    "-2\\alpha  & &-\\gamma &= -3  \n",
    "\\end{cases} \\Rightarrow (\\alpha, \\beta, \\gamma) = (1,1,1) \\Rightarrow [P] = whq\n",
    "$$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 2.a)\n",
    "\n",
    "Applying Lagrange's interpolation, we find\n",
    "\n",
    "$$\n",
    "\\begin{align}\n",
    "f^{''}_i(x) &= k_i + \\frac{k_{i+1} - k_i}{x_{i+1} - x_i}(x-x_i) \\\\\n",
    "&= k_i\\frac{x-x_{i+1}}{x_i-x_{i+1}} + k_{i+1}\\frac{x-x_i}{x_{i+1}-x_i}\\\\\n",
    "&= k_i\\frac{x-x_{i+1}}{x_i-x_{i+1}} + k_{i+1}\\frac{x-x_i}{x_{i+1}-x_i}\\\\\n",
    "&= k_i\\frac{x_{i+1}-x}{h} + k_{i+1}\\frac{x-x_i}{h}\n",
    "\\end{align}\n",
    "$$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 2.b)\n",
    "\n",
    "We integrate once to find\n",
    "\n",
    "$$\n",
    "\\begin{align}\n",
    "f'_i(x) &= \\int f''_i(x) dx \\\\\n",
    "&=\\int \\left(k_i\\frac{x_{i+1}-x}{h} + k_{i+1}\\frac{x-x_i}{h}\\right)dx \\\\\n",
    "&=\\frac{k_i}{h}\\int (x_{i+1}-x)dx +\\frac{k_{i+1}}{h}\\int(x-x_i)dx \\\\\n",
    "&=-\\frac{k_i}{2h}(x_{i+1}-x)^2 +\\frac{k_{i+1}}{2h}(x-x_i)^2 + A_i\n",
    "\\end{align}$$\n",
    "\n",
    "And once more,\n",
    "\n",
    "$$\n",
    "\\begin{align}\n",
    "f_i(x) &= \\int f'_i(x) dx \\\\\n",
    "&=\\int \\left(-\\frac{k_i}{2h}(x_{i+1}-x)^2 +\\frac{k_{i+1}}{2h}(x-x_i)^2 + A_i\\right) dx \\\\\n",
    "&= -\\frac{k_i}{2h} \\int (x_{i+1}-x)^2 dx +\\frac{k_{i+1}}{2h} \\int (x-x_i)^2 dx + \\int A_i dx \\\\\n",
    "&= \\frac{k_i}{6h} (x_{i+1}-x)^3 +\\frac{k_{i+1}}{6h} (x-x_i)^3 + A_ix + B_i\n",
    "\\end{align}$$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 2.c)\n",
    "\n",
    "Since\n",
    "\n",
    "$$\\begin{align}\n",
    "f_i(x_i) &= \\frac{k_i}{6h} (x_{i+1}-x_i)^3 +\\frac{k_{i+1}}{6h} (x_i-x_i)^3 + A_ix_i + B_i \\\\\n",
    "       &= \\frac{h^2k_i}{6} + A_ix_i + B_i\n",
    "\\\\\n",
    "f_i(x_{i+1}) &= \\frac{k_i}{6h} (x_{i+1}-x_{i+1})^3 +\\frac{k_{i+1}}{6h} (x_{i+1}-x_i)^3 + A_ix_{i+1} + B_i \\\\\n",
    "       &= \\frac{h^2k_{i+1}}{6} + A_ix_{i+1} + B_i\n",
    "\\end{align}$$\n",
    "\n",
    "We compose the system of equations\n",
    "\n",
    "$$\\begin{cases}\n",
    "\\frac{h^2k_i}{6} &+ A_ix_i &+ B_i &= y_i\\\\ \n",
    "\\frac{h^2k_{i+1}}{6} &+ A_ix_{i+1} &+ B_i &= y_{i+1}\n",
    "\\end{cases}$$ \n",
    "\n",
    "Taking $(2) - (1)$, we find\n",
    "\n",
    "$$\n",
    "\\frac{h^2k_{i+1}}{6} - \\frac{h^2k_i}{6} + A_ix_{i+1} - A_ix_i = y_{i+1} - y_i\n",
    "$$\n",
    "\n",
    "$$\n",
    "\\frac{h^2(k_{i+1} - k_i)}{6} + A_i(x_{i+1} - x_i) = y_{i+1} - y_i\n",
    "$$\n",
    "\n",
    "$$\n",
    "\\frac{h^2(k_{i+1} - k_i)}{6} + hA_i = y_{i+1} - y_i\n",
    "$$\n",
    "\n",
    "$$\n",
    "\\frac{h(k_{i+1} - k_i)}{6} + A_i = \\frac{y_{i+1} - y_i}{h}\n",
    "$$\n",
    "\n",
    "$$\n",
    "A_i = \\frac{y_{i+1} - y_i}{h} - \\frac{h(k_{i+1} - k_i)}{6}\n",
    "$$\n",
    "\n",
    "Plugging the result back into (1), we find\n",
    "\n",
    "$$\\frac{h^2k_i}{6} + \\left(\\frac{y_{i+1} - y_i}{h} - \\frac{h(k_{i+1} - k_i)}{6}\\right)x_i + B_i = y_i$$\n",
    "\n",
    "$$B_i = y_i -x_i \\left(\\frac{y_{i+1} - y_i}{h} - \\frac{h(k_{i+1} - k_i)}{6}\\right) - \\frac{h^2k_i}{6}$$\n",
    "\n",
    "$$B_i = y_i -x_i \\left(\\frac{y_{i+1} - y_i}{h} - \\frac{hk_{i+1}}{6} + \\frac{hk_i}{6} \\right) - \\frac{(x_{i+1}-x_i)hk_i}{6}$$\n",
    "\n",
    "$$B_i = y_i -x_i \\left(\\frac{y_{i+1} - y_i}{h} - \\frac{hk_{i+1}}{6}\\right) - \\frac{hx_ik_i}{6} + \\frac{hx_ik_i}{6} - \\frac{x_{i+1}hk_i}{6}$$\n",
    "\n",
    "$$B_i = y_i -x_i \\left(\\frac{y_{i+1} - y_i}{h} - \\frac{hk_{i+1}}{6}\\right) - \\frac{x_{i+1}hk_i}{6}$$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 2.d)\n",
    "\n",
    "We have\n",
    "\n",
    "$$f_{i+1}(x) = \\frac{k_{i+1}}{6h} (x_{i+2}-x)^3 +\\frac{k_{i+2}}{6h} (x-x_{i+1})^3 + A_{i+1}x + B_{i+1}$$\n",
    "\n",
    "\n",
    "where\n",
    "\n",
    "$$\n",
    "A_{i+1} = \\frac{y_{i+2} - y_{i+1}}{h} - \\frac{h(k_{i+2} - k_{i+1})}{6}\n",
    "$$\n",
    "\n",
    "$$B_{i+1} = y_{i+1} -x_{i+1} \\left(\\frac{y_{i+2} - y_{i+1}}{h} - \\frac{hk_{i+2}}{6}\\right) - \\frac{x_{i+2}hk_{i+1}}{6}$$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 2.e)\n",
    "\n",
    "Since\n",
    "\n",
    "$$\n",
    "\\begin{align}\n",
    "f'_i(x_{i+1}) &= -\\frac{k_i}{2h}(x_{i+1}-x_{i+1})^2 +\\frac{k_{i+1}}{2h}(x_{i+1}-x_i)^2 + A_i \\\\\n",
    "&= \\frac{hk_{i+1}}{2} + \\frac{y_{i+1} - y_i}{h} - \\frac{h(k_{i+1} - k_i)}{6} \\\\\n",
    "&= \\frac{y_{i+1} - y_i}{h} + \\frac{h(2k_{i+1} + k_i)}{6}\n",
    "\\end{align}$$\n",
    "\n",
    "and\n",
    "\n",
    "$$\n",
    "\\begin{align}\n",
    "f'_{i+1}(x_{i+1}) &= -\\frac{k_{i+1}}{2h}(x_{i+2}-x_{i+1})^2 +\\frac{k_{i+2}}{2h}(x_{i+1}-x_{i+1})^2 + A_{i+1} \\\\\n",
    "&= -\\frac{hk_{i+1}}{2} + \\frac{y_{i+2} - y_{i+1}}{h} - \\frac{h(k_{i+2} - k_{i+1})}{6}\\\\\n",
    "&= \\frac{y_{i+2} - y_{i+1}}{h} - \\frac{h(k_{i+2} +2 k_{i+1})}{6}\n",
    "\\end{align}$$\n",
    "\n",
    "and\n",
    "\n",
    "$$\n",
    "f'_i(x_{i+1}) = f'_{i+1}(x_{i+1})\n",
    ",$$\n",
    "\n",
    "We find\n",
    "\n",
    "$$\\frac{h(k_i+2k_{i+1})}{6} + \\frac{y_{i+1} - y_i}{h} = \\frac{y_{i+2} - y_{i+1}}{h} - \\frac{h(k_{i+2} + 2k_{i+1})}{6}$$\n",
    "\n",
    "$$\\frac{y_{i+2} - 2y_{i+1} + y_i}{h} = \\frac{h(k_i+4k_{i+1}+ k_{i+2})}{6}$$\n",
    "\n",
    "$$6(y_{i+2} - 2y_{i+1} + y_i) = h^2(k_i+4k_{i+1}+ k_{i+2})$$\n",
    "\n",
    "$$ 12\\left( \\frac{y_{i+2} - y_{i+1}}{(x_{i+2}-x_i)(x_{i+1}-x_i)} - \\frac{y_{i+1} - y_i}{(x_{i+2}-x_i)(x_{i+1}-x_i)}\\right) = k_i+4k_{i+1}+ k_{i+2}\\text{, for }i = 2,\\dots,n-2.$$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 2.f)\n",
    "\n",
    "Using a natural spline, we require $f''_1(x_1)=k_1=0$ and $f''_{n-1}(x_n)=k_n=0$. Furthermore, let $d_1 = d_n = 0$ and  \n",
    "\n",
    "$$d_i = 12\\left( \\frac{y_{i+2} - y_{i+1}}{(x_{i+2}-x_i)(x_{i+1}-x_i)} - \\frac{y_{i+1} - y_i}{(x_{i+2}-x_i)(x_{i+1}-x_i)}\\right)\\text{, } i = 2,3,\\dots, n-2$$\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "$$\n",
    "\\begin{bmatrix}\n",
    "1 & 0 & 0 &&&& \\\\\n",
    "&1 & 4  & 1 &&& \\\\\n",
    "&& \\ddots & \\ddots & \\ddots && \\\\\n",
    "&&& \\ddots & \\ddots & \\ddots & \\\\\n",
    "&&&& \\ddots & \\ddots & \\ddots  \\\\\n",
    "&&&&& 1 & 4  \\\\\n",
    "&&&&&& 1  \n",
    "\\end{bmatrix}\n",
    "\\begin{bmatrix}\n",
    "k_1\\\\\n",
    "k_2\\\\\n",
    "\\vdots\\\\\n",
    "\\vdots\\\\\n",
    "\\vdots\\\\\n",
    "k_{n-1}\\\\\n",
    "k_n\n",
    "\\end{bmatrix}=\n",
    "\\begin{bmatrix}\n",
    "d_1\\\\\n",
    "d_2\\\\\n",
    "\\vdots\\\\\n",
    "\\vdots\\\\\n",
    "\\vdots\\\\\n",
    "d_{n-1}\\\\\n",
    "d_n\n",
    "\\end{bmatrix}\n",
    "$$\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "__Reflection:__\n",
    "\n",
    "I will finish my attempt at this question here. I feel that I have made a lot of progress, but the difference in notation between textbook, question, and Wikiversity (n.d.) is confusing. I feel the question pushed us to consider $f_{i+1}(x)$ instead of $f_{i-1}(x)$ which seem to deviate from the convention and prevent me from finding a way of making the spline work at $n-1$. I'm not sure how to fix it and since this is the third time I restart this question from stratch, I feel I have tried hard enough. Thank you for a good challenge! "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 3.a)\n",
    "\n",
    "For the equations in (1) to be dimensionally homogeneous, we require \n",
    "\n",
    "$$[mg\\sin(\\theta)]=[\\frac{1}{2}\\rho v^2C_DA].$$\n",
    "\n",
    "$\\sin(\\theta)$ must be dimensionless. Furthermore, we know the dimensions of the other quantities to be:\n",
    "\n",
    "- $m: [M]$\n",
    "- $g: [L][T]^{-2}$\n",
    "- $\\rho: [M][L]^{-3}$\n",
    "- $v: [L][T]^{-1}$\n",
    "- $A: [L]^2$\n",
    "\n",
    "Therefore,\n",
    "$$[M][L][T]^{-2} = [M][L]^{-3}([L][T]^{-1})^2[C_D][L]^2$$\n",
    "\n",
    "$$[M][L][T]^{-2} = [M][L]^{-3}[L]^2[T]^{-2}[C_D][L]^2$$\n",
    "\n",
    "$$[M][L][T]^{-2} = [M][L][T]^{-2}[C_D]$$\n",
    "\n",
    "$$[C_D] = \\frac{[M][L][T]^{-2}}{[M][L][T]^{-2}} = []$$\n",
    "\n",
    "We could apply the same process to $C_D$ to find that both  coefficients are dimensionless. "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 3.b) \n",
    "\n",
    "We find:\n",
    "\n",
    "- $[M] = m$\n",
    "- $[L] = \\sqrt A$\n",
    "- $[T] = (g^{-2}A)^{0.25}$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 3.c)\n",
    "\n",
    "$$\\bar{\\rho} = \\rho\\times([M][L]^-3)^{-1} = \\rho\\times [M]^{-1}[L]^3 = \\rho m^{-1}A^{1.5}$$\n",
    "\n",
    "$$\\bar{v} = v\\times ([L][T]^{-1})^{-1}= v\\times [L]^{-1}[T]=vA^{-0.5}(g^{-2}A)^{0.25}=vg^{-0.5}A^{-0.25}$$\n",
    "\n",
    "Therefore,\n",
    "\n",
    "$$v = \\bar{v}g^{0.5}A^{0.25}\\text{ and } \\rho=\\bar{\\rho}mA^{-1.5}$$\n",
    "\n",
    "Plugging those in (1), we find\n",
    "\n",
    "$$mg\\sin(\\theta) + \\frac{1}{2}\\rho v^2A = 0$$\n",
    "\n",
    "$$mg\\sin(\\theta) + \\frac{1}{2}\\bar{\\rho}mA^{-1.5}(\\bar{v}g^{0.5}A^{0.25})^2A = 0$$\n",
    "\n",
    "$$mg\\sin(\\theta) + \\frac{1}{2}\\bar{\\rho}mA^{-1.5}\\bar{v}^2gA^{0.5}A = 0$$\n",
    "\n",
    "$$mg\\left( \\sin(\\theta) +\\frac{1}{2}\\bar{\\rho}\\bar{v}^2\\right) = 0$$\n",
    "\n",
    "$$\\sin(\\theta) +\\frac{\\bar{\\rho}\\bar{v}^2}{2} = 0$$\n",
    "\n",
    "We follow the same process for the other equation to find \n",
    "\n",
    "$$\\cos(\\theta) -\\frac{\\bar{\\rho}\\bar{v}^2}{2} = 0$$\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 3.d)\n",
    "\n",
    "Assuming $\\bar{\\rho} = 1$, we have the system of equations\n",
    "\n",
    "$$\n",
    "\\begin{cases}\n",
    "2\\sin\\theta+\\bar{v}^2=0 \\\\\n",
    "2\\cos\\theta-\\bar{v}^2=0\n",
    "\\end{cases}\n",
    "$$\n",
    "\n",
    "We also note\n",
    "\n",
    "$$\n",
    "\\begin{align}\n",
    "f(\\theta, \\bar{v}) &= 2\\sin\\theta+\\bar{v}^2. &f_\\theta(\\theta, \\bar{v}) &= 2\\cos\\theta. &f_{\\bar{v}}(\\theta, \\bar{v}) &= 2\\bar{v}. \\\\\n",
    "g(\\theta, \\bar{v}) &= 2\\cos\\theta-\\bar{v}^2. &g_\\theta(\\theta, \\bar{v}) &= -2\\sin\\theta. &g_{\\bar{v}}(\\theta, \\bar{v}) &= -2\\bar{v}. \\\\\n",
    "\\end{align}\n",
    "$$\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 95,
   "metadata": {},
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "from numpy.linalg import det, inv\n",
    "\n",
    "def solver_stats(fun, *args, **kwargs):\n",
    "    def wrapper(*args, **kwargs):\n",
    "        a, b, c = fun(*args, **kwargs)\n",
    "        x = a[-1].flatten()\n",
    "        y = b[-1].flatten()\n",
    "        print(\"Solution found at {} in {:d} steps. f(x) = {:.2f}, g(x) = {:.2f}.\".format(x, c, y[0], y[1]))\n",
    "        return(a, b, c)\n",
    "    return wrapper\n",
    "\n",
    "def fns(x):\n",
    "    theta, vbar = x\n",
    "    return(np.array([2*np.sin(theta)+vbar**2, 2*np.cos(theta)-vbar**2]).reshape((2,)))\n",
    "\n",
    "def jacobian(x):\n",
    "    theta, vbar = x\n",
    "    return(np.array([[2*np.cos(theta), 2*vbar], [-2*np.sin(theta), -2*vbar]]))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 96,
   "metadata": {},
   "outputs": [],
   "source": [
    "@solver_stats\n",
    "def newton_solver(fns=fns, jacobian=jacobian, x0=(1e-5,1e-5), maxiter=1000, tol = 0.001):\n",
    "    '''Newton's method Solver.\n",
    "       Inputs:\n",
    "       \n",
    "           - fns: System of equations to minimize\n",
    "           - jacobian: the Jacobian of the system of equations \n",
    "           - x0: Initial input\n",
    "           - maxiter: maximum number of iterations\n",
    "           - tol: Tolerance value\n",
    "           \n",
    "       Outputs: \n",
    "       \n",
    "           - x_history: inputs for all iterations\n",
    "           - y_history: y value for all iterations\n",
    "           - k: number of steps\n",
    "    '''\n",
    "    \n",
    "    dim = len(x0)\n",
    "    k = 0  # iterations\n",
    "    x = np.array(x0, dtype=float).reshape((dim, ))\n",
    "    y = fns(x)\n",
    "    x_history = [np.copy(x)]\n",
    "    f_history = [y]\n",
    "\n",
    "    while not np.all(np.abs(y) < tol*np.ones_like(x)) and k < maxiter:\n",
    "        \n",
    "        # Compute the jacobian\n",
    "        jac = jacobian(x)\n",
    "        is_singular = not bool(det(jac))\n",
    "        if is_singular: raise ValueError(\"Jacobian is singular.\")\n",
    "        \n",
    "        # Update\n",
    "        y = fns(x)\n",
    "        x -= np.matmul(inv(jac), y)\n",
    "        k += 1\n",
    "        \n",
    "        # Store\n",
    "        x_history.append(np.copy(x))\n",
    "        f_history.append(y)\n",
    "\n",
    "    return(x_history, f_history, k)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 97,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Solution found at [-0.78539816  1.18920712] in 21 steps. f(x) = 0.00, g(x) = -0.00.\n"
     ]
    }
   ],
   "source": [
    "xhis, yhis, _ = newton_solver()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### HCs\n",
    "\n",
    "\\#modeling: The third problem required the mathematical modeling of a real-world problem through the non-dimensionalization of key variables. Furthermore, the first problem also required a great deal of abstraction to model non-dimensionalization more generally.\n",
    "\n",
    "\\#organization: Hopefully, the document is clearly outlined and the solutions are conveyed appropriately.\n",
    "\n",
    "\\#effort: <- should be an HC.\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### References\n",
    "\n",
    "Wikiversity. (n.d.). Cubic Spline Interpolation. Retrieved from https://en.wikiversity.org/w/index.php?title=Cubic_Spline_Interpolation&oldid=2011779"
   ]
  }
 ],
 "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.7.4"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}