{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {
    "deletable": true,
    "editable": true
   },
   "source": [
    "<center><h2>Interactive Spline Fit with Vpython</h2></center>\n",
    "<p>\n",
    "    Nine points on a plot are fit with cubic splines. The user can slect the number of splines used in the fitting  by moving a slider."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {
    "collapsed": false,
    "deletable": true,
    "editable": true,
    "scrolled": true
   },
   "outputs": [
    {
     "data": {
      "text/html": [
       "<div id=\"glowscript\" class=\"glowscript\"></div>"
      ],
      "text/plain": [
       "<IPython.core.display.HTML object>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "data": {
      "application/javascript": [
       "if (typeof Jupyter !== \"undefined\") { window.__context = { glowscript_container: $(\"#glowscript\").removeAttr(\"id\")};}else{ element.textContent = ' ';}"
      ],
      "text/plain": [
       "<IPython.core.display.Javascript object>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "data": {
      "application/javascript": [
       "if (typeof Jupyter !== \"undefined\") {require.undef(\"nbextensions/vpython_libraries/glow.min\");}else{element.textContent = ' ';}"
      ],
      "text/plain": [
       "<IPython.core.display.Javascript object>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "data": {
      "application/javascript": [
       "if (typeof Jupyter !== \"undefined\") {require.undef(\"nbextensions/vpython_libraries/glowcomm\");}else{element.textContent = ' ';}"
      ],
      "text/plain": [
       "<IPython.core.display.Javascript object>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "data": {
      "application/javascript": [
       "if (typeof Jupyter !== \"undefined\") {require.undef(\"nbextensions/vpython_libraries/jquery-ui.custom.min\");}else{element.textContent = ' ';}"
      ],
      "text/plain": [
       "<IPython.core.display.Javascript object>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "data": {
      "application/javascript": [
       "if (typeof Jupyter !== \"undefined\") {require([\"nbextensions/vpython_libraries/glow.min\"], function(){console.log(\"GLOW LOADED\");});}else{element.textContent = ' ';}"
      ],
      "text/plain": [
       "<IPython.core.display.Javascript object>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "data": {
      "application/javascript": [
       "if (typeof Jupyter !== \"undefined\") {require([\"nbextensions/vpython_libraries/glowcomm\"], function(){console.log(\"GLOWCOMM LOADED\");});}else{element.textContent = ' ';}"
      ],
      "text/plain": [
       "<IPython.core.display.Javascript object>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "data": {
      "application/javascript": [
       "if (typeof Jupyter !== \"undefined\") {require([\"nbextensions/vpython_libraries/jquery-ui.custom.min\"], function(){console.log(\"JQUERY LOADED\");});}else{element.textContent = ' ';}"
      ],
      "text/plain": [
       "<IPython.core.display.Javascript object>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "\"\"\"From \"COMPUTATIONAL PHYSICS\" & \"COMPUTER PROBLEMS in PHYSICS\" \n",
    "by RH Landau, MJ Paez, and CC Bordeianu (deceased). Copyright R Landau, \n",
    "Oregon State Unv, MJ Paez, Univ Antioquia, C Bordeianu (deceased), \n",
    "Univ Bucharest, 2020. Please respect copyright & acknowledge our work.\"\"\"\n",
    "\n",
    "# SplineFitVP.ipynb:  Interactive Spline fit with Vpython slider\n",
    " \n",
    "from vpython import * \n",
    "import numpy as np \n",
    "\n",
    "x = [0., 0.12, 0.25, 0.37, 0.5, 0.62, 0.75, 0.87, 0.99]   # input\n",
    "y = [10.6, 16.0, 45.0, 83.5, 52.8, 19.9, 10.8, 8.25, 4.7]\n",
    "n = 9\n",
    "npt = 15\n",
    "y2 = np.zeros( (n), float);     u = np.zeros( (n), float)\n",
    "graph1 = graph(width=500, height=500, xmin=0,xmax=1, ymin=10,ymax=100,\n",
    "                  title='Spline Fit', xtitle='x', ytitle='y')\n",
    "funct1 = gdots(color=color.blue,radius=5)\n",
    "funct2 = gdots(color = color.red)\n",
    "graph1.visible = 0\n",
    " \n",
    "def S(s):\n",
    "    global Nfit\n",
    "    Nfit=int(s.value)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true,
    "deletable": true,
    "editable": true
   },
   "outputs": [],
   "source": [
    "def update():   \n",
    "    for i in range(0, n):        # Spread out points\n",
    "        funct1.plot(pos = (x[i], y[i]) )\n",
    "        yp1 = (y[1]-y[0]) / (x[1]-x[0]) - (y[2]-y[1])/ \\\n",
    "               (x[2]-x[1])+(y[2]-y[0])/(x[2]-x[0])\n",
    "    ypn = (y[n-1] - y[n-2])/(x[n-1] - x[n-2]) - (y[n-2]-y[n-3])/(x[n-2]-x[n-3]) + (y[n-1]-y[n-3])/(x[n-1]-x[n-3])\n",
    "    if (yp1 > 0.99e30):  y2[0] = 0.;  u[0] = 0.\n",
    "    else:\n",
    "        y2[0] = - 0.5\n",
    "        u[0] = (3./(x[1] - x[0]) )*( (y[1] - y[0])/(x[1] - x[0]) - yp1)\n",
    "    for i in range(1, n - 1):                # Decomp loop\n",
    "        sig = (x[i] - x[i - 1])/(x[i + 1] - x[i - 1]) \n",
    "        p = sig*y2[i - 1] + 2. \n",
    "        y2[i] = (sig - 1.)/p \n",
    "        u[i] = (y[i+1]-y[i])/(x[i+1]-x[i]) - (y[i]-y[i-1])/(x[i]-x[i-1])\n",
    "        u[i] = (6.*u[i]/(x[i + 1] - x[i - 1]) - sig*u[i - 1])/p\n",
    "    if (ypn > 0.99e30):  qn = un = 0.                  # Test for natural\n",
    "    else:\n",
    "        qn = 0.5;\n",
    "        un = (3/(x[n-1]-x[n-2]))*(ypn - (y[n-1]-y[n-2])/(x[n-1]-x[n-2]))\n",
    "    y2[n - 1] = (un - qn*u[n - 2])/(qn*y2[n - 2] + 1.)\n",
    "    for k in range(n - 2, 1,  - 1):  y2[k] = y2[k]*y2[k + 1] + u[k]\n",
    "    for i in range(1, Nfit + 2):                              # Begin fit\n",
    "        xout = x[0] + (x[n - 1] - x[0])*(i - 1)/(Nfit) \n",
    "        klo = 0;    khi = n - 1                         # Bisection algor\n",
    "        while (khi - klo >1):\n",
    "            k = (khi + klo) >> 1\n",
    "            if (x[k] > xout): khi  = k\n",
    "            else: klo = k\n",
    "        h = x[khi] - x[klo] \n",
    "        if (x[k] > xout):  khi = k\n",
    "        else: klo = k \n",
    "        h = x[khi] - x[klo]\n",
    "        a = (x[khi] - xout)/h \n",
    "        b = (xout - x[klo])/h \n",
    "        yout = a*y[klo] + b*y[khi] + ((a*a*a-a)*y2[klo]+(b*b*b-b)*y2[khi])*h*h/6\n",
    "        funct2.plot(pos = (xout, yout) )\n",
    "\n",
    "s = slider( min = 2, max = 100, value=2, bind=S)\n",
    "while 1:\n",
    "    S(s)\n",
    "    rate(50)  \n",
    "    update()\n",
    "    funct2.visible = 0"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 2",
   "language": "python",
   "name": "python2"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 2
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython2",
   "version": "2.7.13"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}