何点かでしか値がわからない関数を、多項式で置き換えて近似するにはどうすればよいか、という話です。点を通る多項式は簡単に作れますが、点をどう並べるかで結果が大きく変わり、等間隔に並べると点を増やすほどかえって悪くなる場合があります。この記事を読み終えると、補間多項式の作り方と誤差の見積もり方、そして点の置き方の良し悪しを数で測る方法がわかります。
相異なる n + 1 n+1 n + 1 点での値を与えると、それを通る n n n 次以下の多項式(補間多項式 )がただ 1 つ存在します。証明は Vandermonde 行列式が 0 0 0 でないことです。
補間多項式は差分商 を係数とする Newton 形で書け、点を 1 つ足すたびに項を 1 つ足すだけで更新できます。差分商は f ( n ) ( ξ ) / n ! f^{(n)}(\xi)/n! f ( n ) ( ξ ) / n ! に等しく、導関数の離散版です。
補間誤差 は f ( x ) − p n ( x ) = f ( n + 1 ) ( ξ ) ( n + 1 ) ! ∏ j ( x − x j ) f(x) - p_n(x) = \dfrac{f^{(n+1)}(\xi)}{(n+1)!}\prod_j (x - x_j) f ( x ) − p n ( x ) = ( n + 1 )! f ( n + 1 ) ( ξ ) ∏ j ( x − x j ) です。証明は Rolle の定理の繰り返しです。値と導関数の両方を合わせる Hermite 補間も同じ形の誤差公式をもち、これが Gauss 求積の誤差評価の土台になります。
補間の安定性は Lebesgue 定数 Λ n \Lambda_n Λ n で測ります。∥ f − p n ∥ ≤ ( 1 + Λ n ) E n ( f ) \|f - p_n\| \le (1 + \Lambda_n) E_n(f) ∥ f − p n ∥ ≤ ( 1 + Λ n ) E n ( f ) で、Λ n \Lambda_n Λ n はデータの誤差が何倍に増幅されるかを表す条件数 でもあります。
等間隔点では Λ n ≥ 2 n / ( 4 n 2 ) \Lambda_n \ge 2^n/(4n^2) Λ n ≥ 2 n / ( 4 n 2 ) と指数的に増え、Runge の関数 1 / ( 1 + 25 x 2 ) 1/(1+25x^2) 1/ ( 1 + 25 x 2 ) では補間が発散します。Chebyshev 点 では Λ n ≤ ( 2 / π ) log ( 2 n + 1 ) + 1 + 4 / π \Lambda_n \le (2/\pi)\log(2n+1) + 1 + 4/\pi Λ n ≤ ( 2/ π ) log ( 2 n + 1 ) + 1 + 4/ π と対数的にしか増えず、解析関数なら誤差は幾何級数的に減ります。
この記事で使う事前知識
平均値の定理・Rolle の定理 : 微分可能な関数が両端で同じ値をとるなら、途中に微分が 0 0 0 になる点がある、という性質。誤差公式の証明で繰り返し使います。→ 平均値の定理とテイラーの定理
行列式 : 正方行列に決まる数で、0 0 0 でないことが連立 1 次方程式がただ 1 つの解をもつ条件になります。→ 行列式
余因子展開 : 行列式を 1 つの行(または列)に沿って小さな行列式の和に分解する方法。→ 余因子
条件数 : 入力の小さな誤差が出力で何倍に拡大されるかを表す数。→ 条件数
一様ノルム : 関数の絶対値の最大値で、関数の大きさを測るもの。→ 一様ノルム
関数 f f f の値がいくつかの点でしかわからないとき、あるいは f f f が複雑で何度も評価したくないとき、f f f を多項式で置き換えて計算します。数値積分(多項式を積分する)、数値微分(多項式を微分する)、常微分方程式の多段法(過去の値を多項式で外挿する)、非線形方程式の割線法(2 点を通る 1 次式の零点)は、どれもこの置き換えの上に立っています。
自然な発想は「点を増やせば近似はよくなるだろう」というものです。ところが 1901 年に C. Runge が示したように、[ − 1 , 1 ] [-1,1] [ − 1 , 1 ] の等間隔点 で 1 / ( 1 + 25 x 2 ) 1/(1+25x^2) 1/ ( 1 + 25 x 2 ) を補間すると、点を増やすほど区間の端の近くで誤差が爆発します。滑らかな(実解析的な)関数でも、点の配置を誤ると補間は収束しないのです。一方、区間の端に点を密に置く Chebyshev 点 を使うと、同じ関数の補間は幾何級数的に収束します。
この違いを説明するのが本章の主題です。誤差公式の中の積 ∏ j ( x − x j ) \prod_j(x - x_j) ∏ j ( x − x j ) の大きさと、Lebesgue 定数という「補間という操作の条件数」の 2 つを調べると、等間隔点と Chebyshev 点の差が定量的にわかります。
以下、P n \mathcal{P}_n P n で実係数の n n n 次以下の多項式全体を表します。
補題 2.1 (Vandermonde 行列式 ) ✓
n ≥ 1 n \ge 1 n ≥ 1 とし、x 0 , … , x n x_0, \ldots, x_n x 0 , … , x n を実数とする。( i , j ) (i, j) ( i , j ) 成分が x i j x_i^j x i j (0 ≤ i , j ≤ n 0 \le i, j \le n 0 ≤ i , j ≤ n )の行列 V V V について
det V = ∏ 0 ≤ i < j ≤ n ( x j − x i ) \det V = \prod_{0 \le i \lt j \le n} (x_j - x_i) det V = 0 ≤ i < j ≤ n ∏ ( x j − x i ) が成り立つ。とくに x 0 , … , x n x_0, \ldots, x_n x 0 , … , x n が相異なれば det V ≠ 0 \det V \ne 0 det V = 0 である。
証明(補題 2.1) n n n についての帰納法で示します。n = 1 n = 1 n = 1 では det V = x 1 − x 0 \det V = x_1 - x_0 det V = x 1 − x 0 で正しいです。n ≥ 2 n \ge 2 n ≥ 2 とし、n n n 個の点についての主張を仮定します。
V V V の第 j j j 列から第 j − 1 j-1 j − 1 列の x 0 x_0 x 0 倍を引く操作を、j = n , n − 1 , … , 1 j = n, n-1, \ldots, 1 j = n , n − 1 , … , 1 の順に行います(後ろから行うので、引く列はつねに元のままです)。この操作で行列式は変わりません(命題 4.3[行列式とその性質] の 3)。操作後の ( i , j ) (i, j) ( i , j ) 成分は、j ≥ 1 j \ge 1 j ≥ 1 で x i j − x 0 x i j − 1 = x i j − 1 ( x i − x 0 ) x_i^j - x_0 x_i^{j-1} = x_i^{j-1}(x_i - x_0) x i j − x 0 x i j − 1 = x i j − 1 ( x i − x 0 ) です。とくに第 0 0 0 行は ( 1 , 0 , … , 0 ) (1, 0, \ldots, 0) ( 1 , 0 , … , 0 ) になります。第 0 0 0 行について余因子展開(定理 5.2[行列式とその性質] )すると、det V \det V det V は第 i i i 行(i = 1 , … , n i = 1, \ldots, n i = 1 , … , n )が ( x i − x 0 ) ( 1 , x i , … , x i n − 1 ) (x_i - x_0)(1, x_i, \ldots, x_i^{n-1}) ( x i − x 0 ) ( 1 , x i , … , x i n − 1 ) の n n n 次の行列式に等しくなります。各行から x i − x 0 x_i - x_0 x i − x 0 をくくり出すと(命題 4.3[行列式とその性質] の 2 を、5 により行に適用)
det V = ∏ i = 1 n ( x i − x 0 ) ⋅ det ( x i j ) 1 ≤ i ≤ n , 0 ≤ j ≤ n − 1 \det V = \prod_{i=1}^{n}(x_i - x_0) \cdot \det\bigl(x_i^{j}\bigr)_{1 \le i \le n,\ 0 \le j \le n-1} det V = i = 1 ∏ n ( x i − x 0 ) ⋅ det ( x i j ) 1 ≤ i ≤ n , 0 ≤ j ≤ n − 1 です。右辺の行列式は x 1 , … , x n x_1, \ldots, x_n x 1 , … , x n の Vandermonde 行列式なので、帰納法の仮定から ∏ 1 ≤ i < j ≤ n ( x j − x i ) \prod_{1 \le i \lt j \le n}(x_j - x_i) ∏ 1 ≤ i < j ≤ n ( x j − x i ) に等しく、主張を得ます。
∎
定理 2.2 (補間多項式の存在と一意性 ) ✓
x 0 , … , x n x_0, \ldots, x_n x 0 , … , x n を相異なる実数、y 0 , … , y n y_0, \ldots, y_n y 0 , … , y n を実数とする。p ( x i ) = y i p(x_i) = y_i p ( x i ) = y i (i = 0 , … , n i = 0, \ldots, n i = 0 , … , n )を満たす p ∈ P n p \in \mathcal{P}_n p ∈ P n がただ 1 つ存在し、
p ( x ) = ∑ i = 0 n y i ℓ i ( x ) , ℓ i ( x ) = ∏ j ≠ i x − x j x i − x j p(x) = \sum_{i=0}^{n} y_i\,\ell_i(x), \qquad \ell_i(x) = \prod_{j \ne i}\frac{x - x_j}{x_i - x_j} p ( x ) = i = 0 ∑ n y i ℓ i ( x ) , ℓ i ( x ) = j = i ∏ x i − x j x − x j と書ける(Lagrange 形 )。
証明(定理 2.2) n = 0 n = 0 n = 0 なら定数 y 0 y_0 y 0 がただ 1 つの解です。n ≥ 1 n \ge 1 n ≥ 1 とします。p ( x ) = ∑ j = 0 n a j x j p(x) = \sum_{j=0}^n a_j x^j p ( x ) = ∑ j = 0 n a j x j とおくと、条件 p ( x i ) = y i p(x_i) = y_i p ( x i ) = y i は係数ベクトル a = ( a 0 , … , a n ) T a = (a_0, \ldots, a_n)^{\mathsf{T}} a = ( a 0 , … , a n ) T についての連立 1 次方程式 V a = y Va = y V a = y です(V V V は 補題 2.1 の行列)。det V ≠ 0 \det V \ne 0 det V = 0 なので、正則性の判定(定理 7.1[行列式とその性質] の 1 ⇒ 5)からただ 1 つの解があります。ℓ i \ell_i ℓ i は n n n 次の多項式で ℓ i ( x k ) = 1 \ell_i(x_k) = 1 ℓ i ( x k ) = 1 (k = i k = i k = i )、0 0 0 (k ≠ i k \ne i k = i )を満たすので、∑ i y i ℓ i \sum_i y_i\ell_i ∑ i y i ℓ i は条件を満たす P n \mathcal{P}_n P n の元であり、一意性からこれが解です。
∎
系 2.3 (零点の個数 ) ✓
q ∈ P n q \in \mathcal{P}_n q ∈ P n が相異なる n + 1 n + 1 n + 1 個の点で 0 0 0 になるならば、q q q は零多項式である。
証明(系 2.3) q q q と零多項式は、どちらもその n + 1 n+1 n + 1 点で値 0 0 0 をとる P n \mathcal{P}_n P n の元です。定理 2.2 の一意性から両者は等しくなります。
∎
以下、f f f の値 y i = f ( x i ) y_i = f(x_i) y i = f ( x i ) を補間する P n \mathcal{P}_n P n の元を p n p_n p n (点を明示したいときは p x 0 , … , x n p_{x_0, \ldots, x_n} p x 0 , … , x n )と書き、ω ( x ) = ∏ j = 0 n ( x − x j ) \omega(x) = \prod_{j=0}^{n}(x - x_j) ω ( x ) = ∏ j = 0 n ( x − x j ) とおきます。Lagrange 形は理論には便利ですが、点を 1 つ足すとすべての ℓ i \ell_i ℓ i が変わります。次の節の Newton 形はこの欠点を解消します。
定義 3.1 (差分商 ) ↗
相異なる点 x 0 , … , x k x_0, \ldots, x_k x 0 , … , x k と関数 f f f について、f f f をこの k + 1 k+1 k + 1 点で補間する P k \mathcal{P}_k P k の元の x k x^k x k の係数を f [ x 0 , … , x k ] f[x_0, \ldots, x_k] f [ x 0 , … , x k ] と書き、f f f の k k k 階の差分商 という。f [ x 0 ] = f ( x 0 ) f[x_0] = f(x_0) f [ x 0 ] = f ( x 0 ) である。
命題 3.2 (差分商の対称性と陽な表示 ) ✓
f [ x 0 , … , x k ] f[x_0, \ldots, x_k] f [ x 0 , … , x k ] は点の並べ方によらず、
f [ x 0 , … , x k ] = ∑ i = 0 k f ( x i ) ∏ j ≠ i ( x i − x j ) f[x_0, \ldots, x_k] = \sum_{i=0}^{k}\frac{f(x_i)}{\prod_{j \ne i}(x_i - x_j)} f [ x 0 , … , x k ] = i = 0 ∑ k ∏ j = i ( x i − x j ) f ( x i ) である。
証明(命題 3.2) 補間多項式は点の集合だけで決まる(定理 2.2 )ので、並べ方によりません。Lagrange 形の ℓ i \ell_i ℓ i の x k x^k x k の係数は 1 / ∏ j ≠ i ( x i − x j ) 1/\prod_{j \ne i}(x_i - x_j) 1/ ∏ j = i ( x i − x j ) なので、∑ i f ( x i ) ℓ i \sum_i f(x_i)\ell_i ∑ i f ( x i ) ℓ i の x k x^k x k の係数は右辺です。
∎
証明(定理 3.3) 1. p k = p x 0 , … , x k p_k = p_{x_0, \ldots, x_k} p k = p x 0 , … , x k とおき、p k = p k − 1 + f [ x 0 , … , x k ] ∏ j < k ( x − x j ) p_k = p_{k-1} + f[x_0, \ldots, x_k]\prod_{j \lt k}(x - x_j) p k = p k − 1 + f [ x 0 , … , x k ] ∏ j < k ( x − x j ) (k ≥ 1 k \ge 1 k ≥ 1 )を示せば、k k k について足し合わせて主張を得ます。r = p k − p k − 1 r = p_k - p_{k-1} r = p k − p k − 1 は P k \mathcal{P}_k P k の元で x 0 , … , x k − 1 x_0, \ldots, x_{k-1} x 0 , … , x k − 1 で 0 0 0 になり、p k − 1 p_{k-1} p k − 1 の次数は k − 1 k-1 k − 1 以下なので x k x^k x k の係数は f [ x 0 , … , x k ] f[x_0, \ldots, x_k] f [ x 0 , … , x k ] です。s = r − f [ x 0 , … , x k ] ∏ j < k ( x − x j ) s = r - f[x_0, \ldots, x_k]\prod_{j \lt k}(x - x_j) s = r − f [ x 0 , … , x k ] ∏ j < k ( x − x j ) は x k x^k x k の項が打ち消し合うので P k − 1 \mathcal{P}_{k-1} P k − 1 の元で、k k k 個の点 x 0 , … , x k − 1 x_0, \ldots, x_{k-1} x 0 , … , x k − 1 で 0 0 0 になります。系 2.3 から s = 0 s = 0 s = 0 です。
2. q 0 = p x 0 , … , x k − 1 q_0 = p_{x_0, \ldots, x_{k-1}} q 0 = p x 0 , … , x k − 1 、q 1 = p x 1 , … , x k q_1 = p_{x_1, \ldots, x_k} q 1 = p x 1 , … , x k (ともに P k − 1 \mathcal{P}_{k-1} P k − 1 の元)とし、
P ( x ) = ( x − x 0 ) q 1 ( x ) − ( x − x k ) q 0 ( x ) x k − x 0 P(x) = \frac{(x - x_0)\,q_1(x) - (x - x_k)\,q_0(x)}{x_k - x_0} P ( x ) = x k − x 0 ( x − x 0 ) q 1 ( x ) − ( x − x k ) q 0 ( x ) とおきます。P ∈ P k P \in \mathcal{P}_k P ∈ P k で、P ( x 0 ) = q 0 ( x 0 ) = f ( x 0 ) P(x_0) = q_0(x_0) = f(x_0) P ( x 0 ) = q 0 ( x 0 ) = f ( x 0 ) 、P ( x k ) = q 1 ( x k ) = f ( x k ) P(x_k) = q_1(x_k) = f(x_k) P ( x k ) = q 1 ( x k ) = f ( x k ) 、0 < i < k 0 \lt i \lt k 0 < i < k では P ( x i ) = ( ( x i − x 0 ) − ( x i − x k ) ) f ( x i ) / ( x k − x 0 ) = f ( x i ) P(x_i) = \bigl((x_i - x_0) - (x_i - x_k)\bigr)f(x_i)/(x_k - x_0) = f(x_i) P ( x i ) = ( ( x i − x 0 ) − ( x i − x k ) ) f ( x i ) / ( x k − x 0 ) = f ( x i ) です。一意性(定理 2.2 )から P = p k P = p_k P = p k で、両辺の x k x^k x k の係数を比べると漸化式を得ます。
∎
漸化式 2 から、差分商は表の形で次々に計算できます(f [ x i , x i + 1 ] f[x_i, x_{i+1}] f [ x i , x i + 1 ] の列、f [ x i , x i + 1 , x i + 2 ] f[x_i, x_{i+1}, x_{i+2}] f [ x i , x i + 1 , x i + 2 ] の列、…)。k = 1 , 2 k = 1, 2 k = 1 , 2 の場合は、非線形方程式の反復解法 で割線法の解析に使った差分商と一致します。そこで示した 2 階の場合(補題 6.1[非線形方程式の反復解法] )は、定理 3.5 の n = 2 n = 2 n = 2 の場合です。
補題 3.4 (Rolle の定理の繰り返し ) ✓
I I I を開区間、k ≥ 1 k \ge 1 k ≥ 1 とし、g : I → R g : I \to \mathbb{R} g : I → R は k k k 回微分可能で、相異なる k + 1 k + 1 k + 1 個の点 t 0 < t 1 < ⋯ < t k t_0 \lt t_1 \lt \cdots \lt t_k t 0 < t 1 < ⋯ < t k (I I I の点)で 0 0 0 になるとする。このとき g ( k ) ( ξ ) = 0 g^{(k)}(\xi) = 0 g ( k ) ( ξ ) = 0 を満たす ξ ∈ ( t 0 , t k ) \xi \in (t_0, t_k) ξ ∈ ( t 0 , t k ) が存在する。
証明(補題 3.4) k k k についての帰納法です。g g g は微分可能なので各 [ t i , t i + 1 ] [t_{i}, t_{i+1}] [ t i , t i + 1 ] で連続かつ内部で微分可能で、両端で 0 0 0 なので、Rolle の定理(定理 3.1[平均値の定理とテイラーの定理] )から g ′ ( s i ) = 0 g'(s_i) = 0 g ′ ( s i ) = 0 となる s i ∈ ( t i , t i + 1 ) s_i \in (t_i, t_{i+1}) s i ∈ ( t i , t i + 1 ) があります(i = 0 , … , k − 1 i = 0, \ldots, k-1 i = 0 , … , k − 1 )。s 0 < ⋯ < s k − 1 s_0 \lt \cdots \lt s_{k-1} s 0 < ⋯ < s k − 1 は相異なり ( t 0 , t k ) (t_0, t_k) ( t 0 , t k ) に入ります。k = 1 k = 1 k = 1 ならこれで終わりです。k ≥ 2 k \ge 2 k ≥ 2 なら g ′ g' g ′ は k − 1 k-1 k − 1 回微分可能で k k k 個の点で 0 0 0 なので、帰納法の仮定から ( g ′ ) ( k − 1 ) = g ( k ) (g')^{(k-1)} = g^{(k)} ( g ′ ) ( k − 1 ) = g ( k ) は ( s 0 , s k − 1 ) ⊂ ( t 0 , t k ) (s_0, s_{k-1}) \subset (t_0, t_k) ( s 0 , s k − 1 ) ⊂ ( t 0 , t k ) の点で 0 0 0 になります。
∎
定理 3.5 (差分商の平均値定理 ) ✓
n ≥ 1 n \ge 1 n ≥ 1 、I I I を開区間、f : I → R f : I \to \mathbb{R} f : I → R を n n n 回微分可能とし、x 0 , … , x n ∈ I x_0, \ldots, x_n \in I x 0 , … , x n ∈ I を相異なる点とする。このとき min i x i < ξ < max i x i \min_i x_i \lt \xi \lt \max_i x_i min i x i < ξ < max i x i を満たす ξ \xi ξ があって
f [ x 0 , … , x n ] = f ( n ) ( ξ ) n ! f[x_0, \ldots, x_n] = \frac{f^{(n)}(\xi)}{n!} f [ x 0 , … , x n ] = n ! f ( n ) ( ξ ) が成り立つ。
証明(定理 3.5) g = f − p n g = f - p_n g = f − p n は I I I 上 n n n 回微分可能で、n + 1 n+1 n + 1 個の点 x 0 , … , x n x_0, \ldots, x_n x 0 , … , x n で 0 0 0 です。補題 3.4 から g ( n ) ( ξ ) = 0 g^{(n)}(\xi) = 0 g ( n ) ( ξ ) = 0 となる ξ \xi ξ が最小の点と最大の点の間にあります。p n p_n p n の x n x^n x n の係数は f [ x 0 , … , x n ] f[x_0, \ldots, x_n] f [ x 0 , … , x n ] なので p n ( n ) = n ! f [ x 0 , … , x n ] p_n^{(n)} = n!\,f[x_0, \ldots, x_n] p n ( n ) = n ! f [ x 0 , … , x n ] で、g ( n ) ( ξ ) = 0 g^{(n)}(\xi) = 0 g ( n ) ( ξ ) = 0 が主張です。
∎
定理 4.1 (補間誤差の公式 ) ✓
I I I を開区間、f : I → R f : I \to \mathbb{R} f : I → R を n + 1 n+1 n + 1 回微分可能とし、x 0 , … , x n ∈ I x_0, \ldots, x_n \in I x 0 , … , x n ∈ I を相異なる点、p n ∈ P n p_n \in \mathcal{P}_n p n ∈ P n をそこでの f f f の補間多項式とする。任意の x ∈ I x \in I x ∈ I について、x , x 0 , … , x n x, x_0, \ldots, x_n x , x 0 , … , x n の最小値と最大値の間の ξ \xi ξ があって
f ( x ) − p n ( x ) = f ( n + 1 ) ( ξ ) ( n + 1 ) ! ω ( x ) , ω ( x ) = ∏ j = 0 n ( x − x j ) f(x) - p_n(x) = \frac{f^{(n+1)}(\xi)}{(n+1)!}\,\omega(x), \qquad \omega(x) = \prod_{j=0}^{n}(x - x_j) f ( x ) − p n ( x ) = ( n + 1 )! f ( n + 1 ) ( ξ ) ω ( x ) , ω ( x ) = j = 0 ∏ n ( x − x j ) が成り立つ。
証明の方針は、誤差を ω ( x ) \omega(x) ω ( x ) の定数倍で打ち消す補助関数 g g g を作り、g g g の零点が多いことから Rolle の定理(補題 3.4 )で g ( n + 1 ) g^{(n+1)} g ( n + 1 ) の零点を取り出す、というものです。なお「モニック」は最高次の係数が 1 1 1 の多項式のことで、ω \omega ω はその例です。
証明(定理 4.1) x x x がどれかの x j x_j x j に等しければ両辺は 0 0 0 です。そうでないとき ω ( x ) ≠ 0 \omega(x) \ne 0 ω ( x ) = 0 なので、c = ( f ( x ) − p n ( x ) ) / ω ( x ) c = (f(x) - p_n(x))/\omega(x) c = ( f ( x ) − p n ( x )) / ω ( x ) とおき、g ( t ) = f ( t ) − p n ( t ) − c ω ( t ) g(t) = f(t) - p_n(t) - c\,\omega(t) g ( t ) = f ( t ) − p n ( t ) − c ω ( t ) を考えます。g g g は n + 1 n+1 n + 1 回微分可能で、n + 2 n+2 n + 2 個の相異なる点 x 0 , … , x n , x x_0, \ldots, x_n, x x 0 , … , x n , x で 0 0 0 です。補題 3.4 (k = n + 1 k = n+1 k = n + 1 )から、それらの最小値と最大値の間の ξ \xi ξ で g ( n + 1 ) ( ξ ) = 0 g^{(n+1)}(\xi) = 0 g ( n + 1 ) ( ξ ) = 0 です。p n ( n + 1 ) = 0 p_n^{(n+1)} = 0 p n ( n + 1 ) = 0 、ω \omega ω はモニックな n + 1 n+1 n + 1 次式なので ω ( n + 1 ) = ( n + 1 ) ! \omega^{(n+1)} = (n+1)! ω ( n + 1 ) = ( n + 1 )! で、0 = f ( n + 1 ) ( ξ ) − c ( n + 1 ) ! 0 = f^{(n+1)}(\xi) - c(n+1)! 0 = f ( n + 1 ) ( ξ ) − c ( n + 1 )! です。c c c の定義に戻せば主張を得ます。
∎
誤差は 2 つの因子の積です。f ( n + 1 ) / ( n + 1 ) ! f^{(n+1)}/(n+1)! f ( n + 1 ) / ( n + 1 )! は関数の性質で、こちらは選べません。ω ( x ) \omega(x) ω ( x ) は点の配置だけで決まり、こちらは選べます。区間の上で max ∣ ω ∣ \max|\omega| max ∣ ω ∣ を最小にする点の選び方が、後の節の Chebyshev 点です(系 6.3 )。
定理 4.2 (Hermite 補間の存在・一意性と誤差 ) ✓
x 0 , … , x n x_0, \ldots, x_n x 0 , … , x n を相異なる実数とする。
任意の実数 y i , z i y_i, z_i y i , z i (i = 0 , … , n i = 0, \ldots, n i = 0 , … , n )について、H ( x i ) = y i H(x_i) = y_i H ( x i ) = y i 、H ′ ( x i ) = z i H'(x_i) = z_i H ′ ( x i ) = z i を満たす H ∈ P 2 n + 1 H \in \mathcal{P}_{2n+1} H ∈ P 2 n + 1 がただ 1 つ存在する。
I I I を x 0 , … , x n x_0, \ldots, x_n x 0 , … , x n を含む開区間、f : I → R f : I \to \mathbb{R} f : I → R を 2 n + 2 2n+2 2 n + 2 回微分可能とし、H H H を y i = f ( x i ) y_i = f(x_i) y i = f ( x i ) 、z i = f ′ ( x i ) z_i = f'(x_i) z i = f ′ ( x i ) についての 1 の多項式とする。任意の x ∈ I x \in I x ∈ I について、x , x 0 , … , x n x, x_0, \ldots, x_n x , x 0 , … , x n の最小値と最大値の間の ξ \xi ξ があって
f ( x ) − H ( x ) = f ( 2 n + 2 ) ( ξ ) ( 2 n + 2 ) ! ω ( x ) 2 f(x) - H(x) = \frac{f^{(2n+2)}(\xi)}{(2n+2)!}\,\omega(x)^2 f ( x ) − H ( x ) = ( 2 n + 2 )! f ( 2 n + 2 ) ( ξ ) ω ( x ) 2
が成り立つ。
証明の方針は、1 では未知数と条件の数がどちらも 2 n + 2 2n+2 2 n + 2 なので、すべて 0 0 0 のデータから H = 0 H = 0 H = 0 しか出ないことを示せば十分だ、というものです。2 は補間誤差の公式と同じ手順で、各 x i x_i x i が値と傾きの 2 回分の零点として数えられるため、零点が 2 n + 2 2n+2 2 n + 2 個になって ω ( x ) 2 \omega(x)^2 ω ( x ) 2 が現れます。
証明(定理 4.2) 1. H ( x ) = ∑ k = 0 2 n + 1 a k x k H(x) = \sum_{k=0}^{2n+1} a_k x^k H ( x ) = ∑ k = 0 2 n + 1 a k x k とおくと、2 n + 2 2n+2 2 n + 2 個の条件は係数 a 0 , … , a 2 n + 1 a_0, \ldots, a_{2n+1} a 0 , … , a 2 n + 1 についての 2 n + 2 2n+2 2 n + 2 元の連立 1 次方程式です。正則性の判定(定理 7.1[行列式とその性質] の 4 ⇒ 5)から、斉次方程式(y i = z i = 0 y_i = z_i = 0 y i = z i = 0 )の解が H = 0 H = 0 H = 0 だけであることを示せば十分です。H ∈ P 2 n + 1 H \in \mathcal{P}_{2n+1} H ∈ P 2 n + 1 が H ( x i ) = H ′ ( x i ) = 0 H(x_i) = H'(x_i) = 0 H ( x i ) = H ′ ( x i ) = 0 を満たすとします。点を小さい順に並べると、Rolle の定理(定理 3.1[平均値の定理とテイラーの定理] )から隣り合う 2 点の間に H ′ H' H ′ の零点が 1 つずつ、計 n n n 個あり、これらは x i x_i x i と異なります。H ′ ( x i ) = 0 H'(x_i) = 0 H ′ ( x i ) = 0 と合わせて、H ′ ∈ P 2 n H' \in \mathcal{P}_{2n} H ′ ∈ P 2 n は相異なる 2 n + 1 2n+1 2 n + 1 個の零点をもつので、系 2.3 から H ′ = 0 H' = 0 H ′ = 0 です。よって H H H は定数で、H ( x 0 ) = 0 H(x_0) = 0 H ( x 0 ) = 0 から H = 0 H = 0 H = 0 です。
2. x x x がどれかの x j x_j x j なら両辺は 0 0 0 です。そうでないとき c = ( f ( x ) − H ( x ) ) / ω ( x ) 2 c = (f(x) - H(x))/\omega(x)^2 c = ( f ( x ) − H ( x )) / ω ( x ) 2 とおき、g ( t ) = f ( t ) − H ( t ) − c ω ( t ) 2 g(t) = f(t) - H(t) - c\,\omega(t)^2 g ( t ) = f ( t ) − H ( t ) − c ω ( t ) 2 を考えます。g g g は n + 2 n + 2 n + 2 個の相異なる点 x 0 , … , x n , x x_0, \ldots, x_n, x x 0 , … , x n , x で 0 0 0 なので、Rolle の定理から、これらを小さい順に並べたときの隣り合う 2 点の間に g ′ g' g ′ の零点が 1 つずつ、計 n + 1 n+1 n + 1 個あります。さらに ω ( x i ) = 0 \omega(x_i) = 0 ω ( x i ) = 0 から ( ω 2 ) ′ ( x i ) = 2 ω ( x i ) ω ′ ( x i ) = 0 (\omega^2)'(x_i) = 2\omega(x_i)\omega'(x_i) = 0 ( ω 2 ) ′ ( x i ) = 2 ω ( x i ) ω ′ ( x i ) = 0 なので g ′ ( x i ) = f ′ ( x i ) − H ′ ( x i ) = 0 g'(x_i) = f'(x_i) - H'(x_i) = 0 g ′ ( x i ) = f ′ ( x i ) − H ′ ( x i ) = 0 です。x i x_i x i は並べた点そのものなので、点の間にある先の n + 1 n+1 n + 1 個とは異なります。よって g ′ g' g ′ は相異なる 2 n + 2 2n+2 2 n + 2 個の零点をもち、補題 3.4 (g ′ g' g ′ に k = 2 n + 1 k = 2n+1 k = 2 n + 1 で適用)から g ( 2 n + 2 ) ( ξ ) = 0 g^{(2n+2)}(\xi) = 0 g ( 2 n + 2 ) ( ξ ) = 0 となる ξ \xi ξ が、x , x 0 , … , x n x, x_0, \ldots, x_n x , x 0 , … , x n の最小値と最大値の間にあります。H ( 2 n + 2 ) = 0 H^{(2n+2)} = 0 H ( 2 n + 2 ) = 0 、ω 2 \omega^2 ω 2 はモニックな 2 n + 2 2n+2 2 n + 2 次式なので ( ω 2 ) ( 2 n + 2 ) = ( 2 n + 2 ) ! (\omega^2)^{(2n+2)} = (2n+2)! ( ω 2 ) ( 2 n + 2 ) = ( 2 n + 2 )! で、f ( 2 n + 2 ) ( ξ ) = c ( 2 n + 2 ) ! f^{(2n+2)}(\xi) = c\,(2n+2)! f ( 2 n + 2 ) ( ξ ) = c ( 2 n + 2 )! です。
∎
Hermite 補間は、n + 1 n+1 n + 1 点の値と傾きという 2 n + 2 2n+2 2 n + 2 個の情報から 2 n + 1 2n+1 2 n + 1 次式を決めます。誤差公式の ω ( x ) 2 \omega(x)^2 ω ( x ) 2 は区間の上で符号を変えないので、積分の平均値の定理と組み合わせると「f ( 2 n + 2 ) ( ξ ) / ( 2 n + 2 ) ! f^{(2n+2)}(\xi)/(2n+2)! f ( 2 n + 2 ) ( ξ ) / ( 2 n + 2 )! × ∫ ω 2 \int\omega^2 ∫ ω 2 」の形の誤差公式が得られます。これが Newton–Cotes 則と Gauss 求積 の章で Gauss 求積の誤差を導く道筋です。
例 4.3 (Runge の関数を等間隔点と Chebyshev 点で補間する )
f ( x ) = 1 / ( 1 + 25 x 2 ) f(x) = 1/(1+25x^2) f ( x ) = 1/ ( 1 + 25 x 2 ) を [ − 1 , 1 ] [-1,1] [ − 1 , 1 ] の n + 1 n+1 n + 1 個の等間隔点と、T n + 1 T_{n+1} T n + 1 の零点 cos ( ( 2 k + 1 ) π / ( 2 n + 2 ) ) \cos\bigl((2k+1)\pi/(2n+2)\bigr) cos ( ( 2 k + 1 ) π / ( 2 n + 2 ) ) (後の節で定義する Chebyshev 点)で補間し、4001 点の格子上の最大誤差を比べます。補間多項式の係数は求めず、SciPy の重心形の補間(注意 6.9 )で値を計算します。
# Python 3.12.10 / NumPy 2.5.2 / SciPy 1.18.1
from numpy.polynomial import chebyshev as C
from scipy.interpolate import BarycentricInterpolator as BI
f = lambda x : 1 / ( 1 + 25 * x ** 2 ) # Runge の関数
xx = np. linspace ( - 1 , 1 , 4001 ) # 誤差を測る点
return np. max ( np. abs ( BI ( nodes , g ( nodes ))( xx ) - g ( xx )))
print ( " n 等間隔点 Chebyshev点 chebinterpolate x^5再現 " )
for n in [ 5 , 10 , 20 , 40 ] :
eq = np. linspace ( - 1 , 1 , n + 1 )
ch = np. cos ( ( 2 * k + 1 ) * np.pi / ( 2 * n + 2 ) ) # T_{n+1} の零点
c = C. chebinterpolate ( f , n ) # 同じ点で別の算法
e2 = np. max ( np. abs ( C. chebval ( xx , c ) - f ( xx )))
ok = err ( ch , lambda x : x ** 5 ) < 1e-12 # 5 次式は誤差 0 のはず
print ( f " {n :3d } { err ( eq ) :.2e } { err ( ch ) :.2e } {e2 :.2e } "
n 等間隔点 Chebyshev点 chebinterpolate x^5再現
5 4.33e-01 5.56e-01 5.56e-01 True
10 1.92e+00 1.09e-01 1.09e-01 True
20 5.98e+01 1.53e-02 1.53e-02 True
40 1.05e+05 2.89e-04 2.89e-04 True
検算. (a) n ≥ 5 n \ge 5 n ≥ 5 の補間は 5 次式を正確に再現するはずで(定理 2.2 の一意性)、最後の列は誤差が 10 − 12 10^{-12} 1 0 − 12 未満であることを示しています。(b) Chebyshev 点の誤差は、重心形の計算と、NumPy の chebinterpolate(同じ点での補間を Chebyshev 級数の係数として求める別の算法)とで表示の桁まで一致します。(c) 等間隔点の誤差は n n n を 2 倍にするごとに大きくなり(Runge 現象)、Chebyshev 点の誤差は n n n を 10 増やすごとにほぼ一定の比で減ります。この比が ρ − 10 \rho^{-10} ρ − 10 (ρ = ( 1 + 26 ) / 5 \rho = (1+\sqrt{26})/5 ρ = ( 1 + 26 ) /5 )であることは、後の例(例 6.8 )で確かめます。
n = 5 n = 5 n = 5 では等間隔点のほうがわずかに良いことにも注意してください。Chebyshev 点が有利なのは点数を増やしたときの振る舞いで、どの n n n でも勝つわけではありません。等間隔点で何が起きているのかを、次の節で補間という操作そのものの性質として捉えます。
この節では、補間という操作が入力のずれをどれだけ増幅するかを 1 つの数で測ります。[ a , b ] [a, b] [ a , b ] 上の連続関数の空間 C [ a , b ] C[a,b] C [ a , b ] に一様ノルム ∥ g ∥ = max [ a , b ] ∣ g ∣ \|g\| = \max_{[a,b]}|g| ∥ g ∥ = max [ a , b ] ∣ g ∣ (関数の絶対値の最大値)を入れます。[ a , b ] [a,b] [ a , b ] の相異なる点 x 0 , … , x n x_0, \ldots, x_n x 0 , … , x n を固定し、f ∈ C [ a , b ] f \in C[a,b] f ∈ C [ a , b ] にその補間多項式を対応させる写像を I n f = ∑ i f ( x i ) ℓ i \mathcal{I}_n f = \sum_i f(x_i)\ell_i I n f = ∑ i f ( x i ) ℓ i と書きます。I n \mathcal{I}_n I n は線形で、定理 2.2 の一意性から q ∈ P n q \in \mathcal{P}_n q ∈ P n について I n q = q \mathcal{I}_n q = q I n q = q です。
定義 5.1 (Lebesgue 関数と Lebesgue 定数 ) ↗
λ n ( x ) = ∑ i = 0 n ∣ ℓ i ( x ) ∣ \lambda_n(x) = \sum_{i=0}^{n}|\ell_i(x)| λ n ( x ) = ∑ i = 0 n ∣ ℓ i ( x ) ∣ を点 x 0 , … , x n x_0, \ldots, x_n x 0 , … , x n の Lebesgue 関数 、Λ n = max x ∈ [ a , b ] λ n ( x ) \Lambda_n = \max_{x \in [a,b]}\lambda_n(x) Λ n = max x ∈ [ a , b ] λ n ( x ) を Lebesgue 定数 という。また f ∈ C [ a , b ] f \in C[a,b] f ∈ C [ a , b ] について E n ( f ) = inf q ∈ P n ∥ f − q ∥ E_n(f) = \inf_{q \in \mathcal{P}_n}\|f - q\| E n ( f ) = inf q ∈ P n ∥ f − q ∥ を、n n n 次以下の多項式による最良一様近似の誤差 という。
λ n \lambda_n λ n は連続なので、最大値は存在します(定理 2.3[平均値の定理とテイラーの定理] )。∑ i ℓ i = 1 \sum_i \ell_i = 1 ∑ i ℓ i = 1 (定数 1 1 1 の補間)なので λ n ≥ 1 \lambda_n \ge 1 λ n ≥ 1 、したがって Λ n ≥ 1 \Lambda_n \ge 1 Λ n ≥ 1 です。
次の定理の 1 の「作用素ノルム」は、∥ f ∥ = 1 \|f\| = 1 ∥ f ∥ = 1 の関数 f f f に補間を施したとき、結果の大きさが最大でどこまで大きくなるかを表す数のことです。1 はそれが Λ n \Lambda_n Λ n に等しいと述べていて、証明の後半は、その最大を実現する折れ線の関数を実際に作る部分です。
定理 5.2 (Lebesgue 定数による誤差の評価 ) ✓
[ a , b ] [a,b] [ a , b ] の相異なる点 x 0 , … , x n x_0, \ldots, x_n x 0 , … , x n と f ∈ C [ a , b ] f \in C[a,b] f ∈ C [ a , b ] について次が成り立つ。
∥ I n f ∥ ≤ Λ n ∥ f ∥ \|\mathcal{I}_n f\| \le \Lambda_n\|f\| ∥ I n f ∥ ≤ Λ n ∥ f ∥ であり、等号を満たす f ≠ 0 f \ne 0 f = 0 が存在する。すなわち Λ n \Lambda_n Λ n は I n \mathcal{I}_n I n の作用素ノルムである。
∥ f − I n f ∥ ≤ ( 1 + Λ n ) E n ( f ) \|f - \mathcal{I}_n f\| \le (1 + \Lambda_n)\,E_n(f) ∥ f − I n f ∥ ≤ ( 1 + Λ n ) E n ( f ) 。
E n ( f ) → 0 E_n(f) \to 0 E n ( f ) → 0 (n → ∞ n \to \infty n → ∞ )。
証明(定理 5.2) 1. ∣ I n f ( x ) ∣ ≤ ∑ i ∣ f ( x i ) ∣ ∣ ℓ i ( x ) ∣ ≤ ∥ f ∥ λ n ( x ) ≤ ∥ f ∥ Λ n |\mathcal{I}_n f(x)| \le \sum_i |f(x_i)||\ell_i(x)| \le \|f\|\lambda_n(x) \le \|f\|\Lambda_n ∣ I n f ( x ) ∣ ≤ ∑ i ∣ f ( x i ) ∣∣ ℓ i ( x ) ∣ ≤ ∥ f ∥ λ n ( x ) ≤ ∥ f ∥ Λ n です。等号のために、λ n ( x ^ ) = Λ n \lambda_n(\hat x) = \Lambda_n λ n ( x ^ ) = Λ n となる x ^ \hat x x ^ を取り、s i s_i s i を ℓ i ( x ^ ) \ell_i(\hat x) ℓ i ( x ^ ) の符号(0 0 0 なら 1 1 1 )とします。点を小さい順に並べ、f f f を各点で s i s_i s i をとる折れ線(両端の点の外では一定)とすると、f f f は連続で ∥ f ∥ = 1 \|f\| = 1 ∥ f ∥ = 1 であり、I n f ( x ^ ) = ∑ i s i ℓ i ( x ^ ) = λ n ( x ^ ) = Λ n \mathcal{I}_n f(\hat x) = \sum_i s_i\ell_i(\hat x) = \lambda_n(\hat x) = \Lambda_n I n f ( x ^ ) = ∑ i s i ℓ i ( x ^ ) = λ n ( x ^ ) = Λ n です。
2. 任意の q ∈ P n q \in \mathcal{P}_n q ∈ P n について I n q = q \mathcal{I}_n q = q I n q = q なので、f − I n f = ( f − q ) − I n ( f − q ) f - \mathcal{I}_n f = (f - q) - \mathcal{I}_n(f - q) f − I n f = ( f − q ) − I n ( f − q ) です。1 から ∥ f − I n f ∥ ≤ ( 1 + Λ n ) ∥ f − q ∥ \|f - \mathcal{I}_n f\| \le (1 + \Lambda_n)\|f - q\| ∥ f − I n f ∥ ≤ ( 1 + Λ n ) ∥ f − q ∥ で、q q q について下限をとれば主張を得ます。
3. ε > 0 \varepsilon \gt 0 ε > 0 について、Weierstrass の多項式近似定理(系 4.4[Arzelà–Ascoli と Stone–Weierstrass] )から ∥ f − P ∥ < ε \|f - P\| \lt \varepsilon ∥ f − P ∥ < ε となる多項式 P P P があります。P P P の次数を N N N とすると、n ≥ N n \ge N n ≥ N で P ∈ P n P \in \mathcal{P}_n P ∈ P n なので E n ( f ) ≤ ε E_n(f) \le \varepsilon E n ( f ) ≤ ε です。
∎
2 と 3 を合わせても、∥ f − I n f ∥ → 0 \|f - \mathcal{I}_n f\| \to 0 ∥ f − I n f ∥ → 0 とは言えません。E n ( f ) E_n(f) E n ( f ) が 0 0 0 に近づく速さより Λ n \Lambda_n Λ n が増える速さのほうが速ければ、補間は発散しえます。Runge の現象はまさにこの状況です。
例 5.3 (Lebesgue 定数は補間の条件数 )
データ f ( x i ) f(x_i) f ( x i ) に大きさ ε \varepsilon ε 以下の誤差(測定誤差や丸め誤差)δ i \delta_i δ i が乗ると、補間多項式は ∑ i δ i ℓ i \sum_i\delta_i\ell_i ∑ i δ i ℓ i だけずれます。その一様ノルムは ε Λ n \varepsilon\Lambda_n ε Λ n 以下で、定理 5.2 の 1 から、最悪の場合(δ i = ε s i \delta_i = \varepsilon s_i δ i = ε s i )にはちょうど ε Λ n \varepsilon\Lambda_n ε Λ n になります。つまり Λ n \Lambda_n Λ n は「データから補間多項式への写像」の(一様ノルムで測った絶対)条件数です。倍精度の丸め誤差 u ≈ 1.1 × 10 − 16 u \approx 1.1 \times 10^{-16} u ≈ 1.1 × 1 0 − 16 を想定すると、Λ n ≈ 10 10 \Lambda_n \approx 10^{10} Λ n ≈ 1 0 10 の点の配置では補間の値は最悪で 6 桁程度しか信用できません。後の表(例 6.6 )では等間隔点の Λ 30 \Lambda_{30} Λ 30 が約 6.6 × 10 6 6.6 \times 10^6 6.6 × 1 0 6 で、n n n を 1 増やすごとにおよそ 2 倍になっていくので、n = 40 n = 40 n = 40 では約 5 × 10 9 5 \times 10^9 5 × 1 0 9 とこの水準に迫ります。
命題 5.4 (等間隔点の Lebesgue 定数の下界 ) ✓
n ≥ 1 n \ge 1 n ≥ 1 、h = ( b − a ) / n h = (b-a)/n h = ( b − a ) / n 、x j = a + j h x_j = a + jh x j = a + j h (j = 0 , … , n j = 0, \ldots, n j = 0 , … , n )とする。このとき
Λ n ≥ λ n ( a + h 2 ) ≥ 2 n 4 n 2 \Lambda_n \ge \lambda_n\Bigl(a + \frac h2\Bigr) \ge \frac{2^n}{4n^2} Λ n ≥ λ n ( a + 2 h ) ≥ 4 n 2 2 n である。
証明の方針は、端の最初の小区間の中点 a + h / 2 a + h/2 a + h /2 だけに注目し、そこでの各 ∣ ℓ j ∣ |\ell_j| ∣ ℓ j ∣ を下から評価して足し合わせる、というものです。足し合わせには二項定理を使うので、式変形が続きますが、最後は ( 2 n n ) \binom{2n}{n} ( n 2 n ) の下からの評価に帰着します。
証明(命題 5.4) x = a + t h x = a + th x = a + t h とおくと ℓ j ( x ) = ∏ k ≠ j ( t − k ) / ( j − k ) \ell_j(x) = \prod_{k \ne j}(t - k)/(j - k) ℓ j ( x ) = ∏ k = j ( t − k ) / ( j − k ) です。t = 1 / 2 t = 1/2 t = 1/2 では
∏ k = 0 n ∣ 1 2 − k ∣ = 1 2 ∏ k = 1 n 2 k − 1 2 = ( 2 n − 1 ) ! ! 2 n + 1 , ∏ k ≠ j ∣ j − k ∣ = j ! ( n − j ) ! \prod_{k=0}^{n}\Bigl|\frac12 - k\Bigr| = \frac12\prod_{k=1}^{n}\frac{2k - 1}{2} = \frac{(2n-1)!!}{2^{n+1}}, \qquad \prod_{k \ne j}|j - k| = j!\,(n-j)! k = 0 ∏ n 2 1 − k = 2 1 k = 1 ∏ n 2 2 k − 1 = 2 n + 1 ( 2 n − 1 )!! , k = j ∏ ∣ j − k ∣ = j ! ( n − j )! なので(( 2 n − 1 ) ! ! = 1 ⋅ 3 ⋯ ( 2 n − 1 ) (2n-1)!! = 1 \cdot 3 \cdots (2n-1) ( 2 n − 1 )!! = 1 ⋅ 3 ⋯ ( 2 n − 1 ) )、
∣ ℓ j ( a + h 2 ) ∣ = ( 2 n − 1 ) ! ! 2 n + 1 ∣ 1 / 2 − j ∣ j ! ( n − j ) ! ≥ ( 2 n − 1 ) ! ! 2 n + 1 n j ! ( n − j ) ! \Bigl|\ell_j\Bigl(a + \frac h2\Bigr)\Bigr| = \frac{(2n-1)!!}{2^{n+1}\,|1/2 - j|\; j!\,(n-j)!} \ge \frac{(2n-1)!!}{2^{n+1}\, n\; j!\,(n-j)!} ℓ j ( a + 2 h ) = 2 n + 1 ∣1/2 − j ∣ j ! ( n − j )! ( 2 n − 1 )!! ≥ 2 n + 1 n j ! ( n − j )! ( 2 n − 1 )!! です(∣ 1 / 2 − j ∣ ≤ n |1/2 - j| \le n ∣1/2 − j ∣ ≤ n を使いました)。二項定理から ∑ j = 0 n 1 / ( j ! ( n − j ) ! ) = 2 n / n ! \sum_{j=0}^{n} 1/(j!(n-j)!) = 2^n/n! ∑ j = 0 n 1/ ( j ! ( n − j )!) = 2 n / n ! なので
λ n ( a + h 2 ) ≥ ( 2 n − 1 ) ! ! 2 n ⋅ n ! = ( 2 n ) ! 2 n + 1 n ( n ! ) 2 = 1 2 n + 1 n ( 2 n n ) \lambda_n\Bigl(a + \frac h2\Bigr) \ge \frac{(2n-1)!!}{2n \cdot n!} = \frac{(2n)!}{2^{n+1}\,n\,(n!)^2} = \frac{1}{2^{n+1}n}\binom{2n}{n} λ n ( a + 2 h ) ≥ 2 n ⋅ n ! ( 2 n − 1 )!! = 2 n + 1 n ( n ! ) 2 ( 2 n )! = 2 n + 1 n 1 ( n 2 n ) です(( 2 n − 1 ) ! ! = ( 2 n ) ! / ( 2 n n ! ) (2n-1)!! = (2n)!/(2^n n!) ( 2 n − 1 )!! = ( 2 n )! / ( 2 n n !) )。最後に ( 2 n n ) ≥ 4 n / ( 2 n ) \binom{2n}{n} \ge 4^n/(2n) ( n 2 n ) ≥ 4 n / ( 2 n ) を示します。( 2 n k + 1 ) / ( 2 n k ) = ( 2 n − k ) / ( k + 1 ) \binom{2n}{k+1}/\binom{2n}{k} = (2n-k)/(k+1) ( k + 1 2 n ) / ( k 2 n ) = ( 2 n − k ) / ( k + 1 ) は k ≤ n − 1 k \le n-1 k ≤ n − 1 で 1 1 1 以上、k ≥ n k \ge n k ≥ n で 1 1 1 未満なので、( 2 n n ) \binom{2n}{n} ( n 2 n ) は ( 2 n k ) \binom{2n}{k} ( k 2 n ) の最大です。二項定理から 4 n = ∑ k = 0 2 n ( 2 n k ) = 2 + ∑ k = 1 2 n − 1 ( 2 n k ) ≤ 2 + ( 2 n − 1 ) ( 2 n n ) ≤ 2 n ( 2 n n ) 4^n = \sum_{k=0}^{2n}\binom{2n}{k} = 2 + \sum_{k=1}^{2n-1}\binom{2n}{k} \le 2 + (2n-1)\binom{2n}{n} \le 2n\binom{2n}{n} 4 n = ∑ k = 0 2 n ( k 2 n ) = 2 + ∑ k = 1 2 n − 1 ( k 2 n ) ≤ 2 + ( 2 n − 1 ) ( n 2 n ) ≤ 2 n ( n 2 n ) (( 2 n n ) ≥ 2 \binom{2n}{n} \ge 2 ( n 2 n ) ≥ 2 を使いました)です。合わせて λ n ( a + h / 2 ) ≥ 4 n / ( 2 n + 1 n ⋅ 2 n ) = 2 n / ( 4 n 2 ) \lambda_n(a + h/2) \ge 4^n/(2^{n+1}n \cdot 2n) = 2^n/(4n^2) λ n ( a + h /2 ) ≥ 4 n / ( 2 n + 1 n ⋅ 2 n ) = 2 n / ( 4 n 2 ) です。
∎
下界は n = 10 n = 10 n = 10 で 2.56 2.56 2.56 、n = 30 n = 30 n = 30 で約 3.0 × 10 5 3.0 \times 10^5 3.0 × 1 0 5 と控えめですが、指数的に増えることが要点です。実際の値はもっと大きく(後の表で n = 30 n = 30 n = 30 で約 6.6 × 10 6 6.6 \times 10^6 6.6 × 1 0 6 )、およそ 2 n + 1 / ( e n log n ) 2^{n+1}/(e\,n\log n) 2 n + 1 / ( e n log n ) の速さで増えることが知られています(Trefethen『Approximation Theory and Approximation Practice』第 15 章)。
Chebyshev 多項式 T k T_k T k は、漸化式 T 0 = 1 T_0 = 1 T 0 = 1 、T 1 = x T_1 = x T 1 = x 、T k + 1 = 2 x T k − T k − 1 T_{k+1} = 2xT_k - T_{k-1} T k + 1 = 2 x T k − T k − 1 で定めます(定義 2.1[誤差保証つきの求積と関数近似] )。
この多項式の要点は、x = cos θ x = \cos\theta x = cos θ とおくと T k ( x ) = cos k θ T_k(x) = \cos k\theta T k ( x ) = cos k θ になることです。角度 θ \theta θ を等間隔に動かすと cos k θ \cos k\theta cos k θ は等しい振れ幅で振動し、[ − 1 , 1 ] [-1,1] [ − 1 , 1 ] 上で ∣ T k ∣ ≤ 1 |T_k| \le 1 ∣ T k ∣ ≤ 1 が保たれます。零点も角度が等間隔の点に対応するので、Chebyshev 点 は、半円周上に等間隔に並べた点を横軸へ落とした点になります。下の図は n = 5 n = 5 n = 5 (6 点)の場合です。
半円周上で角度が等間隔の 6 点(15°, 45°, ..., 165°)を横軸へ真下に落とした点が、Chebyshev 点 cos((2j+1)π/12)。区間の端ほど点が密になる。
補題 6.1 (Chebyshev 多項式の三角関数による表示 ) ✓
k ≥ 0 k \ge 0 k ≥ 0 について次が成り立つ。
すべての実数 θ \theta θ について T k ( cos θ ) = cos k θ T_k(\cos\theta) = \cos k\theta T k ( cos θ ) = cos k θ 。とくに [ − 1 , 1 ] [-1,1] [ − 1 , 1 ] 上で ∣ T k ∣ ≤ 1 |T_k| \le 1 ∣ T k ∣ ≤ 1 。
k ≥ 1 k \ge 1 k ≥ 1 なら T k T_k T k は k k k 次で、最高次の係数は 2 k − 1 2^{k-1} 2 k − 1 である。T k ( − x ) = ( − 1 ) k T k ( x ) T_k(-x) = (-1)^kT_k(x) T k ( − x ) = ( − 1 ) k T k ( x ) 。
k ≥ 1 k \ge 1 k ≥ 1 なら T k T_k T k の零点は cos ( ( 2 j + 1 ) π / ( 2 k ) ) \cos\bigl((2j+1)\pi/(2k)\bigr) cos ( ( 2 j + 1 ) π / ( 2 k ) ) (j = 0 , … , k − 1 j = 0, \ldots, k-1 j = 0 , … , k − 1 )の k k k 個で、すべて ( − 1 , 1 ) (-1, 1) ( − 1 , 1 ) の相異なる点である。また η j = cos ( j π / k ) \eta_j = \cos(j\pi/k) η j = cos ( j π / k ) (j = 0 , … , k j = 0, \ldots, k j = 0 , … , k )で T k ( η j ) = ( − 1 ) j T_k(\eta_j) = (-1)^j T k ( η j ) = ( − 1 ) j である。
証明(補題 6.1) 1. k = 0 , 1 k = 0, 1 k = 0 , 1 では定義どおりです。k − 1 , k k - 1, k k − 1 , k で成り立てば、加法定理 cos ( k + 1 ) θ + cos ( k − 1 ) θ = 2 cos θ cos k θ \cos(k+1)\theta + \cos(k-1)\theta = 2\cos\theta\cos k\theta cos ( k + 1 ) θ + cos ( k − 1 ) θ = 2 cos θ cos k θ から
T k + 1 ( cos θ ) = 2 cos θ T k ( cos θ ) − T k − 1 ( cos θ ) = 2 cos θ cos k θ − cos ( k − 1 ) θ = cos ( k + 1 ) θ T_{k+1}(\cos\theta) = 2\cos\theta\,T_k(\cos\theta) - T_{k-1}(\cos\theta) = 2\cos\theta\cos k\theta - \cos(k-1)\theta = \cos(k+1)\theta T k + 1 ( cos θ ) = 2 cos θ T k ( cos θ ) − T k − 1 ( cos θ ) = 2 cos θ cos k θ − cos ( k − 1 ) θ = cos ( k + 1 ) θ です。x ∈ [ − 1 , 1 ] x \in [-1,1] x ∈ [ − 1 , 1 ] は x = cos θ x = \cos\theta x = cos θ と書けるので ∣ T k ( x ) ∣ ≤ 1 |T_k(x)| \le 1 ∣ T k ( x ) ∣ ≤ 1 です。
2. 漸化式から帰納法で、T k + 1 T_{k+1} T k + 1 の最高次の項は 2 x ⋅ 2 k − 1 x k = 2 k x k + 1 2x \cdot 2^{k-1}x^k = 2^k x^{k+1} 2 x ⋅ 2 k − 1 x k = 2 k x k + 1 です(k = 0 k = 0 k = 0 から 1 1 1 へは T 1 = x T_1 = x T 1 = x )。偶奇も漸化式から帰納法で従います。
3. θ j = ( 2 j + 1 ) π / ( 2 k ) ∈ ( 0 , π ) \theta_j = (2j+1)\pi/(2k) \in (0, \pi) θ j = ( 2 j + 1 ) π / ( 2 k ) ∈ ( 0 , π ) は相異なり、cos \cos cos は [ 0 , π ] [0, \pi] [ 0 , π ] で狭義減少なので cos θ j \cos\theta_j cos θ j は ( − 1 , 1 ) (-1,1) ( − 1 , 1 ) の相異なる k k k 点です。1 から T k ( cos θ j ) = cos ( ( 2 j + 1 ) π / 2 ) = 0 T_k(\cos\theta_j) = \cos\bigl((2j+1)\pi/2\bigr) = 0 T k ( cos θ j ) = cos ( ( 2 j + 1 ) π /2 ) = 0 です。ほかにも零点があれば、T k ∈ P k T_k \in \mathcal{P}_k T k ∈ P k が k + 1 k+1 k + 1 個の相異なる零点をもつので 系 2.3 から T k = 0 T_k = 0 T k = 0 となり、最高次の係数が 2 k − 1 ≠ 0 2^{k-1} \ne 0 2 k − 1 = 0 であることに矛盾します。よってこれがすべての零点です。T k ( cos ( j π / k ) ) = cos j π = ( − 1 ) j T_k(\cos(j\pi/k)) = \cos j\pi = (-1)^j T k ( cos ( j π / k )) = cos j π = ( − 1 ) j も 1 から従います。
∎
定理 6.2 (Chebyshev 多項式の最小最大性 ) ✓
n ≥ 0 n \ge 0 n ≥ 0 とする。
[ − 1 , 1 ] [-1,1] [ − 1 , 1 ] 上で、n + 1 n+1 n + 1 次のモニック多項式 q q q は max [ − 1 , 1 ] ∣ q ∣ ≥ 2 − n \max_{[-1,1]}|q| \ge 2^{-n} max [ − 1 , 1 ] ∣ q ∣ ≥ 2 − n を満たし、q = 2 − n T n + 1 q = 2^{-n}T_{n+1} q = 2 − n T n + 1 で等号が成り立つ。
a < b a \lt b a < b とする。[ a , b ] [a,b] [ a , b ] 上で、n + 1 n+1 n + 1 次のモニック多項式 q q q は max [ a , b ] ∣ q ∣ ≥ 2 ( ( b − a ) / 4 ) n + 1 \max_{[a,b]}|q| \ge 2\bigl((b-a)/4\bigr)^{n+1} max [ a , b ] ∣ q ∣ ≥ 2 ( ( b − a ) /4 ) n + 1 を満たす。
証明(定理 6.2) 1. T ~ = 2 − n T n + 1 \tilde T = 2^{-n}T_{n+1} T ~ = 2 − n T n + 1 は 補題 6.1 の 2 からモニックで、1 から max ∣ T ~ ∣ ≤ 2 − n \max|\tilde T| \le 2^{-n} max ∣ T ~ ∣ ≤ 2 − n 、3 から η j = cos ( j π / ( n + 1 ) ) \eta_j = \cos(j\pi/(n+1)) η j = cos ( j π / ( n + 1 )) (j = 0 , … , n + 1 j = 0, \ldots, n+1 j = 0 , … , n + 1 )で T ~ ( η j ) = ( − 1 ) j 2 − n \tilde T(\eta_j) = (-1)^j2^{-n} T ~ ( η j ) = ( − 1 ) j 2 − n です。よって max ∣ T ~ ∣ = 2 − n \max|\tilde T| = 2^{-n} max ∣ T ~ ∣ = 2 − n です。
max ∣ q ∣ < 2 − n \max|q| \lt 2^{-n} max ∣ q ∣ < 2 − n を満たすモニックな q q q があったと仮定します。r = T ~ − q r = \tilde T - q r = T ~ − q は最高次の項が打ち消し合うので P n \mathcal{P}_n P n の元です。∣ q ( η j ) ∣ < 2 − n |q(\eta_j)| \lt 2^{-n} ∣ q ( η j ) ∣ < 2 − n なので r ( η j ) = ( − 1 ) j 2 − n − q ( η j ) r(\eta_j) = (-1)^j2^{-n} - q(\eta_j) r ( η j ) = ( − 1 ) j 2 − n − q ( η j ) の符号は ( − 1 ) j (-1)^j ( − 1 ) j です。η 0 > η 1 > ⋯ > η n + 1 \eta_0 \gt \eta_1 \gt \cdots \gt \eta_{n+1} η 0 > η 1 > ⋯ > η n + 1 の隣り合う 2 点で r r r の符号は逆なので、中間値の定理(定理 5.7[極限と連続性] )から各区間 ( η j + 1 , η j ) (\eta_{j+1}, \eta_j) ( η j + 1 , η j ) に r r r の零点があり、相異なる零点が n + 1 n+1 n + 1 個得られます。系 2.3 から r = 0 r = 0 r = 0 ですが、r ( η 0 ) > 0 r(\eta_0) \gt 0 r ( η 0 ) > 0 に矛盾します。
2. x = a + b 2 + b − a 2 t x = \frac{a+b}{2} + \frac{b-a}{2}t x = 2 a + b + 2 b − a t とおくと、Q ( t ) = ( 2 / ( b − a ) ) n + 1 q ( a + b 2 + b − a 2 t ) Q(t) = (2/(b-a))^{n+1}q\bigl(\frac{a+b}{2} + \frac{b-a}{2}t\bigr) Q ( t ) = ( 2/ ( b − a ) ) n + 1 q ( 2 a + b + 2 b − a t ) は t t t の n + 1 n+1 n + 1 次のモニック多項式です。1 から max t ∈ [ − 1 , 1 ] ∣ Q ( t ) ∣ ≥ 2 − n \max_{t \in [-1,1]}|Q(t)| \ge 2^{-n} max t ∈ [ − 1 , 1 ] ∣ Q ( t ) ∣ ≥ 2 − n なので、max [ a , b ] ∣ q ∣ ≥ ( ( b − a ) / 2 ) n + 1 2 − n = 2 ( ( b − a ) / 4 ) n + 1 \max_{[a,b]}|q| \ge ((b-a)/2)^{n+1}2^{-n} = 2((b-a)/4)^{n+1} max [ a , b ] ∣ q ∣ ≥ (( b − a ) /2 ) n + 1 2 − n = 2 (( b − a ) /4 ) n + 1 です。
∎
この定理を補間誤差の公式に組み合わせると、max ∣ ω ∣ \max|\omega| max ∣ ω ∣ を最小にする点が決まります。
系 6.3 (Chebyshev 点での補間誤差 ) ✓
a < b a \lt b a < b 、n ≥ 0 n \ge 0 n ≥ 0 とし、Chebyshev 点
x j = a + b 2 + b − a 2 cos ( 2 j + 1 ) π 2 n + 2 ( j = 0 , … , n ) x_j = \frac{a+b}{2} + \frac{b-a}{2}\cos\frac{(2j+1)\pi}{2n+2} \qquad (j = 0, \ldots, n) x j = 2 a + b + 2 b − a cos 2 n + 2 ( 2 j + 1 ) π ( j = 0 , … , n ) をとる。このとき max [ a , b ] ∣ ω ∣ = 2 ( ( b − a ) / 4 ) n + 1 \max_{[a,b]}|\omega| = 2((b-a)/4)^{n+1} max [ a , b ] ∣ ω ∣ = 2 (( b − a ) /4 ) n + 1 で、これは相異なる n + 1 n+1 n + 1 点のすべての選び方の中で最小である。f f f が [ a , b ] [a,b] [ a , b ] を含む開区間で n + 1 n+1 n + 1 回微分可能で ∣ f ( n + 1 ) ∣ ≤ M n + 1 |f^{(n+1)}| \le M_{n+1} ∣ f ( n + 1 ) ∣ ≤ M n + 1 ([ a , b ] [a,b] [ a , b ] 上)ならば
max [ a , b ] ∣ f − p n ∣ ≤ 2 M n + 1 ( n + 1 ) ! ( b − a 4 ) n + 1 \max_{[a,b]}|f - p_n| \le \frac{2M_{n+1}}{(n+1)!}\Bigl(\frac{b-a}{4}\Bigr)^{n+1} [ a , b ] max ∣ f − p n ∣ ≤ ( n + 1 )! 2 M n + 1 ( 4 b − a ) n + 1 である。とくに [ a , b ] = [ − 1 , 1 ] [a,b] = [-1,1] [ a , b ] = [ − 1 , 1 ] では右辺は M n + 1 / ( 2 n ( n + 1 ) ! ) M_{n+1}/(2^n(n+1)!) M n + 1 / ( 2 n ( n + 1 )!) である。
証明(系 6.3) t ( x ) = ( 2 x − a − b ) / ( b − a ) t(x) = (2x - a - b)/(b-a) t ( x ) = ( 2 x − a − b ) / ( b − a ) とおきます。補題 6.1 の 2・3 から、x x x の多項式 T n + 1 ( t ( x ) ) T_{n+1}(t(x)) T n + 1 ( t ( x )) は n + 1 n+1 n + 1 次で、最高次の係数は 2 n ( 2 / ( b − a ) ) n + 1 2^n(2/(b-a))^{n+1} 2 n ( 2/ ( b − a ) ) n + 1 、零点はちょうど x 0 , … , x n x_0, \ldots, x_n x 0 , … , x n です。ω − T n + 1 ( t ( x ) ) / ( 2 n ( 2 / ( b − a ) ) n + 1 ) \omega - T_{n+1}(t(x))/\bigl(2^n(2/(b-a))^{n+1}\bigr) ω − T n + 1 ( t ( x )) / ( 2 n ( 2/ ( b − a ) ) n + 1 ) は P n \mathcal{P}_n P n の元で n + 1 n+1 n + 1 個の点 x j x_j x j で 0 0 0 なので、系 2.3 から ω ( x ) = 2 − n ( ( b − a ) / 2 ) n + 1 T n + 1 ( t ( x ) ) \omega(x) = 2^{-n}((b-a)/2)^{n+1}T_{n+1}(t(x)) ω ( x ) = 2 − n (( b − a ) /2 ) n + 1 T n + 1 ( t ( x )) です。t ( x ) t(x) t ( x ) は [ a , b ] [a,b] [ a , b ] を [ − 1 , 1 ] [-1,1] [ − 1 , 1 ] に写すので、補題 6.1 の 1・3 から max [ a , b ] ∣ ω ∣ = 2 − n ( ( b − a ) / 2 ) n + 1 = 2 ( ( b − a ) / 4 ) n + 1 \max_{[a,b]}|\omega| = 2^{-n}((b-a)/2)^{n+1} = 2((b-a)/4)^{n+1} max [ a , b ] ∣ ω ∣ = 2 − n (( b − a ) /2 ) n + 1 = 2 (( b − a ) /4 ) n + 1 です。他の点の選び方でも ω \omega ω はモニックな n + 1 n+1 n + 1 次式なので、定理 6.2 の 2 からこれ以上です。誤差の評価は、定理 4.1 の ξ \xi ξ が [ a , b ] [a,b] [ a , b ] に入ることから従います。
∎
補題 6.4 (区間の外での Chebyshev 多項式の増大 ) ✓
x ≥ 1 x \ge 1 x ≥ 1 、w = x + x 2 − 1 w = x + \sqrt{x^2 - 1} w = x + x 2 − 1 とすると、k ≥ 0 k \ge 0 k ≥ 0 について
T k ( x ) = w k + w − k 2 , w k 2 ≤ T k ( x ) ≤ w k T_k(x) = \frac{w^k + w^{-k}}{2}, \qquad \frac{w^k}{2} \le T_k(x) \le w^k T k ( x ) = 2 w k + w − k , 2 w k ≤ T k ( x ) ≤ w k である。また x ≤ − 1 x \le -1 x ≤ − 1 では ∣ T k ( x ) ∣ = T k ( ∣ x ∣ ) |T_k(x)| = T_k(|x|) ∣ T k ( x ) ∣ = T k ( ∣ x ∣ ) である。
証明(補題 6.4) w ≥ 1 w \ge 1 w ≥ 1 で w − 1 = x − x 2 − 1 w^{-1} = x - \sqrt{x^2 - 1} w − 1 = x − x 2 − 1 なので w + w − 1 = 2 x w + w^{-1} = 2x w + w − 1 = 2 x です。S k = ( w k + w − k ) / 2 S_k = (w^k + w^{-k})/2 S k = ( w k + w − k ) /2 とおくと S 0 = 1 S_0 = 1 S 0 = 1 、S 1 = x S_1 = x S 1 = x で、
2 x S k − S k − 1 = ( w + w − 1 ) ( w k + w − k ) − w k − 1 − w 1 − k 2 = w k + 1 + w − k − 1 2 = S k + 1 2xS_k - S_{k-1} = \frac{(w + w^{-1})(w^k + w^{-k}) - w^{k-1} - w^{1-k}}{2} = \frac{w^{k+1} + w^{-k-1}}{2} = S_{k+1} 2 x S k − S k − 1 = 2 ( w + w − 1 ) ( w k + w − k ) − w k − 1 − w 1 − k = 2 w k + 1 + w − k − 1 = S k + 1 なので、T k T_k T k と同じ漸化式と初期値を満たし、帰納法で T k ( x ) = S k T_k(x) = S_k T k ( x ) = S k です。0 < w − k ≤ 1 ≤ w k 0 \lt w^{-k} \le 1 \le w^k 0 < w − k ≤ 1 ≤ w k から不等式が従います。最後の主張は 補題 6.1 の 2 の偶奇性です。
∎
この補題は、区間 [ − 1 , 1 ] [-1,1] [ − 1 , 1 ] の中では 1 1 1 以下に抑えられている T k T_k T k が、外に出ると w k w^k w k の速さで増えることを述べています。「区間の上では小さく、区間の外の 1 点では大きい」多項式が必要な場面(共役勾配法 の収束の評価など)で、Chebyshev 多項式が最適な道具になる理由です(演習 7.4 を参照)。
定理 6.5 (Chebyshev 点の Lebesgue 定数の上界 ) ✓
n ≥ 0 n \ge 0 n ≥ 0 とし、[ − 1 , 1 ] [-1,1] [ − 1 , 1 ] の Chebyshev 点 x j = cos ( ( 2 j + 1 ) π / ( 2 n + 2 ) ) x_j = \cos\bigl((2j+1)\pi/(2n+2)\bigr) x j = cos ( ( 2 j + 1 ) π / ( 2 n + 2 ) ) (j = 0 , … , n j = 0, \ldots, n j = 0 , … , n )の Lebesgue 定数を Λ n \Lambda_n Λ n とする。このとき
Λ n ≤ 2 π log ( 2 n + 1 ) + 1 + 4 π \Lambda_n \le \frac{2}{\pi}\log(2n+1) + 1 + \frac{4}{\pi} Λ n ≤ π 2 log ( 2 n + 1 ) + 1 + π 4 である。[ a , b ] [a,b] [ a , b ] の Chebyshev 点(系 6.3 )でも同じ評価が成り立つ。
証明は長いので Appendix に置きます。最後の主張は、ℓ j \ell_j ℓ j が 1 次の変数変換 t ( x ) t(x) t ( x ) で [ − 1 , 1 ] [-1,1] [ − 1 , 1 ] の場合の ℓ j \ell_j ℓ j に移ることから従います。等間隔点の指数的な増大と比べると、Chebyshev 点の Lebesgue 定数は n = 1000 n = 1000 n = 1000 でもこの評価から 7.2 7.2 7.2 以下です。
例 6.6 (Lebesgue 定数を数値で求める )
格子(各小区間を 20 等分した点)の上で Lebesgue 関数の最大を求め、等間隔点では 命題 5.4 の下界と、Chebyshev 点では 定理 6.5 の上界と比べます。等間隔点の値は ℓ j \ell_j ℓ j の積に大きな数と小さな数が混ざるので、丸めの影響でないことを mpmath の 40 桁の計算で確かめます。格子上の最大は Λ n \Lambda_n Λ n を下から近似する値です。
# Python 3.12.10 / NumPy 2.5.2 / mpmath 1.3.0
from mpmath import mp, mpf
def leb_np ( nd , xs ) : # 格子 xs 上の Lebesgue 関数の最大
L = np. ones ( ( len ( nd ) , len ( xs ) ) )
L[j] *= (xs - nd[k]) / (nd[j] - nd[k])
return np. abs ( L ). sum ( axis = 0 ). max ()
def leb_mp ( nd , xs ) : # 同じ量を mpmath の 40 桁で
nd, best = [ mpf ( t ) for t in nd ] , mpf ( 0 )
s = sum ( abs ( mp. fprod ( (x - b) / (a - b) for b in nd if b != a ))
print ( " n 等間隔(倍) 等間隔(40桁) 下界 Cheb 上界 " )
for n in [ 5 , 10 , 20 , 30 ] :
xs = np. linspace ( - 1 , 1 , 20 * n + 1 ) # 各小区間を 20 等分
eq = np. linspace ( - 1 , 1 , n + 1 )
ch = np. cos ( ( 2 * np. arange ( n + 1 ) + 1 ) * np.pi / ( 2 * n + 2 ) )
lo = 2.0 ** n / ( 4 * n * n)
up = 2 / np.pi * np. log ( 2 * n + 1 ) + 1 + 4 / np.pi
print ( f " {n :3d } { leb_np ( eq , xs ) :.4e } { float ( leb_mp ( eq , xs )) :.4e } "
f " {lo :.2e } { leb_np ( ch , xs ) :.3f } {up :.3f } " )
n 等間隔(倍) 等間隔(40桁) 下界 Cheb 上界
5 3.1050e+00 3.1050e+00 3.20e-01 2.104 3.800
10 2.9891e+01 2.9891e+01 2.56e+00 2.489 4.211
20 1.0987e+04 1.0987e+04 6.55e+02 2.901 4.637
30 6.5710e+06 6.5710e+06 2.98e+05 3.149 4.890
検算. (a) n = 2 n = 2 n = 2 の場合は手で計算でき(演習 7.3 )、等間隔点で 5 / 4 5/4 5/4 、Chebyshev 点で 5 / 3 5/3 5/3 です。(b) 等間隔点の値は倍精度と 40 桁で表示の 5 桁まで一致し、丸め誤差の影響ではありません。(c) 等間隔点の値は下界 2 n / ( 4 n 2 ) 2^n/(4n^2) 2 n / ( 4 n 2 ) の上にあり、n n n を 10 増やすごとに数百倍(20 → 30 20 \to 30 20 → 30 で約 600 600 600 倍、1 増やすごとにおよそ 1.9 1.9 1.9 倍)になって、下界の 2 n 2^n 2 n に近い速さで増えています。Chebyshev 点の値は上界の下にあり、n n n を 10 増やしても 0.25 0.25 0.25 〜0.41 0.41 0.41 ずつしか増えません。参考までに、よく引用される評価 ( 2 / π ) log ( n + 1 ) + 1 (2/\pi)\log(n+1) + 1 ( 2/ π ) log ( n + 1 ) + 1 は n = 30 n = 30 n = 30 で 3.186 3.186 3.186 で、表の 3.149 3.149 3.149 はその下にもあります。
Lebesgue 定数と最良近似を組み合わせると、Chebyshev 点での補間が解析関数に対して幾何級数的に収束することが、導関数を使わずに示せます。
次の系の「Bernstein 楕円 E ρ E_\rho E ρ 」は、複素数の平面で [ − 1 , 1 ] [-1,1] [ − 1 , 1 ] を囲む、焦点が ± 1 \pm1 ± 1 の楕円で、ρ \rho ρ が大きいほど大きな楕円です。「正則」は複素数の変数で微分できることで、f f f が E ρ E_\rho E ρ の内側でこの意味でなめらかなら、ρ \rho ρ が大きいほど補間の誤差が速く減る、というのが系の内容です。
系 6.7 (解析関数の Chebyshev 補間 ) ✓
ρ > 1 \rho \gt 1 ρ > 1 とし、f f f は Bernstein 楕円 E ρ E_\rho E ρ で正則、[ − 1 , 1 ] [-1,1] [ − 1 , 1 ] 上で実数値とする。p n p_n p n を [ − 1 , 1 ] [-1,1] [ − 1 , 1 ] の Chebyshev 点での f f f の補間多項式とすると、1 < r < ρ 1 \lt r \lt \rho 1 < r < ρ と M r = max ∂ E r ∣ f ∣ M_r = \max_{\partial E_r}|f| M r = max ∂ E r ∣ f ∣ について
max [ − 1 , 1 ] ∣ f − p n ∣ ≤ ( 2 + 4 π + 2 π log ( 2 n + 1 ) ) 2 M r r − n r − 1 \max_{[-1,1]}|f - p_n| \le \Bigl(2 + \frac4\pi + \frac2\pi\log(2n+1)\Bigr)\frac{2M_r\,r^{-n}}{r - 1} [ − 1 , 1 ] max ∣ f − p n ∣ ≤ ( 2 + π 4 + π 2 log ( 2 n + 1 ) ) r − 1 2 M r r − n が成り立つ。
証明(系 6.7) Chebyshev 展開の打ち切り(系 3.2[誤差保証つきの求積と関数近似] 。係数の減衰は 定理 3.1[誤差保証つきの求積と関数近似] )∑ k = 0 n a k T k \sum_{k=0}^{n}a_kT_k ∑ k = 0 n a k T k は P n \mathcal{P}_n P n の元で、f f f との差は [ − 1 , 1 ] [-1,1] [ − 1 , 1 ] 上で 2 M r r − n / ( r − 1 ) 2M_rr^{-n}/(r-1) 2 M r r − n / ( r − 1 ) 以下です。よって E n ( f ) ≤ 2 M r r − n / ( r − 1 ) E_n(f) \le 2M_rr^{-n}/(r-1) E n ( f ) ≤ 2 M r r − n / ( r − 1 ) です。定理 5.2 の 2 と 定理 6.5 から主張を得ます。
∎
対数の因子は幾何級数の減衰に比べて無視できるので、誤差はおよそ r − n r^{-n} r − n の速さで減ります。r r r は ρ \rho ρ にいくらでも近く取れるので、収束の速さを決めるのは f f f が正則でいられる最大の Bernstein 楕円の大きさ ρ \rho ρ です。Runge の関数は ± i / 5 \pm i/5 ± i /5 に極をもち、ρ = 1 / 5 + 1 + 1 / 25 = ( 1 + 26 ) / 5 ≈ 1.2198 \rho = 1/5 + \sqrt{1 + 1/25} = (1 + \sqrt{26})/5 \approx 1.2198 ρ = 1/5 + 1 + 1/25 = ( 1 + 26 ) /5 ≈ 1.2198 です。
例 6.8 (Chebyshev 補間の誤差の減り方 )
e x e^x e x には 系 6.3 の評価 e / ( 2 n ( n + 1 ) ! ) e/(2^n(n+1)!) e / ( 2 n ( n + 1 )!) を、Runge の関数には 系 6.7 の予言する速さ ρ − n \rho^{-n} ρ − n を当てはめます。
# Python 3.12.10 / NumPy 2.5.2 / SciPy 1.18.1
from scipy.interpolate import BarycentricInterpolator as BI
xx = np. linspace ( - 1 , 1 , 4001 )
def cheb_err ( g , n ) : # T_{n+1} の零点で補間
x = np. cos ( ( 2 * k + 1 ) * np.pi / ( 2 * n + 2 ) )
return np. max ( np. abs ( BI ( x , g ( x ))( xx ) - g ( xx )))
print ( " n e^x の誤差 評価 e/(2^n (n+1)!) 比 " )
for n in [ 2 , 4 , 6 , 8 , 10 ] :
b = math.e / ( 2 ** n * math. factorial ( n + 1 ))
print ( f " {n :3d } {e :.2e } {b :.2e } {e / b :.3f } " )
rho = ( 1 + math. sqrt ( 26 )) / 5 # Bernstein 楕円の大きさ
runge = lambda x : 1 / ( 1 + 25 * x ** 2 )
print ( " n Runge の誤差 誤差×ρ^n " )
for n in [ 10 , 20 , 40 , 80 , 120 ] :
print ( f " {n :3d } {e :.2e } {e * rho ** n :.3f } " )
print ( f "ρ = {rho :.4f } , ρ^(-10) = {rho **- 10 :.4f } " )
n e^x の誤差 評価 e/(2^n (n+1)!) 比
2 5.65e-02 1.13e-01 0.499
4 6.40e-04 1.42e-03 0.452
6 3.62e-06 8.43e-06 0.430
8 1.22e-08 2.93e-08 0.417
10 2.71e-11 6.65e-11 0.408
ρ = 1.2198, ρ^(-10) = 0.1371
検算. (a) e x e^x e x の誤差はすべて評価 e / ( 2 n ( n + 1 ) ! ) e/(2^n(n+1)!) e / ( 2 n ( n + 1 )!) 以下で(比が 1 未満)、比がほぼ一定なので、誤差は評価式と同じ速さで減っています。比が 1 1 1 より小さいのは、評価が ∣ f ( n + 1 ) ∣ |f^{(n+1)}| ∣ f ( n + 1 ) ∣ の最大値 e e e を使っているためです。(b) Runge の関数の n = 10 , 20 , 40 n = 10, 20, 40 n = 10 , 20 , 40 の誤差は、別の点の配置と並べた先の例(例 4.3 )の Chebyshev 点の列と一致します。そこでは別の算法(chebinterpolate)とも一致していました。(c) 誤差に ρ n \rho^n ρ n を掛けた値は 0.82 0.82 0.82 付近で一定になり、誤差がちょうど ρ − n \rho^{-n} ρ − n の速さで減ること、すなわち n n n を 10 増やすごとに ρ − 10 ≈ 0.137 \rho^{-10} \approx 0.137 ρ − 10 ≈ 0.137 倍になることを示しています。これは 系 6.7 の予言する速さ(任意の r < ρ r \lt \rho r < ρ について r − n r^{-n} r − n )と合っています。
演習 7.1 易
f ( x ) = x 3 f(x) = x^3 f ( x ) = x 3 の、点 0 , 1 , 2 , 3 0, 1, 2, 3 0 , 1 , 2 , 3 での差分商の表を定理 3.3 の漸化式で作り、Newton 形が x 3 x^3 x 3 に一致すること、f [ 0 , 1 , 2 , 3 ] f[0,1,2,3] f [ 0 , 1 , 2 , 3 ] が定理 3.5 の f ′ ′ ′ ( ξ ) / 3 ! f'''(\xi)/3! f ′′′ ( ξ ) /3 ! に等しいことを確かめてください。
解答 値は 0 , 1 , 8 , 27 0, 1, 8, 27 0 , 1 , 8 , 27 、1 階の差分商は 1 , 7 , 19 1, 7, 19 1 , 7 , 19 、2 階は ( 7 − 1 ) / 2 = 3 (7-1)/2 = 3 ( 7 − 1 ) /2 = 3 、( 19 − 7 ) / 2 = 6 (19-7)/2 = 6 ( 19 − 7 ) /2 = 6 、3 階は ( 6 − 3 ) / 3 = 1 (6-3)/3 = 1 ( 6 − 3 ) /3 = 1 です。Newton 形は
0 + 1 ⋅ x + 3 x ( x − 1 ) + 1 ⋅ x ( x − 1 ) ( x − 2 ) = x + 3 x 2 − 3 x + x 3 − 3 x 2 + 2 x = x 3 0 + 1\cdot x + 3x(x-1) + 1\cdot x(x-1)(x-2) = x + 3x^2 - 3x + x^3 - 3x^2 + 2x = x^3 0 + 1 ⋅ x + 3 x ( x − 1 ) + 1 ⋅ x ( x − 1 ) ( x − 2 ) = x + 3 x 2 − 3 x + x 3 − 3 x 2 + 2 x = x 3 で一致します。f ′ ′ ′ = 6 f''' = 6 f ′′′ = 6 なので f ′ ′ ′ ( ξ ) / 3 ! = 1 = f [ 0 , 1 , 2 , 3 ] f'''(\xi)/3! = 1 = f[0,1,2,3] f ′′′ ( ξ ) /3 ! = 1 = f [ 0 , 1 , 2 , 3 ] で、ξ \xi ξ はどこでもかまいません。
演習 7.2 標準
[ 0 , 1 ] [0,1] [ 0 , 1 ] の Chebyshev 点で e x e^x e x を補間します。系 6.3 の評価で誤差 10 − 8 10^{-8} 1 0 − 8 以下が保証される最小の n n n を求めてください。
解答 [ 0 , 1 ] [0,1] [ 0 , 1 ] 上 ∣ f ( n + 1 ) ∣ ≤ e |f^{(n+1)}| \le e ∣ f ( n + 1 ) ∣ ≤ e なので、評価は 2 e ( 1 / 4 ) n + 1 / ( n + 1 ) ! 2e(1/4)^{n+1}/(n+1)! 2 e ( 1/4 ) n + 1 / ( n + 1 )! です。n = 6 n = 6 n = 6 では 2 e / ( 4 7 ⋅ 7 ! ) = 2 e / 82,575,360 ≈ 6.6 × 10 − 8 2e/(4^7 \cdot 7!) = 2e/82{,}575{,}360 \approx 6.6 \times 10^{-8} 2 e / ( 4 7 ⋅ 7 !) = 2 e /82 , 575 , 360 ≈ 6.6 × 1 0 − 8 で足りず、n = 7 n = 7 n = 7 では 2 e / ( 4 8 ⋅ 8 ! ) = 2 e / 2,642,411,520 ≈ 2.1 × 10 − 9 2e/(4^8\cdot 8!) = 2e/2{,}642{,}411{,}520 \approx 2.1 \times 10^{-9} 2 e / ( 4 8 ⋅ 8 !) = 2 e /2 , 642 , 411 , 520 ≈ 2.1 × 1 0 − 9 です。よって n = 7 n = 7 n = 7 (8 点)です。
演習 7.3 標準
[ − 1 , 1 ] [-1,1] [ − 1 , 1 ] の 3 点について Lebesgue 定数 Λ 2 \Lambda_2 Λ 2 を求めてください。(i) 等間隔点 − 1 , 0 , 1 -1, 0, 1 − 1 , 0 , 1 。(ii) Chebyshev 点 − 3 / 2 , 0 , 3 / 2 -\sqrt3/2, 0, \sqrt3/2 − 3 /2 , 0 , 3 /2 。どちらが小さいか、その理由も考えてください。
解答 (i) ℓ 0 = x ( x − 1 ) / 2 \ell_0 = x(x-1)/2 ℓ 0 = x ( x − 1 ) /2 、ℓ 1 = 1 − x 2 \ell_1 = 1 - x^2 ℓ 1 = 1 − x 2 、ℓ 2 = x ( x + 1 ) / 2 \ell_2 = x(x+1)/2 ℓ 2 = x ( x + 1 ) /2 です。対称性から [ 0 , 1 ] [0,1] [ 0 , 1 ] で考えると λ 2 = x ( 1 − x ) / 2 + ( 1 − x 2 ) + x ( x + 1 ) / 2 = 1 + x − x 2 \lambda_2 = x(1-x)/2 + (1 - x^2) + x(x+1)/2 = 1 + x - x^2 λ 2 = x ( 1 − x ) /2 + ( 1 − x 2 ) + x ( x + 1 ) /2 = 1 + x − x 2 で、最大は x = 1 / 2 x = 1/2 x = 1/2 での 5 / 4 5/4 5/4 です。
(ii) s = 3 / 2 s = \sqrt3/2 s = 3 /2 とすると ℓ 0 = x ( x − s ) / ( 2 s 2 ) \ell_0 = x(x-s)/(2s^2) ℓ 0 = x ( x − s ) / ( 2 s 2 ) 、ℓ 1 = ( s 2 − x 2 ) / s 2 \ell_1 = (s^2 - x^2)/s^2 ℓ 1 = ( s 2 − x 2 ) / s 2 、ℓ 2 = x ( x + s ) / ( 2 s 2 ) \ell_2 = x(x+s)/(2s^2) ℓ 2 = x ( x + s ) / ( 2 s 2 ) です。0 ≤ x ≤ s 0 \le x \le s 0 ≤ x ≤ s では λ 2 = 1 + ( x s − x 2 ) / s 2 \lambda_2 = 1 + (xs - x^2)/s^2 λ 2 = 1 + ( x s − x 2 ) / s 2 で最大は x = s / 2 x = s/2 x = s /2 での 5 / 4 5/4 5/4 、s ≤ x ≤ 1 s \le x \le 1 s ≤ x ≤ 1 では ℓ 1 ≤ 0 \ell_1 \le 0 ℓ 1 ≤ 0 なので λ 2 = 2 x 2 / s 2 − 1 \lambda_2 = 2x^2/s^2 - 1 λ 2 = 2 x 2 / s 2 − 1 で、最大は x = 1 x = 1 x = 1 での 8 / 3 − 1 = 5 / 3 8/3 - 1 = 5/3 8/3 − 1 = 5/3 です。よって Λ 2 = 5 / 3 \Lambda_2 = 5/3 Λ 2 = 5/3 です。
n = 2 n = 2 n = 2 では等間隔点のほうが小さく、Chebyshev 点(第 1 種)は端点 ± 1 \pm1 ± 1 を含まないので、端の近くで補間が外挿になって大きくなります。Chebyshev 点の利点は n n n が大きいときの対数的な増え方にあり、数値の例(例 6.6 )では n = 5 n = 5 n = 5 ですでに逆転しています。
演習 7.4 難
0 < a < b 0 \lt a \lt b 0 < a < b 、κ = b / a \kappa = b/a κ = b / a とします。p ( 0 ) = 1 p(0) = 1 p ( 0 ) = 1 を満たす k k k 次以下の多項式 p p p で
max x ∈ [ a , b ] ∣ p ( x ) ∣ ≤ 2 ( κ − 1 κ + 1 ) k \max_{x \in [a,b]}|p(x)| \le 2\Bigl(\frac{\sqrt\kappa - 1}{\sqrt\kappa + 1}\Bigr)^k x ∈ [ a , b ] max ∣ p ( x ) ∣ ≤ 2 ( κ + 1 κ − 1 ) k となるものがあることを、補題 6.4 を使って示してください(この評価は共役勾配法の収束の評価の形です)。
解答 p ( x ) = T k ( ( b + a − 2 x ) / ( b − a ) ) / T k ( x 0 ) p(x) = T_k\bigl((b + a - 2x)/(b - a)\bigr)/T_k(x_0) p ( x ) = T k ( ( b + a − 2 x ) / ( b − a ) ) / T k ( x 0 ) 、x 0 = ( b + a ) / ( b − a ) x_0 = (b+a)/(b-a) x 0 = ( b + a ) / ( b − a ) とおきます。x 0 > 1 x_0 \gt 1 x 0 > 1 なので 補題 6.4 から T k ( x 0 ) > 0 T_k(x_0) \gt 0 T k ( x 0 ) > 0 で、p ( 0 ) = 1 p(0) = 1 p ( 0 ) = 1 です。x ∈ [ a , b ] x \in [a,b] x ∈ [ a , b ] なら ( b + a − 2 x ) / ( b − a ) ∈ [ − 1 , 1 ] (b+a-2x)/(b-a) \in [-1,1] ( b + a − 2 x ) / ( b − a ) ∈ [ − 1 , 1 ] なので、補題 6.1 の 1 から分子の絶対値は 1 1 1 以下です。分母について、x 0 = ( κ + 1 ) / ( κ − 1 ) x_0 = (\kappa+1)/(\kappa-1) x 0 = ( κ + 1 ) / ( κ − 1 ) 、x 0 2 − 1 = 4 κ / ( κ − 1 ) 2 x_0^2 - 1 = 4\kappa/(\kappa-1)^2 x 0 2 − 1 = 4 κ / ( κ − 1 ) 2 なので
w = x 0 + x 0 2 − 1 = κ + 1 + 2 κ κ − 1 = ( κ + 1 ) 2 ( κ − 1 ) ( κ + 1 ) = κ + 1 κ − 1 w = x_0 + \sqrt{x_0^2 - 1} = \frac{\kappa + 1 + 2\sqrt\kappa}{\kappa - 1} = \frac{(\sqrt\kappa + 1)^2}{(\sqrt\kappa - 1)(\sqrt\kappa + 1)} = \frac{\sqrt\kappa + 1}{\sqrt\kappa - 1} w = x 0 + x 0 2 − 1 = κ − 1 κ + 1 + 2 κ = ( κ − 1 ) ( κ + 1 ) ( κ + 1 ) 2 = κ − 1 κ + 1 で、補題 6.4 から T k ( x 0 ) ≥ w k / 2 T_k(x_0) \ge w^k/2 T k ( x 0 ) ≥ w k /2 です。よって max [ a , b ] ∣ p ∣ ≤ 2 w − k \max_{[a,b]}|p| \le 2w^{-k} max [ a , b ] ∣ p ∣ ≤ 2 w − k です。
L. N. Trefethen, Approximation Theory and Approximation Practice , SIAM, 2013(Extended Edition, 2019)— Chapters 3–4(Chebyshev 多項式と補間)、8(解析関数の収束)、13(等間隔点と Runge 現象)、15(Lebesgue 定数)。
K. E. Atkinson, An Introduction to Numerical Analysis , 2nd ed., Wiley, 1989 — Chapter 3(補間理論)、Chapter 4(関数近似)。
T. J. Rivlin, Chebyshev Polynomials: From Approximation Theory to Algebra and Number Theory , 2nd ed., Wiley, 1990 — Chapter 1。
山本哲朗『数値解析入門』(増訂版)サイエンス社、2003 — 補間の章。
J.-P. Berrut, L. N. Trefethen, “Barycentric Lagrange interpolation”, SIAM Review 46 (2004), 501–517. doi:10.1137/S0036144502417715
証明の流れは次のとおりです。まず ℓ j \ell_j ℓ j を θ \theta θ で表し、次に λ n \lambda_n λ n を余接(cot \cot cot )の和で上から押さえます。最後に、x x x に最も近い節点の項だけを別に扱い、残りの項を調和級数に似た和 ∑ 1 / ( 2 i − 1 ) \sum 1/(2i-1) ∑ 1/ ( 2 i − 1 ) で評価して、それが対数で抑えられることを使います。
記号. N = n + 1 N = n + 1 N = n + 1 、θ j = ( 2 j + 1 ) π / ( 2 N ) \theta_j = (2j+1)\pi/(2N) θ j = ( 2 j + 1 ) π / ( 2 N ) (j = 0 , … , n j = 0, \ldots, n j = 0 , … , n )とし、x ∈ [ − 1 , 1 ] x \in [-1,1] x ∈ [ − 1 , 1 ] を x = cos θ x = \cos\theta x = cos θ (θ ∈ [ 0 , π ] \theta \in [0,\pi] θ ∈ [ 0 , π ] )と書きます。証明に使う初等的な不等式を先に確かめます。平均値の定理(定理 3.3[平均値の定理とテイラーの定理] )から、0 < u < π / 2 0 \lt u \lt \pi/2 0 < u < π /2 で tan u = u / cos 2 c ≥ u \tan u = u/\cos^2 c \ge u tan u = u / cos 2 c ≥ u (0 < c < u 0 \lt c \lt u 0 < c < u )、したがって cot u ≤ 1 / u \cot u \le 1/u cot u ≤ 1/ u です。cot \cot cot は ( 0 , π / 2 ] (0, \pi/2] ( 0 , π /2 ] で正、導関数 − 1 / sin 2 u < 0 -1/\sin^2 u \lt 0 − 1/ sin 2 u < 0 なので狭義減少です(系 3.4[平均値の定理とテイラーの定理] )。y > 0 y \gt 0 y > 0 で log ( 1 + y ) = y / ( 1 + c ) ≥ y / ( 1 + y ) \log(1+y) = y/(1+c) \ge y/(1+y) log ( 1 + y ) = y / ( 1 + c ) ≥ y / ( 1 + y ) (0 < c < y 0 \lt c \lt y 0 < c < y )です。同じく平均値の定理から、実数 a , b a, b a , b について ∣ cos N a − cos N b ∣ = N ∣ sin N c ∣ ∣ a − b ∣ ≤ N ∣ a − b ∣ |\cos Na - \cos Nb| = N|\sin Nc|\,|a - b| \le N|a - b| ∣ cos N a − cos N b ∣ = N ∣ sin N c ∣ ∣ a − b ∣ ≤ N ∣ a − b ∣ (c c c は a a a と b b b の間)です。
証明(定理 6.5) ℓ j \ell_j ℓ j の表示. 系 6.3 の証明のとおり ω = 2 − n T N \omega = 2^{-n}T_N ω = 2 − n T N で、積の微分から ω ′ ( x j ) = ∏ k ≠ j ( x j − x k ) \omega'(x_j) = \prod_{k \ne j}(x_j - x_k) ω ′ ( x j ) = ∏ k = j ( x j − x k ) なので、x ≠ x j x \ne x_j x = x j で ℓ j ( x ) = ω ( x ) / ( ( x − x j ) ω ′ ( x j ) ) = T N ( x ) / ( ( x − x j ) T N ′ ( x j ) ) \ell_j(x) = \omega(x)/\bigl((x - x_j)\omega'(x_j)\bigr) = T_N(x)/\bigl((x - x_j)T_N'(x_j)\bigr) ℓ j ( x ) = ω ( x ) / ( ( x − x j ) ω ′ ( x j ) ) = T N ( x ) / ( ( x − x j ) T N ′ ( x j ) ) です。補題 6.1 の 1 の両辺を θ \theta θ で微分すると T N ′ ( cos θ ) sin θ = N sin N θ T_N'(\cos\theta)\sin\theta = N\sin N\theta T N ′ ( cos θ ) sin θ = N sin N θ で、sin N θ j = sin ( ( 2 j + 1 ) π / 2 ) = ( − 1 ) j \sin N\theta_j = \sin\bigl((2j+1)\pi/2\bigr) = (-1)^j sin N θ j = sin ( ( 2 j + 1 ) π /2 ) = ( − 1 ) j から ∣ T N ′ ( x j ) ∣ = N / sin θ j |T_N'(x_j)| = N/\sin\theta_j ∣ T N ′ ( x j ) ∣ = N / sin θ j です。よって x = cos θ ≠ x j x = \cos\theta \ne x_j x = cos θ = x j で
∣ ℓ j ( x ) ∣ = ∣ cos N θ ∣ N ⋅ sin θ j ∣ cos θ − cos θ j ∣ |\ell_j(x)| = \frac{|\cos N\theta|}{N}\cdot\frac{\sin\theta_j}{|\cos\theta - \cos\theta_j|} ∣ ℓ j ( x ) ∣ = N ∣ cos N θ ∣ ⋅ ∣ cos θ − cos θ j ∣ sin θ j です。x x x がどれかの節点なら λ n ( x ) = 1 \lambda_n(x) = 1 λ n ( x ) = 1 なので、以下 x x x は節点でないとします。
余接への分解. α = ( θ j + θ ) / 2 \alpha = (\theta_j + \theta)/2 α = ( θ j + θ ) /2 、β = ( θ j − θ ) / 2 \beta = (\theta_j - \theta)/2 β = ( θ j − θ ) /2 とおくと、cos θ − cos θ j = 2 sin α sin β \cos\theta - \cos\theta_j = 2\sin\alpha\sin\beta cos θ − cos θ j = 2 sin α sin β 、sin θ j = sin ( α + β ) = sin α cos β + cos α sin β \sin\theta_j = \sin(\alpha + \beta) = \sin\alpha\cos\beta + \cos\alpha\sin\beta sin θ j = sin ( α + β ) = sin α cos β + cos α sin β なので
sin θ j cos θ − cos θ j = 1 2 ( cot β + cot α ) \frac{\sin\theta_j}{\cos\theta - \cos\theta_j} = \frac12\bigl(\cot\beta + \cot\alpha\bigr) cos θ − cos θ j sin θ j = 2 1 ( cot β + cot α ) です(0 < α < π 0 \lt \alpha \lt \pi 0 < α < π 、β ∉ π Z \beta \notin \pi\mathbb{Z} β ∈ / π Z )。∣ cot ∣ |\cot| ∣ cot ∣ は偶関数で周期 π \pi π なので ∣ cot α ∣ = ∣ cot ( ( θ − ( 2 π − θ j ) ) / 2 ) ∣ |\cot\alpha| = |\cot\bigl((\theta - (2\pi - \theta_j))/2\bigr)| ∣ cot α ∣ = ∣ cot ( ( θ − ( 2 π − θ j )) /2 ) ∣ です。{ θ j } ∪ { 2 π − θ j } \{\theta_j\} \cup \{2\pi - \theta_j\} { θ j } ∪ { 2 π − θ j } はちょうど s m = ( 2 m + 1 ) π / ( 2 N ) s_m = (2m+1)\pi/(2N) s m = ( 2 m + 1 ) π / ( 2 N ) (m = 0 , … , 2 N − 1 m = 0, \ldots, 2N-1 m = 0 , … , 2 N − 1 )の全体なので、v m = ( θ − s m ) / 2 v_m = (\theta - s_m)/2 v m = ( θ − s m ) /2 とおくと
λ n ( x ) ≤ ∣ cos N θ ∣ 2 N ∑ m = 0 2 N − 1 ∣ cot v m ∣ \lambda_n(x) \le \frac{|\cos N\theta|}{2N}\sum_{m=0}^{2N-1}|\cot v_m| λ n ( x ) ≤ 2 N ∣ cos N θ ∣ m = 0 ∑ 2 N − 1 ∣ cot v m ∣ です。
点の並び. 実数 v v v と π Z \pi\mathbb{Z} π Z の距離を δ ( v ) ∈ [ 0 , π / 2 ] \delta(v) \in [0, \pi/2] δ ( v ) ∈ [ 0 , π /2 ] と書くと、∣ cot v ∣ = cot δ ( v ) |\cot v| = \cot\delta(v) ∣ cot v ∣ = cot δ ( v ) です。v m v_m v m は公差 − h -h − h (h = π / ( 2 N ) h = \pi/(2N) h = π / ( 2 N ) 、2 N h = π 2Nh = \pi 2 N h = π )の 2 N 2N 2 N 項の等差数列なので、π \pi π を法として周長 π \pi π の円周上に間隔 h h h で並ぶ 2 N 2N 2 N 点です。δ ( v m 0 ) = d \delta(v_{m_0}) = d δ ( v m 0 ) = d が最小となる m 0 m_0 m 0 を取ると d ≤ h / 2 d \le h/2 d ≤ h /2 です。x x x は節点でないので d > 0 d \gt 0 d > 0 です。m ≠ m 0 m \ne m_0 m = m 0 では v m − v m 0 = k h v_m - v_{m_0} = kh v m − v m 0 = k h で、k k k を 2 N 2N 2 N で割った余り k ′ k' k ′ は 1 , … , 2 N − 1 1, \ldots, 2N-1 1 , … , 2 N − 1 を 1 回ずつとり、δ ( k h ) = min ( k ′ , 2 N − k ′ ) h \delta(kh) = \min(k', 2N - k')\,h δ ( k h ) = min ( k ′ , 2 N − k ′ ) h です。v m 0 = p π + e v_{m_0} = p\pi + e v m 0 = p π + e (p ∈ Z p \in \mathbb{Z} p ∈ Z 、∣ e ∣ = d |e| = d ∣ e ∣ = d )と書くと δ ( v m ) = δ ( e + k h ) ≥ δ ( k h ) − d \delta(v_m) = \delta(e + kh) \ge \delta(kh) - d δ ( v m ) = δ ( e + k h ) ≥ δ ( k h ) − d なので、i = min ( k ′ , 2 N − k ′ ) i = \min(k', 2N - k') i = min ( k ′ , 2 N − k ′ ) として δ ( v m ) ≥ ( i − 1 / 2 ) h \delta(v_m) \ge (i - 1/2)h δ ( v m ) ≥ ( i − 1/2 ) h です。i i i は 1 , … , N − 1 1, \ldots, N-1 1 , … , N − 1 を 2 回ずつ、N N N を 1 回とります。cot \cot cot の減少性から
∑ m ≠ m 0 ∣ cot v m ∣ ≤ 2 ∑ i = 1 N cot ( 2 i − 1 ) π 4 N ≤ 2 ∑ i = 1 N 4 N ( 2 i − 1 ) π \sum_{m \ne m_0}|\cot v_m| \le 2\sum_{i=1}^{N}\cot\frac{(2i-1)\pi}{4N} \le 2\sum_{i=1}^{N}\frac{4N}{(2i-1)\pi} m = m 0 ∑ ∣ cot v m ∣ ≤ 2 i = 1 ∑ N cot 4 N ( 2 i − 1 ) π ≤ 2 i = 1 ∑ N ( 2 i − 1 ) π 4 N です。
最も近い項. cos ( N s m 0 ) = cos ( ( 2 m 0 + 1 ) π / 2 ) = 0 \cos(Ns_{m_0}) = \cos\bigl((2m_0+1)\pi/2\bigr) = 0 cos ( N s m 0 ) = cos ( ( 2 m 0 + 1 ) π /2 ) = 0 で、θ − ( s m 0 + 2 p π ) = 2 e \theta - (s_{m_0} + 2p\pi) = 2e θ − ( s m 0 + 2 p π ) = 2 e なので、∣ cos N θ ∣ = ∣ cos N θ − cos N ( s m 0 + 2 p π ) ∣ ≤ 2 N d |\cos N\theta| = |\cos N\theta - \cos N(s_{m_0} + 2p\pi)| \le 2Nd ∣ cos N θ ∣ = ∣ cos N θ − cos N ( s m 0 + 2 p π ) ∣ ≤ 2 N d です。よって
∣ cos N θ ∣ 2 N ∣ cot v m 0 ∣ ≤ 2 N d 2 N cot d ≤ d ⋅ 1 d = 1 \frac{|\cos N\theta|}{2N}|\cot v_{m_0}| \le \frac{2Nd}{2N}\cot d \le d \cdot \frac1d = 1 2 N ∣ cos N θ ∣ ∣ cot v m 0 ∣ ≤ 2 N 2 N d cot d ≤ d ⋅ d 1 = 1 です。
まとめ. 残りの項は ∣ cos N θ ∣ ≤ 1 |\cos N\theta| \le 1 ∣ cos N θ ∣ ≤ 1 と上の評価から 1 2 N ⋅ 8 N π ∑ i = 1 N 1 2 i − 1 = 4 π ∑ i = 1 N 1 2 i − 1 \frac{1}{2N}\cdot\frac{8N}{\pi}\sum_{i=1}^{N}\frac{1}{2i-1} = \frac4\pi\sum_{i=1}^{N}\frac1{2i-1} 2 N 1 ⋅ π 8 N ∑ i = 1 N 2 i − 1 1 = π 4 ∑ i = 1 N 2 i − 1 1 以下です。i ≥ 2 i \ge 2 i ≥ 2 では y = 2 / ( 2 i − 3 ) y = 2/(2i-3) y = 2/ ( 2 i − 3 ) として log 2 i − 1 2 i − 3 = log ( 1 + y ) ≥ y 1 + y = 2 2 i − 1 \log\frac{2i-1}{2i-3} = \log(1+y) \ge \frac{y}{1+y} = \frac{2}{2i-1} log 2 i − 3 2 i − 1 = log ( 1 + y ) ≥ 1 + y y = 2 i − 1 2 なので、
∑ i = 1 N 1 2 i − 1 ≤ 1 + 1 2 ∑ i = 2 N log 2 i − 1 2 i − 3 = 1 + 1 2 log ( 2 N − 1 ) \sum_{i=1}^{N}\frac1{2i-1} \le 1 + \frac12\sum_{i=2}^{N}\log\frac{2i-1}{2i-3} = 1 + \frac12\log(2N-1) i = 1 ∑ N 2 i − 1 1 ≤ 1 + 2 1 i = 2 ∑ N log 2 i − 3 2 i − 1 = 1 + 2 1 log ( 2 N − 1 ) です。2 N − 1 = 2 n + 1 2N - 1 = 2n+1 2 N − 1 = 2 n + 1 なので λ n ( x ) ≤ 1 + 4 π + 2 π log ( 2 n + 1 ) \lambda_n(x) \le 1 + \frac4\pi + \frac2\pi\log(2n+1) λ n ( x ) ≤ 1 + π 4 + π 2 log ( 2 n + 1 ) で、x x x について最大をとれば主張を得ます。
[ a , b ] [a,b] [ a , b ] の Chebyshev 点 x j = a + b 2 + b − a 2 x ^ j x_j = \frac{a+b}{2} + \frac{b-a}{2}\hat x_j x j = 2 a + b + 2 b − a x ^ j では、t ( x ) = ( 2 x − a − b ) / ( b − a ) t(x) = (2x - a - b)/(b - a) t ( x ) = ( 2 x − a − b ) / ( b − a ) について ℓ j ( x ) = ℓ ^ j ( t ( x ) ) \ell_j(x) = \hat\ell_j(t(x)) ℓ j ( x ) = ℓ ^ j ( t ( x )) (ℓ ^ j \hat\ell_j ℓ ^ j は [ − 1 , 1 ] [-1,1] [ − 1 , 1 ] の点 x ^ j \hat x_j x ^ j の Lagrange 基底)です。x − x k = b − a 2 ( t ( x ) − x ^ k ) x - x_k = \frac{b-a}{2}(t(x) - \hat x_k) x − x k = 2 b − a ( t ( x ) − x ^ k ) の比をとれば因子 b − a 2 \frac{b-a}{2} 2 b − a が約分されるからです。t t t は [ a , b ] [a,b] [ a , b ] を [ − 1 , 1 ] [-1,1] [ − 1 , 1 ] の上へ写すので、Lebesgue 定数は等しくなります。
∎ 執筆: Claude Code(Opus 5.5)、Claude Code(Sonnet 5.5)・監査 2 回(Claude Code)・最終 2026-10-12