## Piecewise Interpolation: Quadratic Spline Interpolation

#### Scheme 1:

In the first scheme, the intervals between the data points are used as intervals on which a quadratic function is defined. On the intermediate nodes, the slope of the parabola on the left is equated to the slope of the parabola on the right (Figure 3). Consider a set of data points , with intervals. We seek to find a quadratic polynomial on the interval . In total, there are coefficients. On each interval, the two equations and ensure that the spline passes through the data points which also ensure the continuity of the splines. These continuity equations will provide equations. In addition, the connection of the quadratic polynomials at intermediate nodes will be assumed to be smooth, and so at every intermediate point the first derivative is assumed to be continuous, that is . These first-derivative continuity equations provide another equations. Finally, the second derivative at the starting point can be set to 0 to provide one additional condition so that the number of equations would be equal to which is equal to the number of unknowns. The following are the resulting equations:

Figure 3. Scheme 1 of piecewise quadratic interpolation with continuity in the interpolating function and its first derivative

The following Mathematica code implements this procedure for a set of data (the procedure will only work if ).

View Mathematica Code
Data = {{-1, 0.038}, {-0.8, 0.058}, {-0.60, 0.10}, {-0.4, 0.20}, {-0.2, 0.50}, {0, 1}, {0.2, 0.5}, {0.4, 0.2}, {0.6, 0.1}, {0.8, 0.058}, {1, 0.038}};
Spline2[Data_] := (k1 = Length[Data] - 1;
y = Table[If[i <= 2 k1,If[EvenQ[i], Data[[i/2 + 1, 2]], Data[[(i + 1)/2, 2]]], 0], {i,1, 3 k1}];
M = Table[0, {i, 1, 3 k1}, {j, 1, 3 k1}];
For[i = 1, i <= k1, i = i + 1,
M[[2 i - 1, 3 (i - 1) + 1]] = M[[2 i, 3 (i - 1) + 1]] = 1;
M[[2 i - 1, 3 (i - 1) + 2]] = Data[[i, 1]];
M[[2 i - 1, 3 (i - 1) + 1 + 2]] = Data[[i, 1]]^2;
M[[2 i, 3 (i - 1) + 1 + 1]] = Data[[i + 1, 1]];
M[[2 i, 3 (i - 1) + 1 + 2]] = Data[[i + 1, 1]]^2];
For[i = 1, i <= k1 - 1, i = i + 1,
M[[2 k1 + i, 3 (i - 1) + 1 + 1]] = 1;
M[[2 k1 + i, 3 (i - 1) + 1 + 4]] = -1;
M[[2 k1 + i, 3 (i - 1) + 1 + 2]] = 2 Data[[i + 1, 1]];
M[[2 k1 + i, 3 (i - 1) + 1 + 5]] = -2 Data[[i + 1, 1]]];
M[[3 k1, 3]] = 2;
Coef = Inverse[M].y;
pf = Table[{Coef[[3 (i - 1) + 1]] + Coef[[3 (i - 1) + 2]]*x + Coef[[3 (i - 1) + 3]]*x^2, Data[[i, 1]] <= x <= Data[[i + 1, 1]]}, {i, 1, k1}];
Piecewise[pf])
y = Spline2[Data]
M // MatrixForm
yactual=1/(1+25x^2);
a = Plot[{y, yactual}, {x, -1, 1},Epilog -> {PointSize[Large], Point[Data]}, AxesLabel -> {"x", "y1"},ImageSize -> Medium, PlotRange -> All,PlotLegends -> {"y2", "yactual"}, PlotLabel ->"Quadratic spline"]
View Python Code
import numpy as np
import matplotlib.pyplot as plt

