{
"cells": [
 {
		   "cell_type": "markdown",
	   "metadata": {},
	   "source": [
       "# Chapter 9: Curve Fitting Interpolation"
	   ]
	},
{
		   "cell_type": "markdown",
		   "metadata": {},
		   "source": [
			"## Example 9.10: Splines.sce"
		   ]
		  },
  {
"cell_type": "code",
	   "execution_count": null,
	   "metadata": {
	    "collapsed": true
	   },
	   "outputs": [],
"source": [
"//Example No. 9_10\n",
"//Splines\n",
"//Pg No. 301\n",
"clear ; close ; clc ;\n",
"\n",
"x = poly(0,'x');\n",
"function isitspline(f1,f2,f3,x0,x1,x2,x3)\n",
"    n1 = degree(f1),n2 = degree(f2),n3 = degree(f3)\n",
"    n = max(n1,n2,n3)\n",
"    f1_x1 = horner(f1,x1)\n",
"    f2_x1 = horner(f2,x1)\n",
"    f2_x2 = horner(f2,x2)\n",
"    f3_x2 = horner(f3,x2)\n",
"    if n ==1 & f1_x1 == f2_x1 & f2_x2 == f3_x2 then\n",
"        printf('The piecewise polynomials are continuous and f(x) is a linear spline')\n",
"    elseif f1_x1 == f2_x1 & f2_x2 == f3_x2\n",
"        for i = 1:n-1\n",
"            df1 = derivat(f1)\n",
"            df2 = derivat(f2)\n",
"            df3 = derivat(f3)\n",
"            df1_x1 = horner(df1,x1)\n",
"            df2_x1 = horner(df2,x1)\n",
"            df2_x2 = horner(df2,x2)\n",
"            df3_x2 = horner(df3,x2)\n",
"            f1 = df1,  f2 = df2,  f3 = df3\n",
"            if df1_x1 ~= df2_x1 | df2_x2 ~= df3_x2 then\n",
"                printf('The %ith derivative of polynomial is not continuours',i)\n",
"                break\n",
"            end\n",
"        end\n",
"        if i == n-1 & df1_x1 == df2_x1 & df2_x2 == df3_x2 then\n",
"            printf('The polynomial is continuous and its derivatives from 1 to %i are continuous, f(x) is a %ith degree polynomial',i,i+1)\n",
"        end\n",
"    else\n",
"            printf('The polynomial is not continuous')\n",
"    end\n",
"    \n",
"endfunction\n",
"n = 4 , x0 = -1 , x1 = 0, x2 = 1 , x3 = 2\n",
"f1 = x+1 ;\n",
"f2 = 2*x + 1 ;\n",
"f3 = 4 - x ;\n",
"disp('case 1')\n",
"isitspline(f1,f2,f3,x0,x1,x2,x3)\n",
"n = 4, x0 = 0 , x1= 1 , x2 = 2 , x3 = 3\n",
"f1 = x^2 + 1 ;\n",
"f2 = 2*x^2 ;\n",
"f3 = 5*x - 2 ;\n",
"disp('case 2')\n",
"isitspline(f1,f2,f3,x0,x1,x2,x3)\n",
"n = 4, x0 = 0, x1 = 1, x2 = 2, x3 = 3\n",
"f1 = x,\n",
"f2 = x^2 - x + 1,\n",
"f3 = 3*x - 3\n",
"disp('case 3')\n",
"isitspline(f1,f2,f3,x0,x1,x2,x3)"
   ]
   }
,
{
		   "cell_type": "markdown",
		   "metadata": {},
		   "source": [
			"## Example 9.11: Cubic_Spline_Interpolation.sce"
		   ]
		  },