** Interpolate using scipy's interpolate lagrange method. ** **
** Interpolate a given N + 1 dataset ($ x_i $, $ y_i $) (i = 0,2, 3, ..., N) by Nth order polynomial ([Lagrange interpolation](https:: //ja.wikipedia.org/wiki/%E3%83%A9%E3%82%B0%E3%83%A9%E3%83%B3%E3%82%B8%E3%83%A5%E8%A3 % 9C% E9% 96% 93)) Implement the program in Python3. ** **
Consider $ y = 1 / (1 + x ^ 2) $ as an example. Sample an 11-point dataset ($ x_i $, $ y_i $) and interpolate it with a 10th-order polynomial. This function does not work with Lagrange polynomial [Runge phenomenon](https://ja.wikipedia.org/wiki/%E3%83%AB%E3%83%B3%E3%82%B2%E7%8F%BE%E8 It is known as an example of causing% B1% A1).
Code 1: Use scipy / numpy. Minimal code (** Click here if you are in a hurry **) Code 2: Calculate the coefficients for interpolation in the code. Click here to study the method.
code(1)scipyを利用してラクをするcode。
from scipy.interpolate import lagrange
import numpy as np
import matplotlib.pyplot as plt
##Main
x =np.linspace(-5,5,num=11) #[-5,5]Divide the range of into 11 equal parts and store in x
y = 1.0/(1.0+x**2) #The function considered in this example, y= 1/(1+x^2)Define
f_Lag=lagrange(x,y) #scipy.interpolate.Lagrange interpolation execution by lagrange
##
#for plot
xnew =np.linspace(-5,5,num=51) # [-5,5]Divide the range of into 51 equal parts and store in xnew
plt.plot(x, y, 'o', xnew, f_Lag(xnew), '-') #Raw data"o"So, Lagrange interpolated line('-')Draw with.
plt.legend(['Raw data','Lagrange'], loc='best') #specification of legend
plt.xlim([-6, 6]) #x-axis plot range
plt.ylim([0, 1.4]) #y-axis plot range
plt.show()
11 data points sampled by blue. The orange line is Lagrange polynomial.
code(2)code中に直接数値計算を実行する
"""
interpolation: ラグランジュinterpolation
example: y = 1/(1+x**2):Section-11 points from 5 to 5 are sampled and Lagrange interpolated.
"""
from math import pi,e, log, factorial
import matplotlib.pyplot as plt
###
def g(i, x): #Calculation of coefficients of Lagrange interpolation. Main.
dum=1.0
for j in range(len(x_lis)):
if j != i :
dum *= (x-x_lis[j])/(x_lis[i]-x_lis[j])
return dum
#
def fLag(x,m): #Lagrange interpolation
dum=0.0
for j in range(m):
dum += y_lis[j]*g(j, x)
return dum
###
##Dataset construction for the example
m = 11 # x= -11-point sampling at regular intervals from 5 to 5
x_lis = []
y_lis = []
def yy(x):
return 1/(1+x**2) #Example function y= 1/(1+x**2)
for k in range(m):
xm = -5.0 + 10.0*k / m
x_lis.append(xm)
y_lis.append(yy(xm))
plt.plot(x_lis,y_lis, 'o',label='Row data')
##
###Lagrange interpolation execution
mm = 5000
y_Laglis = []
xx_lis = []
for k in range(mm):
xm = -5.0 + 10.0*k / mm
xx_lis.append(xm)
y_lis_exact=[]
for j in range(mm):
y_Laglis.append(fLag(xx_lis[j],m))
y_lis_exact.append(yy(xx_lis[j]))
#plot
plt.grid(True)
plt.xlabel('x',fontsize=24)
plt.ylabel('f(x)',fontsize=24)
plt.xticks(fontsize=18)
plt.yticks(fontsize=18)
plt.plot(xx_lis,y_Laglis, color='Red',label='Lagrange')
plt.plot(xx_lis,y_lis_exact, color='Black',label='Exact')
plt.legend(loc='upper left')
plt.show()¥
11 data points sampled by blue. The red line is Lagrange polynomial. The black line represents the exact solution.
** As the order (data points) is increased, the approximation of the central part becomes better, but the approximation of the part near both ends becomes worse. It is known as an example in which the error increases as the order is increased and finally diverges to infinity. This is because equinox equinoxes are used, and it is known that it can be avoided if the coordination of the equinoxes is rough near the center and finely divided at both ends [1]. ** **
** For interval [-1,1]: When the domain of $ f (x) $ is analyzed and connected to the complex plane, its singularity $ z $ is **
** Does not meet **.
In this example, the singularity z of $ f (z) = 1 / (1 + z ^ 2) $ is $ ± i $.
The right side of the inequality is $ 4e ^ \ pi = 92.562770531 ...
If f (x) is a continuous function in the interval [a, b], the interpolation polynomial approaches f (x) at the limit where N is infinite if the equinox $ x_n $ is selected as follows [2]. ..
In addition, polynomial interpolation that selects the zero point of an orthogonal polynomial (for example, Chebisif polynomial) as the equinox does not cause the Rung phenomenon as seen above [1].
[1] Masatake Mori, ["Numerical Analysis Second Edition"](https://www.amazon.co.jp/%E6%95%B0%E5%80%A4%E8%A7%A3%E6%9E % 90-% E5% 85% B1% E7% AB% 8B% E6% 95% B0% E5% AD% A6% E8% AC% 9B% E5% BA% A7-% E6% A3% AE-% E6% AD% A3% E6% AD% A6 / dp / 4320017013), Kyoritsu Shuppan, 2002. [2] Masao Iri and Kazutake Fujino, ["Common knowledge of numerical calculation"](https://www.amazon.co.jp/%E6%95%B0%E5%80%A4%E8%A8%88 % E7% AE% 97% E3% 81% AE% E5% B8% B8% E8% AD% 98-% E4% BC% 8A% E7% 90% 86-% E6% AD% A3% E5% A4% AB / dp / 4320013433), Kyoritsu Shuppan, 1985.