Data = [[-1, 0.038], [-0.8, 0.058], [-0.60, 0.10], [-0.4, 0.20], [-0.2, 0.50], [0, 1], [0.2, 0.5], [0.4, 0.2], [0.6, 0.1], [0.8, 0.058], [1, 0.038]]
def Spline2(x, Data):
k1 = len(Data) - 1
y = [(Data[i//2][1] if i%2==0 else Data[(i)//2][1]) if i <= 2*k1 else 0 for i in range(1, 3*k1 + 1)]
M = np.zeros([3*k1, 3*k1])
for i in range(k1):
M[2*(i + 1) - 2][3*i] = M[2*(i + 1) - 1][3*i] = 1
M[2*(i + 1) - 2][3*i + 1] = Data[i][0]
M[2*(i + 1) - 2][3*i + 2] = Data[i][0]**2
M[2*(i + 1) - 1][3*i + 1] = Data[i + 1][0]
M[2*(i + 1) - 1][3*i + 2] = Data[i + 1][0]**2
for i in range(k1 - 1):
M[2*k1 + i][3*i + 1] = 1
M[2*k1 + i][3*i + 4] = -1
M[2*k1 + i][3*i + 2] = 2*Data[i + 1][0]
M[2*k1 + i][3*i + 5] = -2*Data[i + 1][0]
M[3*k1 - 1][2] = 2
print("M Matrix\n", M)
Coef = np.matmul(np.linalg.inv(M),y)
display([['{:30}'.format(str(round(Coef[3*i],3)) + ' + ' + str(round(Coef[3*i + 1],3)) + 'x + ' + str(round(Coef[3*i + 2],3)) + 'x**2'),
'{:15}'.format(str(round(Data[i][0],3))+' <=x<= '+str(round(Data[i + 1][0],3)))] for i in range(k1)])
return np.piecewise(x, [(x >= Data[i][0])&(x <= Data[i + 1][0]) for i in range(k1)],
[lambda x, j=i: Coef[3*j] + Coef[3*j + 1]*x + Coef[3*j + 2]*x**2 for i in range(k1)])
x = np.arange(-1,1,0.01)
y = Spline2(x, Data)
yactual = 1/(1 + 25*x**2)
plt.plot(x,yactual, label="yactual")
plt.plot(x, y, label="y2")
plt.scatter([point[0] for point in Data],[point[1] for point in Data], c='k')
plt.xlabel("x"); plt.ylabel("y2")
plt.legend(); plt.grid(); plt.show()

The following MATLAB code applies scheme 1 to data extracted from the Runge function:

##### Example

Use scheme 1 of quadratic interpolation to find the interpolating function given the following 4 data points (-1, 0.038) (-0.8, 0.058), (-0.60, 0.10), (-0.4, 0.20)

##### Solution

The four data points () with three intervals (), therefore, there are unknowns. The following are the nine equations according to scheme 1 to be solved to find the unknowns:

Solving the above equations produces the following interpolation function:

Figure 4 shows the resulting interpolation function along with the four data points. The Mathematica code is provided below the figure.

Figure 4. Scheme 1 quadratic interpolation applied to four data points

View Mathematica Code
Data = {{-1, 0.038}, {-0.8, 0.058}, {-0.60, 0.10}, {-0.4, 0.20}};
Spline2[Data_] := (k1 = Length[Data] - 1;
y = Table[
If[i <= 2 k1,
If[EvenQ[i], Data[[i/2 + 1, 2]], Data[[(i + 1)/2, 2]]], 0], {i,
1, 3 k1}];
M = Table[0, {i, 1, 3 k1}, {j, 1, 3 k1}];
For[i = 1, i <= k1, i = i + 1,
M[[2 i - 1, 3 (i - 1) + 1]] = M[[2 i, 3 (i - 1) + 1]] = 1;
M[[2 i - 1, 3 (i - 1) + 2]] = Data[[i, 1]];
M[[2 i - 1, 3 (i - 1) + 1 + 2]] = Data[[i, 1]]^2;
M[[2 i, 3 (i - 1) + 1 + 1]] = Data[[i + 1, 1]];
M[[2 i, 3 (i - 1) + 1 + 2]] = Data[[i + 1, 1]]^2];
For[i = 1, i <= k1 - 1, i = i + 1,
M[[2 k1 + i, 3 (i - 1) + 1 + 1]] = 1;
M[[2 k1 + i, 3 (i - 1) + 1 + 4]] = -1;
M[[2 k1 + i, 3 (i - 1) + 1 + 2]] = 2 Data[[i + 1, 1]];
M[[2 k1 + i, 3 (i - 1) + 1 + 5]] = -2 Data[[i + 1, 1]]];
M[[3 k1, 3]] = 2;
Coef = Inverse[M].y;
pf = Table[{Coef[[3 (i - 1) + 1]] + Coef[[3 (i - 1) + 2]]*x + Coef[[3 (i - 1) + 3]]*x^2, Data[[i, 1]] <= x <= Data[[i + 1, 1]]}, {i, 1, k1}];
Piecewise[pf])
y = Spline2[Data]
M // MatrixForm
a = Plot[{y}, {x, -1, -0.4}, Epilog -> {PointSize[Large], Point[Data]}, AxesLabel -> {"x", "y2"},ImageSize -> Medium, PlotRange -> All, PlotLegends -> {"y2"}, PlotLabel -> "Quadratic spline"]
View Python Code
import numpy as np
import matplotlib.pyplot as plt

Data = [[-1, 0.038], [-0.8, 0.058], [-0.60, 0.10], [-0.4, 0.20]]
def Spline2(x, Data):
k1 = len(Data) - 1
y = [(Data[i//2][1] if i%2==0 else Data[(i)//2][1]) if i <= 2*k1 else 0 for i in range(1, 3*k1 + 1)]
M = np.zeros([3*k1, 3*k1])
for i in range(k1):
M[2*(i + 1) - 2][3*i] = M[2*(i + 1) - 1][3*i] = 1
M[2*(i + 1) - 2][3*i + 1] = Data[i][0]
M[2*(i + 1) - 2][3*i + 2] = Data[i][0]**2
M[2*(i + 1) - 1][3*i + 1] = Data[i + 1][0]
M[2*(i + 1) - 1][3*i + 2] = Data[i + 1][0]**2
for i in range(k1 - 1):
M[2*k1 + i][3*i + 1] = 1
M[2*k1 + i][3*i + 4] = -1
M[2*k1 + i][3*i + 2] = 2*Data[i + 1][0]
M[2*k1 + i][3*i + 5] = -2*Data[i + 1][0]
M[3*k1 - 1][2] = 2
print("M Matrix\n", M)
Coef = np.matmul(np.linalg.inv(M),y)
display([['{:30}'.format(str(round(Coef[3*i],3)) + ' + ' + str(round(Coef[3*i + 1],3)) + 'x + ' + str(round(Coef[3*i + 2],3)) + 'x**2'),
'{:15}'.format(str(round(Data[i][0],3))+' <=x<= '+str(round(Data[i + 1][0],3)))] for i in range(k1)])
return np.piecewise(x, [(x >= Data[i][0])&(x <= Data[i + 1][0]) for i in range(k1)],
[lambda x, j=i: Coef[3*j] + Coef[3*j + 1]*x + Coef[3*j + 2]*x**2 for i in range(k1)])
x = np.arange(-1,-0.4,0.01)
y = Spline2(x, Data)
plt.plot(x, y, label="y2")
plt.scatter([point[0] for point in Data],[point[1] for point in Data], c='k')
plt.xlabel("x"); plt.ylabel("y2")
plt.legend(); plt.grid(); plt.show()

The following MATLAB code applies scheme 1 to the data from this example:

In certain cases, such as fluctuating data, this scheme might lead to some oscillations within the intervals between the data points. As an example, when this scheme is used to fit through the 11 data points of the Runge function given in the previous section, a nice smooth interpolating function is produced for the first few intervals (Figure 5). The condition of having zero second derivative on the initial interval causes the interpolating in the first interval to follow the line connecting the first two data points. However, in the last five intervals, the interpolating function oscillates and deviates away from the expected behaviour.

Figure 5. Behaviour of scheme 1 for quadratic interpolation of the Runge function data points.

#### Scheme 2

For the second scheme, the intervals where the quadratic interpolation functions are defined, are chosen different from the intervals between the data points in a manner to ensure symmetry of the quadratic interpolation. The new intervals where the piecewise quadratic functions are defined are bounded by the points (Figure 6). Consider a set of data points , with intervals. On each of the intervals bounded by a quadratic polynomial is defined. In total, there are coefficients. First,the quadratic polynomials have to pass through the data poitns resulting in equations. These quadratic polynomials have to be continuous and differentiable at the intermediate points that are the bounds of the intervals resulting in equations. Therefore, in total, there are equations which is equal to the number of unknowns. The following are the resulting equations:

Figure 6. Scheme 2 of piecewise quadratic interpolation with continuity in the interpolating function and its first derivative

The following Mathematica code implements this procedure for a set of data (the procedure will only work if ). As shown in Figure 7, the oscillations produced by scheme 1 disappear if this scheme is used for interpolating the data points of the Runge function given in the previous section.

View Mathematica Code
Data = {{-1, 0.038}, {-0.8, 0.058}, {-0.60, 0.10}, {-0.4,0.20}, {-0.2, 0.50}, {0, 1}, {0.2, 0.5}, {0.4, 0.2}, {0.6,0.1}, {0.8, 0.058}, {1, 0.038}};
Spline2[Data_] := (k = Length[Data];
NewData = Table[0, {i, 1, k - 3}];
Do[NewData[[i]] = (Data[[i + 1]] + Data[[i + 2]])/2, {i, 1, k - 3}];
NNewData = Prepend[Append[NewData, Data[[k]]], Data[[1]]];
y = Table[If[i <= k - 1 + 1, Data[[i, 2]], 0], {i, 1, 3 k - 6}];
M = Table[0, {i, 1, 3 k - 6}, {j, 1, 3 k - 6}];
M[[1, 1]] = M[[2, 1]] = 1;
M[[1, 2]] = Data[[1, 1]];
M[[1, 3]] = Data[[1, 1]]^2;
M[[2, 2]] = Data[[2, 1]];
M[[2, 3]] = Data[[2, 1]]^2;
M[[k - 1, 3 k - 8]] = M[[k, 3 k - 8]] = 1;
M[[k - 1, 3 k - 7]] = Data[[k - 1, 1]];
M[[k - 1, 3 k - 6]] = Data[[k - 1, 1]]^2;
M[[k, 3 k - 7]] = Data[[k, 1]];
M[[k, 3 k - 6]] = Data[[k, 1]]^2;
For[i = 3, i <= k - 2, i = i + 1, M[[i, 3 (i - 2) + 1]] = 1;
M[[i, 3 (i - 2) + 2]] = Data[[i, 1]];
M[[i, 3 (i - 2) + 3]] = Data[[i, 1]]^2];
For[i = 1, i <= k - 3, i = i + 1,
M[[k + i, 3 (i - 1) + 1]] = 1; M[[k + i, 3 (i - 1) + 4]] = -1;
M[[k + i, 3 (i - 1) + 2]] = NewData[[i, 1]];
M[[k + i, 3 (i - 1) + 5]] = -NewData[[i, 1]];
M[[k + i, 3 (i - 1) + 3]] = NewData[[i, 1]]^2;
M[[k + i, 3 (i - 1) + 6]] = -1*NewData[[i, 1]]^2;
M[[k + i + k - 3, 3 (i - 1) + 2]] = 1;
M[[k + i + k - 3, 3 (i - 1) + 5]] = -1;
M[[k + i + k - 3, 3 (i - 1) + 3]] = 2*NewData[[i, 1]];
M[[k + i + k - 3, 3 (i - 1) + 6]] = -2*NewData[[i, 1]]];
Coef = Inverse[M].y;
pf = Table[{Coef[[3 (i - 1) + 1]] + Coef[[3 (i - 1) + 2]]*x + Coef[[3 (i - 1) + 3]]*x^2,NNewData[[i, 1]] <= x <= NNewData[[i + 1, 1]]}, {i, 1, k - 2}];
Piecewise[pf])
y = Spline2[Data]
yactual = 1/(1 + 25 x^2);
M//MatrixForm
a = Plot[{y, yactual}, {x, -1, 1}, Epilog -> {PointSize[Large], Point[Data]}, AxesLabel -> {"x", "y2"},ImageSize -> Medium, PlotRange -> All, PlotLegends -> {"y2", "yactual"}, PlotLabel -> "Quadratic spline"]
View Python Code
import numpy as np
import matplotlib.pyplot as plt

Data = [[-1, 0.038], [-0.8, 0.058], [-0.60, 0.10], [-0.4, 0.20], [-0.2, 0.50], [0, 1], [0.2, 0.5], [0.4, 0.2], [0.6,0.1], [0.8, 0.058], [1, 0.038]]
def Spline2(x, Data):
k = len(Data)
Data = np.array(Data)
NewData = np.array([(Data[i + 1] + Data[i + 2])/2 for i in range(k-3)])
NNewData = np.vstack([Data[0], NewData, Data[k-1]])
y = [Data[i][1] if i < k else 0 for i in range(3*k - 6)]
M = np.zeros([3*k - 6, 3*k - 6])
M[0][0] = M[1][0] = 1
M[0][1] = Data[0][0]
M[0][2] = Data[0][0]**2
M[1][1] = Data[1][0]
M[1][2] = Data[1][0]**2
M[k - 2][3*k - 9] = M[k - 1][3*k - 9] = 1
M[k - 2][3*k - 8] = Data[k - 2][0]
M[k - 2][3*k - 7] = Data[k - 2][0]**2
M[k - 1][3*k - 8] = Data[k - 1][0]
M[k - 1][3*k - 7] = Data[k - 1][0]**2
for i in range(2, k-2):
M[i][3*(i - 1)] = 1
M[i][3*(i - 1) + 1] = Data[i][0]
M[i][3*(i - 1) + 2] = Data[i][0]**2
for i in range(k - 3):
M[k + i][3*i] = 1
M[k + i][3*i + 3] = -1
M[k + i][3*i + 1] = NewData[i][0]
M[k + i][3*i + 4] = -NewData[i][0]
M[k + i][3*i + 2] = NewData[i][0]**2
M[k + i][3*i + 5] = -1*NewData[i][0]**2
M[k + i + k - 3][3*i + 1] = 1
M[k + i + k - 3][3*i + 4] = -1
M[k + i + k - 3][3*i + 2] = 2*NewData[i][0]
M[k + i + k - 3][3*i + 5] = -2*NewData[i][0]
print("M Matrix\n", M)
Coef = np.matmul(np.linalg.inv(M),y)
display([['{:30}'.format(str(round(Coef[3*i],3)) + ' + ' + str(round(Coef[3*i + 1],3)) + 'x + ' + str(round(Coef[3*i + 2],3)) + 'x**2'),
'{:15}'.format(str(round(NNewData[i][0],3))+' <=x<= '+str(round(NNewData[i + 1][0],3)))] for i in range(k - 2)])
return np.piecewise(x, [(x >= NNewData[i][0])&(x <= NNewData[i + 1][0]) for i in range(k - 2)],
[lambda x, j=i: Coef[3*j] + Coef[3*j + 1]*x + Coef[3*j + 2]*x**2 for i in range(k - 2)])
x = np.arange(-1,1,0.01)
yactual = 1/(1 + 25*x**2)
y = Spline2(x, Data)
plt.plot(x, yactual, label="yactual")
plt.plot(x, y, label="y2")
plt.scatter([point[0] for point in Data],[point[1] for point in Data], c='k')
plt.xlabel("x"); plt.ylabel("y2")
plt.legend(); plt.grid(); plt.show()

The following MATLAB code applies scheme 2 to data extracted from the Runge function:

MATLAB file

Figure 7. Behaviour of scheme 2 of quadratic interpolation when applied to the Runge function data points

##### Example

Use scheme 2 of quadratic interpolation to find the interpolating function given the following 5 data points (-1, 0.038) (-0.8, 0.058), (-0.60, 0.10), (-0.4, 0.20), (-0.2, 0.5)

##### Solution

The five data points () with four intervals (). intervals for the quadratic interpolation functions are defined between the points . Therefore, there are unknowns. The following are the nine equations according to scheme 2 to be solved to find the unknowns:

Solving the above equations produces the following interpolation function:

Figure 8 shows the resulting interpolation function along with the five data points. The Mathematica code is provided below the figure.

Figure 8. Scheme 2 quadratic interpolation applied to five data points

View Mathematica Code
Data = {{-1, 0.038}, {-0.8, 0.058}, {-0.60, 0.10}, {-0.4,0.20}, {-0.2, 0.50}};
Spline2[Data_] := (k = Length[Data];
NewData = Table[0, {i, 1, k - 3}];
Do[NewData[[i]] = (Data[[i + 1]] + Data[[i + 2]])/2, {i, 1, k - 3}];
NNewData = Prepend[Append[NewData, Data[[k]]], Data[[1]]];
y = Table[If[i <= k - 1 + 1, Data[[i, 2]], 0], {i, 1, 3 k - 6}];
M = Table[0, {i, 1, 3 k - 6}, {j, 1, 3 k - 6}];
M[[1, 1]] = M[[2, 1]] = 1;
M[[1, 2]] = Data[[1, 1]];
M[[1, 3]] = Data[[1, 1]]^2;
M[[2, 2]] = Data[[2, 1]];
M[[2, 3]] = Data[[2, 1]]^2;
M[[k - 1, 3 k - 8]] = M[[k, 3 k - 8]] = 1;
M[[k - 1, 3 k - 7]] = Data[[k - 1, 1]];
M[[k - 1, 3 k - 6]] = Data[[k - 1, 1]]^2;
M[[k, 3 k - 7]] = Data[[k, 1]];
M[[k, 3 k - 6]] = Data[[k, 1]]^2;
For[i = 3, i <= k - 2, i = i + 1, M[[i, 3 (i - 2) + 1]] = 1;
M[[i, 3 (i - 2) + 2]] = Data[[i, 1]];
M[[i, 3 (i - 2) + 3]] = Data[[i, 1]]^2];
For[i = 1, i <= k - 3, i = i + 1, M[[k + i, 3 (i - 1) + 1]] = 1;
M[[k + i, 3 (i - 1) + 4]] = -1;
M[[k + i, 3 (i - 1) + 2]] = NewData[[i, 1]];
M[[k + i, 3 (i - 1) + 5]] = -NewData[[i, 1]];
M[[k + i, 3 (i - 1) + 3]] = NewData[[i, 1]]^2;
M[[k + i, 3 (i - 1) + 6]] = -1*NewData[[i, 1]]^2;
M[[k + i + k - 3, 3 (i - 1) + 2]] = 1;
M[[k + i + k - 3, 3 (i - 1) + 5]] = -1;
M[[k + i + k - 3, 3 (i - 1) + 3]] = 2*NewData[[i, 1]];
M[[k + i + k - 3, 3 (i - 1) + 6]] = -2*NewData[[i, 1]]];
Coef = Inverse[M].y;
pf = Table[{Coef[[3 (i - 1) + 1]] + Coef[[3 (i - 1) + 2]]*x + Coef[[3 (i - 1) + 3]]*x^2, NNewData[[i, 1]] <= x <= NNewData[[i + 1, 1]]}, {i, 1, k - 2}];
Piecewise[pf])
y = Spline2[Data]
M // MatrixForm
a = Plot[{y}, {x, -1, -0.2}, Epilog -> {PointSize[Large], Point[Data]}, AxesLabel -> {"x", "y2"}, ImageSize -> Medium, PlotRange -> All, PlotLegends -> {"y2"},  PlotLabel -> "Quadratic spline"]
View Python Code
import numpy as np
import matplotlib.pyplot as plt

Data = [[-1, 0.038], [-0.8, 0.058], [-0.60, 0.10], [-0.4,0.20], [-0.2, 0.50]]
def Spline2(x, Data):
k = len(Data)
Data = np.array(Data)
NewData = np.array([(Data[i + 1] + Data[i + 2])/2 for i in range(k-3)])
NNewData = np.vstack([Data[0], NewData, Data[k-1]])
y = [Data[i][1] if i < k else 0 for i in range(3*k - 6)]
M = np.zeros([3*k - 6, 3*k - 6])
M[0][0] = M[1][0] = 1
M[0][1] = Data[0][0]
M[0][2] = Data[0][0]**2
M[1][1] = Data[1][0]
M[1][2] = Data[1][0]**2
M[k - 2][3*k - 9] = M[k - 1][3*k - 9] = 1
M[k - 2][3*k - 8] = Data[k - 2][0]
M[k - 2][3*k - 7] = Data[k - 2][0]**2
M[k - 1][3*k - 8] = Data[k - 1][0]
M[k - 1][3*k - 7] = Data[k - 1][0]**2
for i in range(2, k-2):
M[i][3*(i - 1)] = 1
M[i][3*(i - 1) + 1] = Data[i][0]
M[i][3*(i - 1) + 2] = Data[i][0]**2
for i in range(k - 3):
M[k + i][3*i] = 1
M[k + i][3*i + 3] = -1
M[k + i][3*i + 1] = NewData[i][0]
M[k + i][3*i + 4] = -NewData[i][0]
M[k + i][3*i + 2] = NewData[i][0]**2
M[k + i][3*i + 5] = -1*NewData[i][0]**2
M[k + i + k - 3][3*i + 1] = 1
M[k + i + k - 3][3*i + 4] = -1
M[k + i + k - 3][3*i + 2] = 2*NewData[i][0]
M[k + i + k - 3][3*i + 5] = -2*NewData[i][0]
print("M Matrix\n", M)
Coef = np.matmul(np.linalg.inv(M),y)
display([['{:30}'.format(str(round(Coef[3*i],3)) + ' + ' + str(round(Coef[3*i + 1],3)) + 'x + ' + str(round(Coef[3*i + 2],3)) + 'x**2'),
'{:15}'.format(str(round(NNewData[i][0],3))+' <=x<= '+str(round(NNewData[i + 1][0],3)))] for i in range(k - 2)])
return np.piecewise(x, [(x >= NNewData[i][0])&(x <= NNewData[i + 1][0]) for i in range(k - 2)],
[lambda x, j=i: Coef[3*j] + Coef[3*j + 1]*x + Coef[3*j + 2]*x**2 for i in range(k - 2)])
x = np.arange(-1,-0.2,0.01)
y = Spline2(x, Data)
plt.plot(x, y, label="y2")
plt.scatter([point[0] for point in Data],[point[1] for point in Data], c='k')