i tried now for hours to find a function polynomial_interpolation that interpolates π data points with a polynomial of degree πβ1 by solving the linear system for the polynomial's coefficients. The function takes the interpolation points xs,ys and returns the interpolated values y corresponding to x. To solve the linear system use the Gaussian algorithm.

```
xs= np.array([-3,2,-1,0,1,2,3])
n = len(xs)
x =np.zeros(n)
ys= np.array([6, -2, 0, 1, -4, 2, 0])
def polynomial_interpolation(x, xs, ys):
"""
Interpolate a 1-D function using naive polynomial interpolation.
x : a float or an 1d-array
xs : (N,) array_like
A 1-D array of values to interpolate.
ys : (N,) array_like
A 1-D array of real values. The length of y along the
interpolation axis must be equal to the length of x.
returns : (N,) array_like
A 1-D array of the interpolated values at x.
"""
if len(ys) != len(xs):
raise ValueError("xs and ys are not of the same length.")
am = np.zeros((n,n))
forig = np.zeros(n)
forig = ys.copy()
for i in range(0, n):
for j in rnage(0,n):
am[i,j] = x[i]**j
return am
for i in range(0,n):
for j in range(0,n):
ys[j] = ys[j] - ys[i]* am[j,i]/ am[i,i]
am[j] = am[j] -am[i]* am[j,i]/ am[i,i]
print('obere rechte Dreiecksform')
for i in range(n-1,-1,-1):
a[i] = ys[i]
for i in range(i+1, n):
x[i] = x[i] - am[i,j] * x[j]
x[i] = x[i] - am[i,j]
return x
xplot = np.arange(xs[0]-1.1, xs[n-1]+1.1, (x[n-1]-x[0])/100)
ysplot = 0
for i in range(0,n):
ysplot = ysplot+x[i]*xplot**i
pyplot.plot(xplot, ysplot, "r", linewidth = 2)
pyplot.plot(x, forig, "o")
pyplot.show()
```

Read more here: https://stackoverflow.com/questions/65078528/python-linear-system-for-the-coefficients-of-a-polynomial-polynomial-interpol

### Content Attribution

This content was originally published by Chaox_stack at Recent Questions - Stack Overflow, and is syndicated here via their RSS feed. You can read the original post over there.