24 files changed, 672 insertions, 0 deletions
diff --git a/845/CH6/EX6.1/Ex6_1.sce b/845/CH6/EX6.1/Ex6_1.sce
new file mode 100755
index 000000000..70f9a206d
--- /dev/null
+++ b/845/CH6/EX6.1/Ex6_1.sce
@@ -0,0 +1,20 @@
+//Example 6.1
+
+clc
+clear
+
+x = 0.1:0.2:1.3;
+y = [0.003 0.067 0.148 0.248 0.37 0.518 0.697];
+
+n = length(x);
+del = %nan*ones(n,6);
+del(:,1) = y';
+for j = 2:6
+    for i = 1:n-j+1
+        del(i,j) = del(i+1,j-1) - del(i,j-1);
+    end
+end
+del = [x' del];
+del = round(del*10^3)/10^3;
+mprintf("%5s %7s %8s %9s %8s %8s %8s",'x','y','dy','d2y','d3y','d4y','d5y')
+disp(del)
diff --git a/845/CH6/EX6.10/Ex6_10.sce b/845/CH6/EX6.10/Ex6_10.sce
new file mode 100755
index 000000000..44ebd653d
--- /dev/null
+++ b/845/CH6/EX6.10/Ex6_10.sce
@@ -0,0 +1,25 @@
+//Example 6.10
+
+clc
+clear
+
+x = 1:5;
+Y = poly(0, "Y");
+y = [2 5 7 Y 32];
+
+n = length(x);
+del = %nan*ones(n,5);
+del(:,1) = y';
+for j = 2:5
+    for i = 1:n-j+1
+        del(i,j) = del(i+1,j-1) - del(i,j-1);
+    end
+end
+del(:,1) = [];
+
+// del4 = 0
+
+y0 = del(:,4);
+y0(isnan(y0)) = [];
+Y = roots(y0)
+disp(Y,"Missing value f(x3) = ")
diff --git a/845/CH6/EX6.11/Ex6_11.sce b/845/CH6/EX6.11/Ex6_11.sce
new file mode 100755
index 000000000..4159cc4ee
--- /dev/null
+++ b/845/CH6/EX6.11/Ex6_11.sce
@@ -0,0 +1,35 @@
+//Example 6.11
+
+clc
+clear
+
+x = 0:5;
+y = [-3 3 11 27 57 107];
+
+n = length(x);
+del = %nan*ones(n,4);
+del(:,1) = y';
+for j = 2:4
+    for i = 1:n-j+1
+        del(i,j) = del(i+1,j-1) - del(i,j-1);
+    end
+end
+del(:,1) = [];
+
+X = poly(0, "x");
+h = x(2) - x(1);
+p = (X-x(1)) / h;
+x0 = x(1);
+y0 = y(1);
+dely0 = del(1,:);
+
+Y = y0;
+
+for i = 1:length(dely0)
+    t = 1;
+    for j = 1:i
+        t = t * (p-j+1);
+    end
+    Y = Y + t*dely0(i)/factorial(i);
+end
+disp(Y,"Required cubic polynomial:")
diff --git a/845/CH6/EX6.12/Ex6_12.sce b/845/CH6/EX6.12/Ex6_12.sce
new file mode 100755
index 000000000..2c33c5220
--- /dev/null
+++ b/845/CH6/EX6.12/Ex6_12.sce
@@ -0,0 +1,35 @@
+//Example 6.12
+
+clc
+clear
+
+x = 1:8;
+y = x^3;
+
+n = length(x);
+del = %nan*ones(n,4);
+del(:,1) = y';
+for j = 2:4
+    for i = 1:n-j+1
+        del(i+j-1,j) = del(i+j-1,j-1) - del(i+j-2,j-1);
+    end
+end
+
+X = 7.5;
+h = x(2) - x(1);
+p = (X-x(n)) / h;
+xn = x(n);
+yn = y(n);
+delyn = del(n,:);
+
+Y = 0;
+
+for i = 0:length(delyn)-1
+    t = 1;
+    for j = 1:i
+        t = t * (p+j-1);
+    end
+    Y = Y + t*delyn(i+1)/factorial(i);
+end
+disp(Y,"y(7.5) = ")
+
diff --git a/845/CH6/EX6.13/Ex6_13.sce b/845/CH6/EX6.13/Ex6_13.sce
new file mode 100755
index 000000000..b7c5b84e5
--- /dev/null
+++ b/845/CH6/EX6.13/Ex6_13.sce
@@ -0,0 +1,35 @@
+//Example 6.13
+
+clc
+clear
+
+x = 1974:2:1982;
+y = [40 43 48 52 57];
+
+n = length(x);
+del = %nan*ones(n,5);
+del(:,1) = y';
+for j = 2:5
+    for i = 1:n-j+1
+        del(i+j-1,j) = del(i+j-1,j-1) - del(i+j-2,j-1);
+    end
+end
+
+X = 1979;
+h = x(2) - x(1);
+p = (X-x(n)) / h;
+xn = x(n);
+yn = y(n);
+delyn = del(n,:);
+
+Y = 0;
+
+for i = 0:length(delyn)-1
+    t = 1;
+    for j = 1:i
+        t = t * (p+j-1);
+    end
+    Y = Y + t*delyn(i+1)/factorial(i);
+end
+Y = round(Y*10^4)/10^4;
+disp(Y,"Estimated sales for the year 1979: ")
diff --git a/845/CH6/EX6.14/Ex6_14.sce b/845/CH6/EX6.14/Ex6_14.sce
new file mode 100755
index 000000000..dd0d0a424
--- /dev/null
+++ b/845/CH6/EX6.14/Ex6_14.sce
@@ -0,0 +1,23 @@
+//Example 6.14
+
+clc
+clear
+
+x = [1 3 4 6];
+y = [-3 0 30 132];
+
+n = length(x);
+Y = 0;
+X = poly(0, "X");
+//X = 5;
+for i = 1:n
+    t = x;
+    t(i) = [];
+    p = 1;
+    for j = 1:length(t)
+        p = p * (X-t(j))/(x(i)-t(j));
+    end
+    Y = Y + p*y(i);
+end
+Y5 = horner(Y,5);
+disp(Y5,"y(5) = ")
diff --git a/845/CH6/EX6.15/Ex6_15.sce b/845/CH6/EX6.15/Ex6_15.sce
new file mode 100755
index 000000000..05dca459d
--- /dev/null
+++ b/845/CH6/EX6.15/Ex6_15.sce
@@ -0,0 +1,23 @@
+//Example 6.15
+
+clc
+clear
+
+x = [1 2 5];
+y = [1 4 10];
+
+n = length(x);
+Y = 0;
+X = poly(0, "X");
+//X = 5;
+for i = 1:n
+    t = x;
+    t(i) = [];
+    p = 1;
+    for j = 1:length(t)
+        p = p * (X-t(j))/(x(i)-t(j));
+    end
+    Y = Y + p*y(i);
+end
+Y5 = horner(Y,3);
+disp(Y5,"f(3) = ")
diff --git a/845/CH6/EX6.16/Ex6_16.sce b/845/CH6/EX6.16/Ex6_16.sce
new file mode 100755
index 000000000..5d3fb665a
--- /dev/null
+++ b/845/CH6/EX6.16/Ex6_16.sce
@@ -0,0 +1,30 @@
+//Example 6.16
+
+clc
+clear
+
+x = [0 1 2 4];
+y = [1 1 2 5];
+
+n = length(x);
+del = %nan*ones(n,4);
+del(:,1) = y';
+for j = 2:4
+    for i = 1:n-j+1
+        del(i,j) = (del(i+1,j-1) - del(i,j-1)) / (x(i+j-1) - x(i));
+    end
+end
+del(:,1) = [];
+
+Y = 0;
+X = poly(0, "X");
+for i = 1:n
+    t = x;
+    t(i) = [];
+    p = 1;
+    for j = 1:length(t)
+        p = p * (X-t(j))/(x(i)-t(j));
+    end
+    Y = Y + p*y(i);
+end
+disp(round(Y*10^4)/10^4,"Interpolating polynomial:")
diff --git a/845/CH6/EX6.17/Ex6_17.sce b/845/CH6/EX6.17/Ex6_17.sce
new file mode 100755
index 000000000..b3a35ca08
--- /dev/null
+++ b/845/CH6/EX6.17/Ex6_17.sce
@@ -0,0 +1,30 @@
+//Example 6.17
+
+clc
+clear
+
+x = [0 1 4];
+y = [2 1 4];
+
+n = length(x);
+del = %nan*ones(n,3);
+del(:,1) = y';
+for j = 2:3
+    for i = 1:n-j+1
+        del(i,j) = (del(i+1,j-1) - del(i,j-1)) / (x(i+j-1) - x(i));
+    end
+end
+del(:,1) = [];
+
+Y = 0;
+X = 2;
+for i = 1:n
+    t = x;
+    t(i) = [];
+    p = 1;
+    for j = 1:length(t)
+        p = p * (X-t(j))/(x(i)-t(j));
+    end
+    Y = Y + p*y(i);
+end
+disp(Y,"y(2) = ")
diff --git a/845/CH6/EX6.18/Ex6_18.sce b/845/CH6/EX6.18/Ex6_18.sce
new file mode 100755
index 000000000..46bde6a7b
--- /dev/null
+++ b/845/CH6/EX6.18/Ex6_18.sce
@@ -0,0 +1,2 @@
+// Example 6.18
+// This is an analytical problem and need not be coded.
diff --git a/845/CH6/EX6.19/Ex6_19.sce b/845/CH6/EX6.19/Ex6_19.sce
new file mode 100755
index 000000000..11d5d818a
--- /dev/null
+++ b/845/CH6/EX6.19/Ex6_19.sce
@@ -0,0 +1,104 @@
+//Example 6.19
+
+clc
+clear
+
+x = 0:4;
+n = length(x);
+f = "X^2 + Y^2 - Y";
+tab = %nan*ones(n,5);
+
+for j = 0:4
+    fj = strsubst(f,'Y','j');
+    for i = 1:n
+        tab(i,j+1) = eval(strsubst(fj,'X','x(i)'));
+    end
+end
+//tab(:,1) = [];
+mprintf("%4s %6s %6s %6s %6s %6s",'x','y=0','y=1','y=2','y=3','y=4')
+disp([(0:4)' tab])
+tab2 = tab(2:4,2:4)';
+n1 = length(tab2(:,1));
+y = 2:4;
+
+del1 = %nan*ones(n1,3);
+del1(:,1) = tab2(:,1);
+for j = 2:4
+    for i = 1:n1-j+1
+        del1(i,j) = del1(i+1,j-1) - del1(i,j-1);
+    end
+end
+
+del2 = %nan*ones(n1,3);
+del2(:,1) = tab2(:,2);
+for j = 2:4
+    for i = 1:n1-j+1
+        del2(i,j) = del2(i+1,j-1) - del2(i,j-1);
+    end
+end
+
+del3 = %nan*ones(n1,3);
+del3(:,1) = tab2(:,3);
+for j = 2:4
+    for i = 1:n1-j+1
+        del3(i,j) = del3(i+1,j-1) - del3(i,j-1);
+    end
+end
+
+y0 = y(1);
+Y = 3.5;
+hy = y(2) - y(1);
+py = (Y-y0)/hy;
+
+f1y = 0;
+del1y0 = del1(1,:);
+for i = 0:length(del1y0)-1
+    t = 1;
+    for j = 1:i
+        t = t * (py-j+1);
+    end
+    f1y = f1y + t*del1y0(i+1)/factorial(i);
+end
+
+f2y = 0;
+del2y0 = del2(1,:);
+for i = 0:length(del2y0)-1
+    t = 1;
+    for j = 1:i
+        t = t * (py-j+1);
+    end
+    f2y = f2y + t*del2y0(i+1)/factorial(i);
+end
+
+f3y = 0;
+del3y0 = del3(1,:);
+for i = 0:length(del3y0)-1
+    t = 1;
+    for j = 1:i
+        t = t * (py-j+1);
+    end
+    f3y = f3y + t*del3y0(i+1)/factorial(i);
+end
+
+del = %nan*ones(n1,3);
+del(:,1) = [f1y; f2y; f3y];
+for j = 2:4
+    for i = 1:n1-j+1
+        del(i,j) = del(i+1,j-1) - del(i,j-1);
+    end
+end
+
+f = 0;
+X = 2.5;
+x0 = x(2);
+hx = x(2) - x(1);
+px = (X-x0)/hx;
+del0 = del(1,:)
+for i = 0:length(del0)-1
+    t = 1;
+    for j = 1:i
+        t = t * (px-j+1);
+    end
+    f = f + t*del0(i+1)/factorial(i);
+end
+disp(f,"f(2.5,3.5) = ")
diff --git a/845/CH6/EX6.2/Ex6_2.sce b/845/CH6/EX6.2/Ex6_2.sce
new file mode 100755
index 000000000..eba668f6c
--- /dev/null
+++ b/845/CH6/EX6.2/Ex6_2.sce
@@ -0,0 +1,2 @@
+// Example 6.2
+// This is an analytical problem and need not be coded.
diff --git a/845/CH6/EX6.20/Ex6_20.sce b/845/CH6/EX6.20/Ex6_20.sce
new file mode 100755
index 000000000..6f7718043
--- /dev/null
+++ b/845/CH6/EX6.20/Ex6_20.sce
@@ -0,0 +1,53 @@
+//Example 6.20
+
+clc
+clear
+
+function [p] = cubicsplin(x,y)
+// Fits point data to cubic spline fit
+
+n = length(x);
+a = y(1:n-1);   // Spline Initials
+
+M1 = zeros(3*(n-1));
+M2 = zeros(3*(n-1),1);
+// Point Substitutions
+for i = 1:n-1
+    M1(i,i) = x(i+1) - x(i);
+    M1(i,i+n-1) = (x(i+1) - x(i))^2;
+    M1(i,i+2*(n-1)) = (x(i+1) - x(i))^3;
+    M2(i) = y(i+1) - y(i);
+end
+
+// Knot equations
+for i = 1:n-2
+    // Derivative (S') continuity
+    M1(i+n-1,i) = 1;
+    M1(i+n-1,i+1) = -1;
+    M1(i+n-1,i+n-1) = 2*(x(i+1)-x(i));
+    M1(i+n-1,i+2*(n-1)) = 3*(x(i+1)-x(i))^2;
+    // S'' continuity
+    M1(i+2*n-3,i+n-1) = 2;
+    M1(i+2*n-3,i+n) = -2;
+    M1(i+2*n-3,i+2*(n-1)) = 6*(x(i+1)-x(i));
+end
+// Given BC
+M1(3*n-4,n) = 1;
+M1(3*n-3,2*n-2) = 1;
+M1(3*n-3,3*n-3) = 3*(3-2);
+
+var = M1\M2;
+var = round(var);
+b = var(1:n-1);
+c = var(n:2*(n-1));
+d = var(2*(n-1)+1:3*(n-1));
+p = [d c b a(:)];
+endfunction
+
+x = 0:3;
+y = [1 4 0 -2];
+p = cubicsplin(x,y);
+for i = 1:length(p(:,1))
+    disp(strcat(["S",string(i-1),"(x) ="]))
+    disp(poly(p(i,:),"X",["coeff"]))
+end
diff --git a/845/CH6/EX6.21/Ex6_21.sce b/845/CH6/EX6.21/Ex6_21.sce
new file mode 100755
index 000000000..07aa74cc1
--- /dev/null
+++ b/845/CH6/EX6.21/Ex6_21.sce
@@ -0,0 +1,57 @@
+//Example 6.21
+
+clc
+clear
+
+function [p] = cubicsplin(x,y)
+// Fits point data to cubic spline fit
+
+n = length(x);
+a = y(1:n-1);   // Spline Initials
+
+M1 = zeros(3*(n-1));
+M2 = zeros(3*(n-1),1);
+// Point Substitutions
+for i = 1:n-1
+    M1(i,i) = x(i+1) - x(i);
+    M1(i,i+n-1) = (x(i+1) - x(i))^2;
+    M1(i,i+2*(n-1)) = (x(i+1) - x(i))^3;
+    M2(i) = y(i+1) - y(i);
+end
+
+// Knot equations
+for i = 1:n-2
+    // Derivative (S') continuity
+    M1(i+n-1,i) = 1;
+    M1(i+n-1,i+1) = -1;
+    M1(i+n-1,i+n-1) = 2*(x(i+1)-x(i));
+    M1(i+n-1,i+2*(n-1)) = 3*(x(i+1)-x(i))^2;
+    // S'' continuity
+    M1(i+2*n-3,i+n-1) = 2;
+    M1(i+2*n-3,i+n) = -2;
+    M1(i+2*n-3,i+2*(n-1)) = 6*(x(i+1)-x(i));
+end
+// Given BC
+M1(3*n-4,1) = 1;
+M1(3*n-3,n-1) = 1;
+M1(3*n-3,2*n-2) = 2*(3-2);
+M1(3*n-3,3*n-3) = 3*(3-2)^2;
+M2(3*n-4) = 2;
+M2(3*n-3) = 2;
+
+var = M1\M2;
+var = round(var);
+b = var(1:n-1);
+c = var(n:2*(n-1));
+d = var(2*(n-1)+1:3*(n-1));
+
+p = [a(:) b c d];
+endfunction
+
+x = 0:3;
+y = [1 4 0 -2];
+p = cubicsplin(x,y);
+for i=1:length(p(:,1))
+    disp(strcat(["S",string(i-1),"(x) = "]))
+    disp(poly(p(i,:),"x",["coeff"]))
+end
diff --git a/845/CH6/EX6.22/Ex6_22.sce b/845/CH6/EX6.22/Ex6_22.sce
new file mode 100755
index 000000000..e1bf6e2bb
--- /dev/null
+++ b/845/CH6/EX6.22/Ex6_22.sce
@@ -0,0 +1,42 @@
+//Example 6.22
+
+clc
+clear
+
+x = 3:8;
+y = [0.205 0.24 0.259 0.262 0.25 0.224];
+
+n = length(x);
+del = %nan*ones(n,5);
+del(:,1) = y';
+for j = 2:5
+    for i = 1:n-j+1
+        del(i,j) = del(i+1,j-1) - del(i,j-1);
+    end
+end
+
+X = poly(0, "X");
+x0 = x(1);
+y0 = y(1);
+h = x(2) - x(1);
+p = (X-x0)/h;
+del0 = del(1,:);
+del0 = round(del0*10^4)/10^4;
+del0 = del0(1:find(del0==0)-1);
+
+Y = 0;
+for i = 0:length(del0)-1
+    t = 1;
+    for j = 1:i
+        t = t * (p-j+1);
+    end
+    Y = Y + t*del0(i+1)/factorial(i);
+end
+disp(Y,"y = ")
+
+dydx = derivat(Y);
+minx = roots(dydx);
+miny = round(horner(Y,minx)*10^5)/10^5;
+disp(minx,"min_x = ")
+disp(miny,"min_y = ")
+//min_y value is incorrectly displayed in textbook as 0.25425 instead of 0.26278
diff --git a/845/CH6/EX6.23/Ex6_23.sce b/845/CH6/EX6.23/Ex6_23.sce
new file mode 100755
index 000000000..bbf324893
--- /dev/null
+++ b/845/CH6/EX6.23/Ex6_23.sce
@@ -0,0 +1,36 @@
+//Example 6.23
+
+clc
+clear
+
+x = [-1 1 2 3];
+y = [-21 15 12 3];
+
+n = length(x);
+X = poly(0, "X");
+Y = 0;
+for i = 1:n
+    t = x;
+    t(i) = [];
+    p = 1;
+    for j = 1:length(t)
+        p = p * (X-t(j))/(x(i)-t(j));
+    end
+    Y = Y + p*y(i);
+end
+
+dydx = derivat(Y);
+extx = real(roots(dydx));
+extx = round(extx*10^4)/10^4;
+d2ydx = derivat(dydx);
+
+if  horner(d2ydx,extx(1)) < 0 then
+    maxx = extx(1);
+    maxy = horner(Y,maxx);
+else
+    maxx = extx(2);
+    maxy = horner(Y,maxx);
+end
+maxy = round(maxy*10^4)/10^4;
+disp(maxx,"max_x = ")
+disp(maxy,"max_y = ")
diff --git a/845/CH6/EX6.24/Ex6_24.sce b/845/CH6/EX6.24/Ex6_24.sce
new file mode 100755
index 000000000..02a9092e4
--- /dev/null
+++ b/845/CH6/EX6.24/Ex6_24.sce
@@ -0,0 +1,33 @@
+//Example 6.24
+
+clc
+clear
+
+x = 1:3:10;
+F = [500426 329240 175212 40365];
+
+n = length(x);
+del = %nan*ones(n,4);
+del(:,1) = F';
+for j = 2:4
+    for i = 1:n-j+1
+        del(i,j) = del(i+1,j-1) - del(i,j-1);
+    end
+end
+
+del0 = del(1,:);
+X = 2;
+x0 = x(1);
+h = x(2) - x(1);
+p = (X-x0) / h;
+F2 = 0;
+for i = 0:length(del0)-1
+    t = 1;
+    for j = 1:i
+        t = t * (p-j+1);
+    end
+    F2 = F2 + t*del0(i+1)/factorial(i);
+end
+
+f2 = F(1) - F2;
+disp(f2,"f(2) = ")
diff --git a/845/CH6/EX6.3/Ex6_3.sce b/845/CH6/EX6.3/Ex6_3.sce
new file mode 100755
index 000000000..318045c03
--- /dev/null
+++ b/845/CH6/EX6.3/Ex6_3.sce
@@ -0,0 +1,2 @@
+// Example 6.3
+// This is an analytical problem and need not be coded.
diff --git a/845/CH6/EX6.4/Ex6_4.sce b/845/CH6/EX6.4/Ex6_4.sce
new file mode 100755
index 000000000..9bfb8eccb
--- /dev/null
+++ b/845/CH6/EX6.4/Ex6_4.sce
@@ -0,0 +1,2 @@
+// Example 6.4
+// This is an analytical problem and need not be coded.
diff --git a/845/CH6/EX6.5/Ex6_5.sce b/845/CH6/EX6.5/Ex6_5.sce
new file mode 100755
index 000000000..157ff3344
--- /dev/null
+++ b/845/CH6/EX6.5/Ex6_5.sce
@@ -0,0 +1,2 @@
+// Example 6.5
+// This is an analytical problem and need not be coded.
diff --git a/845/CH6/EX6.6/Ex6_6.sce b/845/CH6/EX6.6/Ex6_6.sce
new file mode 100755
index 000000000..84c95049e
--- /dev/null
+++ b/845/CH6/EX6.6/Ex6_6.sce
@@ -0,0 +1,2 @@
+// Example 6.6
+// This is an analytical problem and need not be coded.
diff --git a/845/CH6/EX6.7/Ex6_7.sce b/845/CH6/EX6.7/Ex6_7.sce
new file mode 100755
index 000000000..88f5b9247
--- /dev/null
+++ b/845/CH6/EX6.7/Ex6_7.sce
@@ -0,0 +1,2 @@
+// Example 6.7
+// This is an analytical problem and need not be coded.
diff --git a/845/CH6/EX6.8/Ex6_8.sce b/845/CH6/EX6.8/Ex6_8.sce
new file mode 100755
index 000000000..5e4ae1d5a
--- /dev/null
+++ b/845/CH6/EX6.8/Ex6_8.sce
@@ -0,0 +1,41 @@
+//Example 6.8
+
+clc
+clear
+
+x = 10:10:50;
+y = [46 66 81 93 101];
+
+n = length(x);
+del = %nan*ones(n,5);
+del(:,1) = y';
+for j = 2:5
+    for i = 1:n-j+1
+        del(i,j) = del(i+1,j-1) - del(i,j-1);
+    end
+end
+del(:,1) = [];
+
+X = 15; //input
+for i = 1:n
+    if X>x(i) then
+        h = x(i+1) - x(i);
+        p = (X-x(i)) / h;
+        x0 = x(i);
+        y0 = y(i);
+        dely0 = del(i,:);
+        dely0(isnan(y0)) = [];
+    end
+end
+
+Y = y0;
+
+for i = 1:length(dely0)
+    t = 1;
+    for j = 1:i
+        t = t * (p-j+1);
+    end
+    Y = Y + t*dely0(i)/factorial(i);
+end
+Y = round(Y*10^4)/10^4;
+disp(Y,"f(15) = ")
diff --git a/845/CH6/EX6.9/Ex6_9.sce b/845/CH6/EX6.9/Ex6_9.sce
new file mode 100755
index 000000000..0f47285e1
--- /dev/null
+++ b/845/CH6/EX6.9/Ex6_9.sce
@@ -0,0 +1,36 @@
+//Example 6.9
+
+clc
+clear
+
+x = 0.1:0.1:0.5;
+y = [1.4 1.56 1.76 2 2.28];
+
+n = length(x);
+del = %nan*ones(n,5);
+del(:,1) = y';
+for j = 2:5
+    for i = 1:n-j+1
+        del(i,j) = del(i+1,j-1) - del(i,j-1);
+    end
+end
+del(:,1) = [];
+
+X = poly(0, "X");
+h = x(2) - x(1);
+p = (X-x(1)) / h;
+x0 = x(1);
+y0 = y(1);
+dely0 = del(1,:);
+
+Y = y0;
+
+for i = 1:length(dely0)
+    t = 1;
+    for j = 1:i
+        t = t * (p-j+1);
+    end
+    Y = Y + t*dely0(i)/factorial(i);
+end
+Y = round(Y*10^2)/10^2;
+disp(Y,"Required Newton''s Interpolating Polynomial:")