소개

scipy의 interpolate의 lagrange 메소드를 이용한 보간을 실시한다.

주어진 N + 1 개의 데이터 세트 ($ x_i $, $ y_i $) (i = 0,2, 3, ..., N)를 N 차 다항식으로 보간 (라그랑주 보간)하는 프로그램을 Python3에서 구현한다.

예를 들면 $y=1/(1+x^2)$를 생각한다. 11점의 데이터 세트($x_i$, $y_i$)를 샘플하고, 그것을 10차 다항식으로 보간한다. 이 함수는 라그랑주 보간이 잘되지 않는 Runge 현상을 일으키는 예로 알려져 있습니다.

내용

코드 1: scipy·numpy를 이용. 최소 코드 (급한 사람은 여기)
코드 2: 보간을 위한 계수를 코드에서 계산. 수법의 공부를 한다면 이쪽.

계산 코드 (1) : scipy · numpy 사용

코드 (1) scipy를 이용하여 편하게 하는 코드.

from scipy.interpolate import lagrange
import numpy as np
import matplotlib.pyplot as plt

##  メイン
x =np.linspace(-5,5,num=11) #[-5,5]の範囲を11等分して xに格納
y = 1.0/(1.0+x**2)          #本例題で考えている関数である y = 1/(1+x^2)を定義
f_Lag=lagrange(x,y) #scipy.interpolate.lagrangeによるラグランジュ補間実行
##

#for plot
xnew =np.linspace(-5,5,num=51) # [-5,5]の範囲を51等分して xnewに格納
plt.plot(x, y, 'o',  xnew, f_Lag(xnew), '-')   # 生データを"o"で， ラグランジュ補間したものを線('-')で描く。
plt.legend(['Raw data','Lagrange'], loc='best')  # legendの指定
plt.xlim([-6, 6])  # x軸のプロット範囲 
plt.ylim([0, 1.4]) # y軸のプロット範囲
plt.show()

결과(1)

청인이 샘플한 11개의 데이터점. 오렌지선이 라그랑주 보간한 것.

계산 코드 (2) 보간을위한 계수를 파이썬으로 계산

코드 (2) 코드 중에 직접 수치 계산을 실행합니다.

"""
補間: ラグランジュ補間
例題: y = 1/(1+x**2): を区間 -5から5まで11点サンプリングし，それをラグランジュ補間する。 
"""
from math import pi,e, log, factorial
import matplotlib.pyplot as plt


### 
def g(i, x): # ラグランジュ補間の係数の計算。メイン。
    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): #ラグランジュ補間
    dum=0.0
    for j in range(m):
        dum += y_lis[j]*g(j, x)

    return dum
###

## 例題のためのデータセット構築
m = 11  #  x= -5 から5まで等間隔に11点サンプリング
x_lis = []
y_lis = []
def yy(x):
    return 1/(1+x**2) # 例題の関数 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')
##

### ラグランジュ補間実行
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()¥

결과(2)

청인이 샘플한 11개의 데이터점. 적선이 라그랑주 보간한 것. 검은 선은 엄밀한 해를 나타냅니다.

차수(데이터점)를 늘려 가면 중앙 부분의 근사는 좋아지지만, 양단에 가까운 부분의 근사가 나빠진다. 차수를 올리면 오차가 커지고, 마침내 무한대로 발산하는 예로 알려져 있다. 이것은 등간격 분점을 이용했기 때문에, 분점의 배위를 중앙 부근에서 거칠게 양단에서 세세하게 하면 피할 수 있는 것이 알려져 있다[1].

부록 : 등간격 라그랑주 보간이 점수를 늘리면 원래의 함수 f (x)에 접근하기위한 조건

구간 [-1,1]의 경우: $f(x)$의 정의 영역을 복소 평면에 해석 접속했을 때, 그 특이점 $z$가

$|(1+z)^{1+z}|\leqq 4 e^{\pi |Im(z)|}$

를 채우지 않는 것이다.
본 예제에서 말하면, $f(z) = 1/(1+z^2)$의 특이점 z는 $±i$이고,
부등식의 오른쪽은 $4e^\pi= 92.562770531... $,
왼쪽은 $|(1+i)^{1+i}| = 0.644793...$
그러므로 위의 부등식을 만족시켜 버리기 때문에 등간격 분점을 이용한 라그랑주 보간(다항식 보간)은 분점수의 증화와 함께 오차가 증가한다.

f (x)가 구간 [a, b]에서 연속 함수이면 분점 $ x_n $를 다음과 같이 선택하면 N이 무한대의 극한이고 보간 다항식이 f (x)에 접근한다 [2] .

$x_n = (a+b)/2 + (b-a)(sin(-\pi/2+n\pi/N))/2,$
$(n = 0, 1, 2, ...,N) $

또, 직교 다항식(예를 들면 체비시프 다항식)의 영점을 분점으로 선택하는 다항식 보간에서는 위에서 본 것 같은 Rung 현상은 생기지 않는다[1].

참고문헌

[1]모리 마사타케, 「수치 해석 제2판」 , 공립 출판, 2002.
[2] 이리 마사오 · 후지노 카즈키, 「수치 계산의 상식」, 공립 출판, 1985.