Skip to content

Instantly share code, notes, and snippets.

@ti-nspire
Last active Feb 26, 2018
Embed
What would you like to do?
チェビシェフ補間で多項式近似式を求める
# 参考: 『インターフェース』 CQ 出版, 2017 年 9 月号, p.109-111, 三上直樹
# ある函数 f(x) が下のような N 次の多項式 g(x) で表現できると假定し、そこそこ精度のある近似式を求める。
# f(x) ≈ g(x)
# g(x) ≈ c[0] * T[0] + c[1] * T[1] + ... + c[N] * T[N]
# T[N] はチェビシェフ多項式である。
# ここでは例として f(x) = Sin(Pi * x / 2) の近似式を具体的に計算してみる。区間は -1 ≦ x ≦ 1 とする。
# 区間が a ≦ x ≦ b である場合は次式 u(a, b) で変数を変換する。
# u = lambda a, b: (2 * x - (b + a)) / (b - a)
import sympy as sym
import numpy as np
def cheby(f, N):
x = sym.Symbol("x")
k = np.arange(N+1)
s = [2 / (N + 1), 2 * k + 1, np.pi / (2 * (N + 1)), (1 / (N + 1))]
p = s[1] * s[2]
# 第 1 種チェビシェフ多項式を求める。
T = [sym.chebyshevt_poly(n, x) for n in k]
# 数列 r を求める。
r = np.cos(p)
# 数列 r を使って係数 c[0] を求める。
fVals = f(r)
c = np.zeros(N+1)
c[0] = s[3] * np.sum(fVals)
# 係数 c[1] ~ c[N] を求める。
c[1:N+1] = [s[0] * np.sum(fVals * np.cos(p * n)) for n in np.arange(1, N+1)]
# 係数 c とチェビシェフ多項式 T とで近似式 g(x) を組み立てる。
g = np.sum([c[n] * T[n] for n in k])
# 0 次 ~ N 次の係数を取り出す。
return np.array([round(g.coeff(x, n), 8) for n in k])
########
# test #
########
if __name__ == "__main__":
f = lambda x : np.sin(np.pi * x / 2) # この函数の近似式をチェビシェフ補間による多項式で求める。
N = 5 # N 次まで近似してみる。
cApprx = cheby(f, 5)
print(cApprx)
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment