Problem 4

Application Problem 4.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {
    "hideCode": true,
    "hidePrompt": true
   },
   "source": [
    "Civil Engineering Problem\n",
    "\n",
    "CEE384\n",
    "\n",
    "Adam Badahdah"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "button": false,
    "deletable": true,
    "hideCode": true,
    "hidePrompt": true,
    "new_sheet": false,
    "run_control": {
     "read_only": false
    }
   },
   "source": [
    "The vertical stress under the corner of a rectangular footing subjected to a uniform load of intensity q is given by the solution of Boussinesq’s equation,\n",
    "\n",
    "\n",
    "$$\\sigma = \\frac{q}{4\\pi}\\left[ \n",
    "\\frac{2mn\\sqrt{m^2+n^2+1}m^2+n^2+2}{m^2+n^2+1+m^2n^2+1}\n",
    "+\\sin^{-1}\\left( \\frac{2mn\\sqrt{m^2+n^2+1}}{m^2+n^2+1+m^2n^2}\\right)\n",
    "\\right]    $$\n",
    "\n",
    "\n",
    "Where m and n are dimensionless ratios with $m=a/z$ and $n=b/z$ , where $a=4.6$ and $b = 14$ are the dimensions of the footing and z is the depth below the footing. If a=4.6 and b=14, use a third-order interpolating polynomial to compute the stress, $\\sigma$, at a depth $z = 10m$ below the corner of a footing that subjected to a total load of 100 tones. Note that you will need to use Boussinesq’s equation to generate enough data points over different depths to create an interpolated polynomial."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {
    "collapsed": true,
    "hideCode": true,
    "hideOutput": false,
    "hidePrompt": true
   },
   "outputs": [],
   "source": [
    "#To display latex in output cells\n",
    "from IPython.display import Latex\n",
    "\n",
    "#Import math operators\n",
    "from math import *\n",
    "import numpy as np"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {
    "button": "none",
    "collapsed": false,
    "deletable": true,
    "hideCode": true,
    "hidePrompt": true,
    "new_sheet": false,
    "run_control": {
     "read_only": false
    }
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "11.0053183378\n"
     ]
    }
   ],
   "source": [
    "q = 100\n",
    "a = 4.6;\n",
    "b = 14;\n",
    "\n",
    "def func(a, b, z, q):\n",
    "    \n",
    "    #define m and n\n",
    "    m = a/z;\n",
    "    n = b/z;\n",
    "    \n",
    "    #define some variables\n",
    "    var1 = 100/(4*pi); \n",
    "    var2 = 2*m*n*sqrt(m**2+m**2+1)*m**2+n**2+2;\n",
    "    var3 = m**2+n**2+1+m**2*n**2*m**2+n**2+1;\n",
    "    var4 = np.arcsin((2*m*n*sqrt(m**2+n**2+1))/(m**2+n**2+1+m**2*n**2));\n",
    "    sigma = var1*(var2/var3+var4);\n",
    "    \n",
    "    return sigma\n",
    "\n",
    "\n",
    "sigma = func(a,b,10, q)\n",
    "print(sigma)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {
    "button": false,
    "collapsed": true,
    "deletable": true,
    "hideCode": true,
    "hidePrompt": true,
    "new_sheet": false,
    "run_control": {
     "read_only": false
    }
   },
   "outputs": [],
   "source": [
    "#Import commands for MATLAB-like plotting\n",
    "from numpy import linspace, array\n",
    "from matplotlib.pylab import *"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {
    "collapsed": false,
    "hideCode": true,
    "hideOutput": false,
    "hidePrompt": true,
    "run_control": {}
   },
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAY8AAAEWCAYAAACe8xtsAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAGBhJREFUeJzt3X+U3XWd3/HnyxCX0bKOQGA3ExC20hQqStY5rL+rxTVA\nPZKlqxWx4pFzKKe623pO45Lj6Xr29LSg6R5bXT0sVQ62lfhjDRELGkBw0VNYnRgwQQiyiDoT10Rg\n/HEcNMR3/7h3YDLOJPOdzP0xmefjnDlz53M/33tf+RJ48f0x95OqQpKkJp7R6wCSpMXH8pAkNWZ5\nSJIaszwkSY1ZHpKkxiwPSVJjlofUJ5IMJPlCkp8k+Wx7bEWSB5IMHGLbP0ny/u4kleCoXgeQ9JQ/\nBk4EjquqJ9tjVwDXVdXEIbb9n8BDSf6yqvZ0MqQEHnlI/eR5wIOTxZHkt4BLgP9zqA2r6gngi8Db\nOppQarM8pLYk/zrJz6d8/TLJV2aZNzJt7N1Jbmw/Pj/Jt5P8LMlYkv84h/f+C+DPgckMlwJ/AIxX\n1Wh7zkun5XsiySNTXuYrwL+c5x9fasTTVlJbVX0a+DRAkt8G/g7YNMPULwAfS3JaVX2nPfYW4C/b\njz8OvKmqvprkucCpc3jv9yUp4PlV9dZ2hncCu6bMuQv4R+3nlgO3Af9vysvcD7xojn9c6bB45CFN\nk+QZwPXAV6rqr6c/X1W/AD4PXNSefxrwT4Eb21P2AWck+e2qeryqvjnPKIPAz2Z57kPt5947Zexn\nwHPm+V5SI5aH9Jv+C3AM8KcHmXM97fKgddSxpV0qAP8KOB/4XpK/TfLSeeZ4vJ3jAEn+LfBq4C1V\n9espTx0D/GSe7yU1YnlIUyR5M61S+OOq2neQqbcCK5Kc1Z5//eQTVfWNqroAOAHYAnxmnnG+BfyT\nafleCfxn4IKq+um0+acD987zvaRGLA+pLcka4MPAuqrae7C57WL5LLAROJZWmZDkmUkuTvKc9pyf\nAr+e/ZUO6uvAYJKh9mufRKuI3lZVD84w/5/TuuNK6jjLQ3raBcBzga9NuaPpYP8xvh54LfDZKb+X\nAfBvgEeS/BS4HLgYIMnJ7dc8eS5hqupXwHXAW9tD59D6PZC/mZLvvvZrH03rVNkn5vhnlQ5LXAxK\n6l9JVgBfBdYc7BcFk/wJcFJVvadr4bSkWR6SpMY8bSVJaszykCQ1ZnlIkho7Yj+e5Pjjj69TTjml\n1zEkadHYtm3bj6tqxVzmHrHlccoppzAyMnLoiZIkAJJ8b65zPW0lSWrM8pAkNWZ5SJIaO2KveUhS\nr+3bt4/R0VGeeOKJXkc5wNFHH82qVatYvnz5vF/D8pCkDhkdHeWYY47hlFNOIUmv4wBQVTz66KOM\njo5y6qmHXKdsVl05bZXk2iR7kuycMrYxyQNJvpXkhiSDs2z7SJIdSe6ZvvSnJPWzJ554guOOO65v\nigMgCccdd9xhHw1165rHdcC508ZuBV5QVS8EHgQ2HGT711TVWVU13KF8AGzZPsbLr7qdU6+4iZdf\ndTtbto918u0kLQH9VByTFiJTV8qjqu4EHps2dsuUj7G+G1jVjSyz2bJ9jA2bdzA2PkEBY+MTbNi8\nwwKRpBn0y91W72D2RWwKuC3JtiSXHexFklyWZCTJyN69B13L5zds3LqLiX37Dxib2LefjVt3NXod\nSeonH/rQhzj99NO5+OKLF/R1e37BPMl7gSeBT84y5RVVNZbkBODWJA+0j2R+Q1VdA1wDMDw83Oiz\n5nePz7xUwmzjkrQYfPSjH+W2225j1aqFPbnT0yOPJG8HXg9cXLMsLFJVY+3ve4AbgLM7kWXl4ECj\ncUlaaAt93fXyyy/n4Ycf5rzzzuODH/zgAqVs6Vl5JDkXeA/whqr6xSxznp3kmMnHwOuAnTPNPVzr\n165mYPmyA8YGli9j/drVnXg7STpAJ667Xn311axcuZI77riDd7/73QsXlu7dqrsJuAtYnWQ0yaXA\nXwHH0DoVdU+Sq9tzVya5ub3pibTWk74X+DpwU1V9qRMZ160Z4soLz2RocIAAQ4MDXHnhmaxbM9SJ\nt5OkAyy2665dueZRVRfNMPzxWebuBs5vP34YeFEHox1g3Zohy0JSTyy26679creVJC1pi+26q+Uh\nSX1gsV137fmtupIknjplvnHrLnaPT7BycID1a1cf9qn0Rx55ZAHS/SbLQ5L6xGK67uppK0lSY5aH\nJHXQLL//3FMLkcnykKQOOfroo3n00Uf7qkAm1/M4+uijD+t1vOYhSR2yatUqRkdHafpBrZ02uZLg\n4bA8JKlDli9fflir9fUzT1tJkhqzPCRJjVkekqTGLA9JUmOWhySpMctDktSY5SFJaszykCQ1ZnlI\nkhqzPCRJjVkekqTGLA9JUmOWhySpMctDktSY5SFJaszykCQ1ZnlIkhqzPCRJjVkekqTGXMN8Ediy\nfYyNW3exe3yClYMDrF+7mnVrhnodS9ISZnn0uS3bx9iweQcT+/YDMDY+wYbNOwAsEEk942mrPrdx\n666nimPSxL79bNy6q0eJJMny6Hu7xycajUtSN3SlPJJcm2RPkp1TxjYmeSDJt5LckGRwlm3PTbIr\nyUNJruhG3n6ycnCg0bgkdUO3jjyuA86dNnYr8IKqeiHwILBh+kZJlgEfAc4DzgAuSnJGZ6P2l/Vr\nVzOwfNkBYwPLl7F+7eoeJZKkLpVHVd0JPDZt7JaqerL9493Aqhk2PRt4qKoerqpfAZ8CLuho2D6z\nbs0QV154JkODAwQYGhzgygvP9GK5pJ7ql7ut3gF8eobxIeAHU34eBf5gthdJchlwGcDJJ5+8kPl6\nat2aIctCUl/p+QXzJO8FngQ+ebivVVXXVNVwVQ2vWLHi8MNJkmbU0yOPJG8HXg+cU1U1w5Qx4KQp\nP69qj0mSeqhnRx5JzgXeA7yhqn4xy7RvAKclOTXJM4E3Azd2K6MkaWbdulV3E3AXsDrJaJJLgb8C\njgFuTXJPkqvbc1cmuRmgfUH9XcBW4H7gM1V1XzcyS5Jml5nPFi1+w8PDNTIy0usYkrRoJNlWVcNz\nmdvzC+aSpMXH8pAkNWZ5SJIaszwkSY1ZHpKkxiwPSVJjlockqTHLQ5LUmOUhSWrM8pAkNWZ5SJIa\nszwkSY1ZHpKkxiwPSVJjlockqTHLQ5LUmOUhSWrM8pAkNWZ5SJIaO6rXAbR4bdk+xsatu9g9PsHK\nwQHWr13NujVDvY4lqQssD83Llu1jbNi8g4l9+wEYG59gw+YdABaItAR42krzsnHrrqeKY9LEvv1s\n3LqrR4kkdZPloXnZPT7RaFzSkcXy0LysHBxoNC7pyGJ5aF7Wr13NwPJlB4wNLF/G+rWre5RIUjd5\nwVzzMnlR3LutpKXJ8tC8rVszZFlIS5SnrSRJjVkekqTGLA9JUmOWhySpsa6UR5Jrk+xJsnPK2BuT\n3Jfk10mGD7LtI0l2JLknyUg38kqSDq5bRx7XAedOG9sJXAjcOYftX1NVZ1XVrCUjSeqertyqW1V3\nJjll2tj9AEm6EUGStIAWwzWPAm5Lsi3JZQebmOSyJCNJRvbu3duleJK09CyG8nhFVZ0FnAe8M8mr\nZptYVddU1XBVDa9YsaJ7CSVpien78qiqsfb3PcANwNm9TSRJ6uvySPLsJMdMPgZeR+tCuySph7p1\nq+4m4C5gdZLRJJcm+aMko8BLgZuSbG3PXZnk5vamJwJfS3Iv8HXgpqr6UjcyS5Jm1627rS6a5akb\nZpi7Gzi//fhh4EUdjCZJmoe+Pm0lSepPlockqTHLQ5LUmOUhSWqscXm0b59dduiZkqQj1SHLI8kz\nkrwlyU1J9gAPAD9M8u0kG5M8v/MxJUn9ZC5HHncA/xjYAPxOVZ1UVScArwDuBt6f5K0dzChJ6jNz\n+T2P11bVvumDVfUY8Dngc0mWL3gyaZ62bB9j49Zd7B6fYOXgAOvXrmbdmqFex5KOKIc88pgsjiRf\nTnL+1OeSXDN1jtRrW7aPsWHzDsbGJyhgbHyCDZt3sGX7WK+jSUeUJhfMTwX+LMn7poy5OJP6ysat\nu5jYt/+AsYl9+9m4dVePEklHpiblMQ6cA5yY5AtJntOhTNK87R6faDQuaX6alEeq6smq+ne0rnV8\nDTihM7Gk+Vk5ONBoXNL8NCmPqycfVNV1wNuBWxY4j3RY1q9dzcDyA38NaWD5MtavXd2jRNKR6ZB3\nWyVJtfz11PGq2ga8Y+qcDmWU5mzyrirvtpI6ay636t6R5HPA56vq+5ODSZ5J63c9LqH1uyDXdSSh\n1NC6NUOWhdRhcymPc2kdYWxK8nvA48AArVNetwD/vaq2dy6iJKnfHLI8quoJ4KNJTgD+K3AcMFFV\n450OJ0nqT01WEvxPtI44jgW+mWSTBSJJS1PTT9V9AtgKnATcleSshY8kSep3TY48Hqiqyd8u/5sk\n19G6ffdfLHgqSVJfa3Lk8eMkL578oaoeBFYsfCRJUr9rcuTxp8CnkmwDdgAvBL7bkVSSpL425yOP\nqroXOAvY1B66A7ioE6EkSf2tyZEHVfVL4Kb2lyRpiWq8hrkkSZaHJKkxy0OS1JjlIUlqzPKQJDVm\neUiSGrM8JEmNWR6SpMa6Uh5Jrk2yJ8nOKWNvTHJfkl8nGT7Itucm2ZXkoSRXdCOv1Albto/x8qtu\n59QrbuLlV93Olu1jvY4kzVu3jjyuo7Ui4VQ7gQuBO2fbKMky4CPAecAZwEVJzuhQRqljtmwfY8Pm\nHYyNT1DA2PgEGzbvsEC0aHWlPKrqTuCxaWP3V9WuQ2x6NvBQVT1cVb8CPgVc0KGYUsds3LqLiX37\nDxib2LefjVsP9a+A1J/6/ZrHEPCDKT+PtsdmlOSyJCNJRvbu3dvxcNJc7R6faDQu9bt+L49Gquqa\nqhququEVK1xqRP1j5eBAo3Gp3/V7eYzRWvJ20qr2mLSorF+7moHlyw4YG1i+jPVrV/cokXR4Gn0k\new98Azgtyam0SuPNwFt6G0lqbt2a1tnWjVt3sXt8gpWDA6xfu/qpcWmx6Up5JNkEvBo4Psko8D5a\nF9A/TGsp25uS3FNVa5OsBD5WVedX1ZNJ3gVsBZYB11bVfd3ILC20dWuGLAsdMVJVvc7QEcPDwzUy\nMtLrGJK0aCTZVlWz/t7dVP1+zUOS1IcsD0lSY5aHJKkxy0OS1JjlIUlqzPKQJDVmeUiSGrM8JEmN\nWR6SpMYsD0lSY5aHJKkxy0OS1JjlIUlqrN/X85DUYVu2j7nOiBqzPKQlbMv2MTZs3sHEvv0AjI1P\nsGHzDgALRAflaStpCdu4dddTxTFpYt9+Nm7d1aNEWiwsD2kJ2z0+0WhcmmR5SEvYysGBRuPSJMtD\nWsLWr13NwPJlB4wNLF/G+rWre5RIi4UXzKUlbPKiuHdbqSnLQ1ri1q0ZsizUmKetJEmNWR6SpMYs\nD0lSY5aHJKkxy0OS1JjlIUlqzPKQJDVmeUiSGrM8JEmNWR6SpMYsD0lSY135bKsk1wKvB/ZU1Qva\nY8cCnwZOAR4B3lRVj8+w7SPAz4D9wJNVNdyNzJJ6y+Vx+1u3jjyuA86dNnYF8OWqOg34cvvn2bym\nqs6yOKSlYXJ53LHxCYqnl8fdsn2s19HU1pXyqKo7gcemDV8AfKL9+BPAum5kkdT/XB63//XymseJ\nVfXD9uN/AE6cZV4BtyXZluSyg71gksuSjCQZ2bt370JmldRFLo/b//rignlVFa2SmMkrquos4Dzg\nnUledZDXuaaqhqtqeMWKFZ2IKqkLXB63//WyPH6U5HcB2t/3zDSpqsba3/cANwBndy2hpJ5wedz+\n18vyuBG4pP34EuDz0yckeXaSYyYfA68DdnYtoaSeWLdmiCsvPJOhwQECDA0OcOWFZ3q3VR/p1q26\nm4BXA8cnGQXeB1wFfCbJpcD3gDe1564EPlZV59O6DnJDksms11fVl7qRWVJvuTxuf+tKeVTVRbM8\ndc4Mc3cD57cfPwy8qIPRJEnz0BcXzCVJi4vlIUlqzPKQJDVmeUiSGrM8JEmNWR6SpMYsD0lSY5aH\nJKkxy0OS1FhXfsNcko4UrnDYYnlI0hxNrnA4uVDV5AqHwJIrEE9bSdIcucLh0ywPSZojVzh8muUh\nSXPkCodPszwkaY5c4fBpXjCXpDmavCju3VaWhyQ14gqHLZ62kiQ1ZnlIkhqzPCRJjVkekqTGLA9J\nUmOWhySpMctDktSY5SFJasxfEpSkRa4Xa4xYHpK0iPVqjRFPW0nSItarNUYsD0laxHq1xojlIUmL\nWK/WGLE8JGkR69UaI10pjyTXJtmTZOeUsWOT3JrkO+3vz51l23OT7EryUJIrupFXkhaLdWuGuPLC\nMxkaHCDA0OAAV154ZsfvtkpVdfQNAJK8Cvg58L+q6gXtsQ8Aj1XVVe1SeG5V/dm07ZYBDwJ/CIwC\n3wAuqqpvH+o9h4eHa2RkZIH/JJJ05EqyraqG5zK3K0ceVXUn8Ni04QuAT7QffwJYN8OmZwMPVdXD\nVfUr4FPt7SRJPdTLax4nVtUP24//AThxhjlDwA+m/DzaHpMk9VBfXDCv1rmzwz5/luSyJCNJRvbu\n3bsAySRJM+llefwoye8CtL/vmWHOGHDSlJ9XtcdmVFXXVNVwVQ2vWLFiQcNKkp7Wy/K4Ebik/fgS\n4PMzzPkGcFqSU5M8E3hzeztJUg9161bdTcBdwOoko0kuBa4C/jDJd4DXtn8mycokNwNU1ZPAu4Ct\nwP3AZ6rqvm5kliTNriu36vZCkr3A9+a5+fHAjxcwzkIxVzPmasZczfRjrsPN9LyqmtM5/yO2PA5H\nkpG53uvcTeZqxlzNmKuZfszVzUx9cbeVJGlxsTwkSY1ZHjO7ptcBZmGuZszVjLma6cdcXcvkNQ9J\nUmMeeUiSGrM8JEmNLenySPLuJPcl2ZlkU5Kjpz2fJB9qryXyrSS/3ye5Xp3kJ0nuaX/9eZdy/ft2\npvuS/IcZnu/V/jpUrq7sr35dt+Ywcz2SZEd7vy3oGgez5Hpj+5/jr5PMestpD/bXXHN1ZH/Nkmlj\nkgfa/67dkGRwlm07s6+qakl+0fp03u8CA+2fPwO8fdqc84EvAgFeAvxdn+R6NfB/u7y/XgDsBJ4F\nHAXcBjy/D/bXXHJ1ZX8BrwJ+H9g5ZewDwBXtx1cA759hu2XA3wO/BzwTuBc4o9e52s89Ahzfxf11\nOrAa+AowPMt2vdhfh8zVyf01S6bXAUe1H7+/23+3lvSRB63/2AwkOYrWf3x2T3v+AloLWFVV3Q0M\nTn6YY49z9cLptMrgF9X62Ji/BS6cNqcX+2suubqi+nTdmsPI1VEz5aqq+6tq1yE27fr+mmOujpkl\n0y3tv/MAd9P64NjpOravlmx5VNUY8N+A7wM/BH5SVbdMm9b19UTmmAvgZe3D1S8m+WedzNS2E3hl\nkuOSPIvWUcZJ0+b0Yv2VueSC7u+vSf26bs1cckFrqYTbkmxLclmHM81VP6/z06v99Q5aR/3TdWxf\nLdnyaJ/jvQA4FVgJPDvJW3ubas65vgmcXFUvBD4MbOl0rqq6n9ah8S3Al4B7gP2dft9DmWOuru+v\nmVTrPELf3Rt/iFyvqKqzgPOAd6a1pLRm1/X9leS9wJPAJzv9XlMt2fKg9Um+362qvVW1D9gMvGza\nnEbriXQrV1X9tKp+3n58M7A8yfEdzkVVfbyqXlxVrwIep7W+/FS92F+HzNWr/dW24OvWdDHX5JEw\nVbUHuIHWaZBe68nfs7no9v5K8nbg9cDF7f8JmK5j+2opl8f3gZckeVaSAOfQ+tj3qW4E3ta+i+gl\ntE4h/XD6C3U7V5LfaT9HkrNp/XN8tMO5SHJC+/vJtK4rXD9tSi/21yFz9Wp/tfXrujWHzJXk2UmO\nmXxM6wLtzunzeqAv1/np9v5Kci7wHuANVfWLWaZ1bl8t9F0Bi+kL+AvgAVr/gP838FvA5cDl7ecD\nfITW3Qo7OMhdFl3O9S7gPlp3TtwNvKxLub4KfLv9vue0x/phfx0qV1f2F7CJ1nWqfbTOLV8KHAd8\nGfgOrTvBjm3PXQncPGXb82kdMf098N5+yEXrDp1721/3dSnXH7Uf/xL4EbC1T/bXIXN1cn/Nkukh\nWtcz7ml/Xd3NfeXHk0iSGlvKp60kSfNkeUiSGrM8JEmNWR6SpMYsD0lSY5aHJKkxy0OS1JjlIXVJ\nksvz9Joi301yR68zSfPlLwlKXZZkOXA78IGq+kKv80jz4ZGH1H3/A7jd4tBidlSvA0hLSftTUJ9H\n6/O2pEXL01ZSlyR5Ma1V+15ZVY/3Oo90ODxtJXXPu4BjgTvaF80/1utA0nx55CFJaswjD0lSY5aH\nJKkxy0OS1JjlIUlqzPKQJDVmeUiSGrM8JEmN/X/U+iLJwcBhkAAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<matplotlib.figure.Figure at 0x10d2c6da0>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "z = array(linspace(8, 12, 9));\n",
    "\n",
    "xlabel('z'); ylabel('$\\sigma$(z)'); title('z vs. f(z)');\n",
    "scatter(z, func(a, b, z, q));\n",
    "legend('f(x)')\n",
    "show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "hideCode": true,
    "hidePrompt": true
   },
   "source": [
    "## Select Points to interpolate with. "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {
    "collapsed": false,
    "hideCode": true,
    "hidePrompt": true
   },
   "outputs": [],
   "source": [
    "z = np.array([9.5, 10, 10.5, 11])\n",
    "T = func(a, b, z, q)\n",
    "\n",
    "z = np.array([-9, -8, -7, -6])\n",
    "T = np.array([9.9, 11.7, 17.6, 18.2])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "hideCode": true,
    "hidePrompt": true
   },
   "source": [
    "Calculate $b_0$"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {
    "collapsed": false,
    "hideCode": true,
    "hidePrompt": true
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "9.9\n"
     ]
    }
   ],
   "source": [
    "b0 = T[0]\n",
    "print(b0)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "hideCode": true,
    "hidePrompt": true
   },
   "source": [
    "Calculate $b_1$"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {
    "collapsed": false,
    "hideCode": true,
    "hidePrompt": true
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "1.8\n"
     ]
    }
   ],
   "source": [
    "b1 = (T[1]-T[0])/(z[1]-z[0])\n",
    "print(b1)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "hideCode": true,
    "hidePrompt": true
   },
   "source": [
    "Calculate $b_2$"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {
    "collapsed": false,
    "hideCode": true,
    "hidePrompt": true
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "2.05\n"
     ]
    }
   ],
   "source": [
    "#T[z2, z1]\n",
    "var1 = (T[2]-T[1])/(z[2]-z[1]) #RIGHT\n",
    "\n",
    "#T[z1, z0]\n",
    "var2 = (T[1]-T[0])/(z[1]-z[0]) #RIGHT\n",
    "\n",
    "b2 = (var1-var2)/(z[2]-z[0])\n",
    "print(b2)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "hideCode": true,
    "hidePrompt": true
   },
   "source": [
    "Calculate $b_3$"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {
    "collapsed": false,
    "hideCode": true,
    "hidePrompt": true
   },
   "outputs": [
    {
     "data": {
      "text/latex": [
       "Therefore $b_3$ is: -1.56666666667"
      ],
      "text/plain": [
       "<IPython.core.display.Latex object>"
      ]
     },
     "execution_count": 9,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "#T[z3 ,z2, z1]\n",
    "var1 = (T[3]-T[2])/(z[3]-z[2])\n",
    "var2 = (T[2]-T[1])/(z[2]-z[1])\n",
    "\n",
    "temp1 = (var1-var2)/(z[3]-z[1])\n",
    "\n",
    "#T[z2 ,z1, z0]\n",
    "var1 = (T[2]-T[1])/(z[2]-z[1])\n",
    "var2 = (T[1]-T[0])/(z[1]-z[0])\n",
    "\n",
    "temp2 = (var1-var2)/(z[2]-z[0])\n",
    "\n",
    "b3 = (temp1-temp2)/(z[3]-z[0])\n",
    "Latex(r\"Therefore $b_3$ is: \" + str(b3))"
   ]
  }
 ],
 "metadata": {
  "celltoolbar": "Hide code",
  "hide_code_all_hidden": true,
  "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.6.0"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}