
  {
"cell_type": "code",
	   "execution_count": null,
	   "metadata": {
	    "collapsed": true
	   },
	   "outputs": [],
"source": [
"//Example No. 9_11\n",
"//Cubic Spline Interpolation\n",
"//Pg No. 306\n",
"clear ; close ; clc ;\n",
"\n",
"X = [ 4 9 16]\n",
"Fx = [ 2 3 4]\n",
"n = length(X)\n",
"h = diff(X)\n",
"disp(h)\n",
"x = poly(0,'x');\n",
"A(1) = 0;\n",
"A(n) = 0;\n",
"\n",
"//Since we do not know only a1 = A(2) we just have one equation which can be solved directly without solving tridiagonal matrix\n",
"A(2) = 6*( ( Fx(3) - Fx(2) )/h(2) - ( Fx(2) - Fx(1) )/h(1) )/( 2*( h(1) + h(2) ) );\n",
"disp(A(2),'a1 = ');\n",
"xuk = 7;\n",
"for i = 1:n-1\n",
"    if xuk <= X(i+1) then\n",
"        break;\n",
"    end\n",
"end\n",
"u = x*ones(1,2) - X(i:i+1)\n",
"s = (  A(2)*( u(i)^3 - ( h(i)^2 )*u(i))/6*h(i)  ) + ( Fx(i+1)*u(i)- Fx(i)*u(i+1) )/h(i);\n",
"disp(s,'s(x) = ');\n",
"s_7 = horner(s,xuk);\n",
"disp(s_7,'s(7)')"
   ]
   }
,
{
		   "cell_type": "markdown",
		   "metadata": {},
		   "source": [
			"## Example 9.12: Cubic_Spline_Interpolation.sce"
		   ]
		  },
  {
"cell_type": "code",
	   "execution_count": null,
	   "metadata": {
	    "collapsed": true
	   },
	   "outputs": [],
"source": [
"//Example No. 9_12\n",
"//Cubic Spline Interpolation\n",
"//Pg No. 313\n",
"clear ; close ;clc ;\n",
"\n",
"X = 1:4 ;\n",
"Fx = [ 0.5 0.3333 0.25 0.2]\n",
"n = length(X)\n",
"h = diff(X)\n",
"disp(h)\n",
"x = poly(0,'x');\n",
"A(1) = 0;\n",
"A(n) = 0;\n",
"//Forming Tridiagonal Matrix\n",
"//take make diagonal below main diagonal be 1 , main diagonal is 2 and diagonal above main diagonal is 3\n",
"diag1 = h(2:n-2);\n",
"diag2 = 2*(h(1:n-2)+h(2:n-1));\n",
"diag3 = h(2:n-2);\n",
"TridiagMat = diag(diag1,-1)+diag(diag2)+diag(diag3,1)\n",
"disp(TridiagMat);\n",
"D = diff(Fx);\n",
"D = 6*diff(D./h);\n",
"disp(D) \n",
"A(2:n-1) = TridiagMat\D' \n",
"disp(A)\n",
"xuk = 2.5;\n",
"for i = 1:n-1\n",
"    if xuk <= X(i+1) then\n",
"        break;\n",
"    end\n",
"end\n",
"u = x*ones(1,2) - X(i:i+1)\n",
"s = (  A(i)*( h(i+1)^2*u(2) - u(2)^2 )/( 6*h(i+1) ) ) + (  A(i+1)*( u(1)^3 - ( h(i)^2 )*u(1))/6*h(i)  ) + ( Fx(i+1)*u(1)- Fx(i)*u(2) )/h(i);\n",
"disp(s,'s(x) = ');\n",
"s2_5 = horner(s,xuk);\n",
"disp(s2_5,'s(2.5)') "
   ]
   }
,
{
		   "cell_type": "markdown",
		   "metadata": {},
		   "source": [
			"## Example 9.1: Polynomial_Forms.sce"
		   ]
		  },
  {
"cell_type": "code",
	   "execution_count": null,
	   "metadata": {
	    "collapsed": true
	   },
	   "outputs": [],
"source": [
"//Example No. 9_01\n",
"//Pg No.277\n",
"clear ; close ; clc ;\n",
"\n",
"printf('solving linear equations \n a0 + 100a1 = 3/7 \n a0 + 101a1 = -4/7 \n we get,\n');\n",
"C = [ 1 100 ; 1 101] \n",
"p = [ 3/7 ; -4/7] \n",
"a = C\p \n",
"printf('\n a0 = %f \n a1 = %f \n',a(1),a(2));\n",
"x = poly(0,'x') ;\n",
"px = a(1) + a(2)*x \n",
"p100 = horner(px,100)\n",
"p101 = horner(px,101)\n",
"printf('\n p(100) = %f \n p(101) = %f\n',p100,p101)"
   ]
   }
,