SHOW:
|
|
- or go back to the newest paste.
| 1 | #!/usr/bin/env python | |
| 2 | ||
| 3 | import numpy as np | |
| 4 | from matplotlib import pyplot as plt | |
| 5 | ||
| 6 | data_fname = 'knot_points.csv' | |
| 7 | # x1,y1 | |
| 8 | # x2,y2 | |
| 9 | # ... | |
| 10 | ||
| 11 | def read_data(fname): | |
| 12 | X = [] | |
| 13 | Y = [] | |
| 14 | with open(fname, 'r') as f: | |
| 15 | for line in f: | |
| 16 | (x, y) = line.split(',')
| |
| 17 | X.append(float(x)) | |
| 18 | Y.append(float(y)) | |
| 19 | return (X, Y) | |
| 20 | ||
| 21 | def langrange_polynomial(X, Y): | |
| 22 | def L(i): | |
| 23 | return lambda x: np.prod([(x-X[j])/(X[i]-X[j]) for j in range(len(X)) if i != j]) * Y[i] | |
| 24 | Sx = [L(i) for i in range(len(X))] # summands | |
| 25 | return lambda x: np.sum([s(x) for s in Sx]) | |
| 26 | ||
| 27 | # F = langrange_polynomial(*read_data(data_fname)) | |
| 28 | (X, Y) = read_data(data_fname) | |
| 29 | F = langrange_polynomial(X, Y) | |
| 30 | ||
| 31 | x_range = np.linspace(X[0], X[-1], 100) | |
| 32 | plt.plot(X, Y, 'ro') | |
| 33 | plt.plot(x_range, map(F, x_range)) | |
| 34 | plt.xlabel(r'$x$') | |
| 35 | plt.ylabel(r'$F(x)$') | |
| 36 | plt.title('Lagrange polynomial interpolation')
| |
| 37 | plt.grid(True) | |
| 38 | plt.show() |