
{
		   "cell_type": "markdown",
		   "metadata": {},
		   "source": [
			"## Example 9.2: Shifted_Power_form.sce"
		   ]
		  },
  {
"cell_type": "code",
	   "execution_count": null,
	   "metadata": {
	    "collapsed": true
	   },
	   "outputs": [],
"source": [
"//Example No. 9_02\n",
"//Page No. 278\n",
"clear ; close ; clc ;\n",
"\n",
"C = [ 1 100-100 ; 1 101-100] \n",
"p = [ 3/7 ; -4/7] \n",
"a = C\p \n",
"printf('\n a0 = %f \n a1 = %f \n',a(1),a(2));\n",
"x = poly(0,'x') ;\n",
"px = a(1) + a(2)*(x - 100) \n",
"p100 = horner(px,100)\n",
"p101 = horner(px,101)\n",
"printf('\n p(100) = %f \n p(101) = %f\n',p100,p101)"
   ]
   }
,
{
		   "cell_type": "markdown",
		   "metadata": {},
		   "source": [
			"## Example 9.3: Linear_Interpolation.sce"
		   ]
		  },
  {
"cell_type": "code",
	   "execution_count": null,
	   "metadata": {
	    "collapsed": true
	   },
	   "outputs": [],
"source": [
"//Example No. 9_03\n",
"//Page No. 280\n",
"clear ; close ; clc ;\n",
"\n",
"x = 1:5\n",
"f = [1 1.4142 1.7321 2 2.2361]\n",
"n = 2.5\n",
"for i = 1:5\n",
"    if n <= x(i) then\n",
"        break ;\n",
"    end\n",
"end\n",
"printf('%f lies between points %i and %i',n,x(i-1),x(i))\n",
"f2_5 = f(i-1) + ( n - x(i-1) )*( f(i) - f(i-1) )/( x(i) - x(i-1) )\n",
"err1 = 1.5811 - f2_5\n",
"disp(f2_5,'f(2.5) = ')\n",
"disp(err1,'error1 = ')\n",
"disp('The correct answer is 1.5811.The difference between results is due to use of a linear model to a nonlinear use')\n",
"disp('repeating the procedure using x1 = 2 and x2 = 4')\n",
"x1 = 2\n",
"x2 = 4\n",
"f2_5 = f(x1) + ( 2.5 - x1 )*( f(x2) - f(x1) )/( x2 - x1 )\n",
"err2 = 1.5811 - f2_5\n",
"disp(err2,'error2 = ')\n",
"disp(f2_5,'f(2.5) = ')\n",
"disp('NOTE- The increase in error due to the increase in the interval between the interpolating data')"
   ]
   }
,
{
		   "cell_type": "markdown",
		   "metadata": {},
		   "source": [
			"## Example 9.4: Lagrange_Interpolation.sce"
		   ]
		  },
  {
"cell_type": "code",
	   "execution_count": null,
	   "metadata": {
	    "collapsed": true
	   },
	   "outputs": [],
"source": [
"//Example No. 9_04\n",
"//Lagrange Interpolation- Second order\n",
"//Pg No. 282\n",
"clear ; close ; clc ;\n",
"\n",
"X = [ 1 2 3 4 5]\n",
"Fx = [ 1 1.4142 1.7321 2 2.2361];\n",
"X = X(2:4)\n",
"Fx = Fx(2:4)\n",
"x0 = 2.5 \n",
"x = poly(0,'x')\n",
"p = 0 \n",
"for i = 1:3\n",
"    L(i) = 1\n",
"    for j = 1:3\n",
"        if j == i then\n",
"            continue ;\n",
"        else\n",
"            L(i) = L(i)*( x - X(j) )/( X(i) - X(j) ) ;\n",
"        end\n",
"    end\n",
"    p = p + Fx(i)*L(i) \n",
"end\n",
"L0 = horner(L(1),2.5) ;\n",
"L1 = horner(L(2),2.5) ;\n",
"L2 = horner(L(3),2.5) ;\n",
"p2_5 = horner(p,2.5) ;\n",
"printf('For x = 2.5 we have,\n L0(2.5) = %f \n L1(2.5) = %f \n L2(2.5) = %f \n p(2.5) = %f \n',L0,L1,L2,p2_5)\n",
"\n",
"err = sqrt(2.5) - p2_5 ;\n",
"printf('The error is %f',err);\n",
""
   ]
   }
,
{
		   "cell_type": "markdown",
		   "metadata": {},
		   "source": [
			"## Example 9.5: Lagrange_Interpolation.sce"
		   ]
		  },
  {
"cell_type": "code",
	   "execution_count": null,
	   "metadata": {
	    "collapsed": true
	   },
	   "outputs": [],
"source": [
"//Example No. 9_05\n",
"//Lagrange Interpolation\n",
"//Pg No. 283\n",
"clear ; close ; clc ;\n",
"\n",
"i = [ 0 1 2 3 ]\n",
"X = [ 0 1 2 3 ]\n",
"Fx = [ 0 1.7183 6.3891 19.0855 ]\n",
"x = poly(0,'x');\n",
"n = 3 //order of lagrange polynomial \n",
"p = 0 \n",
"for i = 1:n+1\n",
"    L(i) = 1\n",
"    for j = 1:n+1\n",
"        if j == i then\n",
"            continue ;\n",
"        else\n",
"            L(i) = L(i)*( x - X(j) )/( X(i) - X(j) ) ;\n",
"        end\n",
"    end\n",
"    p = p + Fx(i)*L(i) \n",
"end\n",
"disp('The Lagrange basis polynomials are')\n",
"for i = 1:4\n",
"        disp(string(L(i)))\n",
"end\n",
"disp('The interpolation polynomial is ')\n",
"disp(string(p))\n",
"disp('The interpolation value at x = 1.5 is ' )\n",
"p1_5 = horner(p,1.5);\n",
"e1_5 = p1_5 + 1 ;\n",
"disp(e1_5,'e^1.5 = ',p1_5);\n",
"\n",
" "
   ]
   }
,
{
		   "cell_type": "markdown",
		   "metadata": {},
		   "source": [
			"## Example 9.6: Newton_Interpolation.sce"
		   ]
		  },
  {
"cell_type": "code",
	   "execution_count": null,
	   "metadata": {
	    "collapsed": true
	   },
	   "outputs": [],
"source": [
"//Example No. 9_06\n",
"//Newton Interpolation - Second order\n",
"//Pg No. 288\n",
"clear ; close ; clc ;\n",
"\n",
"i = [ 0 1 2 3]\n",
"X = 1:4 \n",
"Fx = [ 0 0.3010 0.4771 0.6021] \n",
"X = 1:3\n",
"Fx = Fx(1:3)\n",
"x = poly(0,'x');\n",
"A = Fx'\n",
"for i = 2:3\n",
"    for j = 1:4-i\n",
"        A(j,i) = ( A(j+1,i-1)-A(j,i-1) )/(X(j+i-1)-X(j)) ;\n",
"    end\n",
"end\n",
"printf('The coefficients of the polynomial are,\n a0 = %.4G \n a1 = %.4G \n a2 = %.4G \n',A(1,1),A(1,2),A(1,3))\n",
"p = A(1,1);\n",
"for i = 2:3\n",
"    p = p +A(1,i)* prod(x*ones(1,i-1) - X(1:i-1));\n",
"end\n",
"disp(string(p))\n",
"p2_5 = horner(p,2.5)\n",
"printf('p(2.5) = %.4G \n',p2_5)"
   ]
   }
,
{
		   "cell_type": "markdown",
		   "metadata": {},
		   "source": [
			"## Example 9.7: Newton_Divided_Difference_Interpolation.sce"
		   ]
		  },
  {
"cell_type": "code",
	   "execution_count": null,
	   "metadata": {
	    "collapsed": true
	   },
	   "outputs": [],
"source": [
"//Example No. 9_07\n",
"//Newton Divided Difference Interpolation\n",
"//Pg No. 291\n",
"clear ; close ; clc ;\n",
"\n",
"i = 0:4\n",
"X = 1:5\n",
"Fx = [ 0 7 26 63 124];\n",
"x = poly(0,'x');\n",
"A = [ i' X' Fx']\n",
"for i = 4:7\n",
"     for j = 1:8-i\n",
"        A(j,i) = ( A(j+1,i-1)-A(j,i-1) )/(X(j+i-3)-X(j)) ;\n",
"    end\n",
"end\n",
"disp(A)\n",
"p = A(1,3);\n",
"p1_5(1) = p ;\n",
"for i = 4:7\n",
"    p = p +A(1,i)* prod(x*ones(1,i-3) - X(1:i-3));\n",
"    p1_5(i-2) = horner(p,1.5);\n",
"end\n",
"printf('p0(1.5) = %f \n p1(1.5) = %f \n p2(1.5) = %f \n p3(1.5) = %f \n p4(1.5) = %f \n',p1_5(1),p1_5(2),p1_5(3),p1_5(4),p1_5(5));\n",
"disp(p1_5(5),'The function value at x = 1.5 is')"
   ]
   }
,
{
		   "cell_type": "markdown",
		   "metadata": {},
		   "source": [
			"## Example 9.8: Newton_Gregory_Forward_Difference_Formula.sce"
		   ]
		  },
  {
"cell_type": "code",
	   "execution_count": null,
	   "metadata": {
	    "collapsed": true
	   },
	   "outputs": [],
"source": [
"//Example No. 9_08\n",
"//Newton-Gregory forward difference formula\n",
"//Pg No. 297\n",
"clear ; close ; clc ;\n",
"\n",
"X = [ 10 20 30 40 50]\n",
"Fx = [ 0.1736 0.3420 0.5000 0.6428 0.7660]\n",
"x = poly(0,'x');\n",
"A = [X' Fx'];\n",
"for i = 3:6\n",
"     A(1:7-i,i) = diff(A(1:8-i,i-1))\n",
"end\n",
"disp(A)\n",
"x0 = X(1);\n",
"h = X(2) - X(1) ;\n",
"x1 = 25\n",
"s = (x1 - x0)/h ;\n",
"p(1) = Fx(1); \n",
"for j = 1:4\n",
"    p(j+1) = p(j) + prod(s*ones(1,j)-[0:j-1])*A(1,j+2)/factorial(j)\n",
"end\n",
"printf('p1(s) = %.4G \n p2(s) = %.4G \n p3(s) = %.4G \n p4(s) = %.4G \n',p(2),p(3),p(4),p(5)) \n",
"printf(' Thus sin(%d) = %.4G \n ',x1,p(5)) \n",
""
   ]
   }
,
{
		   "cell_type": "markdown",
		   "metadata": {},
		   "source": [
			"## Example 9.9: Newton_Backward_Difference_Formula.sce"
		   ]
		  },
  {
"cell_type": "code",
	   "execution_count": null,
	   "metadata": {
	    "collapsed": true
	   },
	   "outputs": [],
"source": [
"//Example No. 9_09\n",
"//Newton Backward difference formula\n",
"//Pg No. 299\n",
"clear ;close ;clc ;\n",
"\n",
"X = [ 10 20 30 40 50]\n",
"Fx = [ 0.1736 0.3420 0.5000 0.6428 0.7660]\n",
"x = poly(0,'x');\n",
"A = [X' Fx'];\n",
"for i = 3:6\n",
"     A(i-1:5,i) = diff(A(i-2:5,i-1))\n",
"end\n",
"disp(A);\n",
"xn = X(5);\n",
"h = 10 ;\n",
"xuk = 25;\n",
"s = (xuk - xn)/h ;\n",
"disp(s,'s = ');\n",
"p(1) = Fx(5)\n",
"for j = 1:4\n",
"    p(j+1) = p(j) + prod(s*ones(1,j)-[0:j-1])*A(5,j+2)/factorial(j)\n",
"end\n",
"printf('\n\n p1(s) = %.4G \n p2(s) = %.4G \n p3(s) = %.4G \n p4(s) = %.4G \n',p(2),p(3),p(4),p(5)) \n",
"printf(' Thus sin(%d) = %.4G \n ',xuk,p(5)) "
   ]
   }
],
"metadata": {
		  "kernelspec": {
		   "display_name": "Scilab",
		   "language": "scilab",
		   "name": "scilab"
		  },
		  "language_info": {
		   "file_extension": ".sce",
		   "help_links": [
			{
			 "text": "MetaKernel Magics",
			 "url": "https://github.com/calysto/metakernel/blob/master/metakernel/magics/README.md"
			}
		   ],
		   "mimetype": "text/x-octave",
		   "name": "scilab",
		   "version": "0.7.1"
		  }
		 },
		 "nbformat": 4,
		 "nbformat_minor": 0
}