УказательРазделыОбозначенияАвторО проекте


§

Вспомогательная страница к разделу ИНТЕРПОЛЯЦИЯ


Интерполяция по избыточной таблице

Предположим, что задана таблица

\begin{array}{c||c|c|c|c} x & x_1 & x_2 & \ldots & x_N \\ \hline y & y_1 & y_2 & \ldots & y_N \end{array} \quad \quad \{x_j,y_j\}_{j=1}^N \subset \mathbb R,

узлы \{x_j\}_{j=1}^N — все различны. Эта таблица однозначно определяет полином f(x) \in \mathbb R[x] такой, что \left\{f(x_j)=y_j \right\}_{j=1}^N и \deg f \le N-1. Если же оказывается, что n=\deg f < N-1, то получается, что таблица избыточна: тот же полином f(x) можно построить — например, в формах Лагранжа или Ньютона — по любой выборке из этой таблицы, состоящей из (n+1)-го узла.

Обозначим

W(x)=\prod_{j=1}^N (x-x_j)

и вычислим две последовательности:

\tau_k = \sum_{j=1}^{N} y_j \frac{x_j^{k}}{W^{\prime}(x_j)} \quad npu \ k \in \{0,\dots, 2\,n-2 \}

и, при дополнительном предположении \{ y_j\ne 0\}_{j=1}^N,

\widetilde \tau_k = \sum_{j=1}^{N} \frac{1}{y_j} \frac{x_j^{k}}{W^{\prime}(x_j)} \quad npu \ k \in \{0,\dots, 2\,n-2 \} \, .

Каждое из этих выражений является симметрической функцией пар значений (x_1,y_1),\dots, (x_N, y_N) в том смысле, что значение выражения не меняется при произвольной перестановке этих пар. Каждая из последовательностей генерирует последовательность ганкелевых полиномов

\mathcal H_1(x; \{\tau\})=\left|\begin{array}{cc} \tau_0 & \tau_1 \\ 1 & x \end{array} \right| ,\mathcal H_2(x; \{\tau\}) = \left|\begin{array}{ccc} \tau_0 & \tau_1 & \tau_2 \\ \tau_1 & \tau_2 & \tau_3 \\ 1 & x & x^2 \end{array} \right|, \dots
\mathcal H_1(x; \{\widetilde \tau\})=\left|\begin{array}{cc} \widetilde \tau_0 & \widetilde \tau_1 \\ 1 & x \end{array} \right| ,\mathcal H_2(x; \{\widetilde \tau\}) = \left|\begin{array}{ccc} \widetilde \tau_0 & \widetilde \tau_1 & \widetilde \tau_2 \\ \widetilde \tau_1 & \widetilde \tau_2 & \widetilde \tau_3 \\ 1 & x & x^2 \end{array} \right|, \dots

Посмотрим как ведут себя эти полиномы при различных таблицах значений. Начнем со второй последовательности.

Т

Теорема. Пусть таблица

\begin{array}{c||c|c|c|c} x & x_1 & x_2 & \ldots & x_N \\ \hline y & y_1 & y_2 & \ldots & y_N \end{array}

состоит из значений полинома f(x) степени n \le N-1 и, вдобавок, \{ y_j\ne 0\}_{j=1}^N. Тогда этот полином может быть представлен в виде

f(x)\equiv (-1)^{n(n+1)/2} \left(\prod_{j=1}^N y_j \right) \mathcal H_{n}(x; \{ \widetilde \tau \}) \equiv (-1)^{n(n+1)/2} \left(\prod_{j=1}^N y_j \right) \left| \begin{array}{lllll} \widetilde \tau_0 & \widetilde \tau_1 & \widetilde \tau_2 & \dots & \widetilde \tau_{n} \\ \widetilde \tau_1 & \widetilde \tau_2 & \widetilde \tau_3 & \dots & \widetilde \tau_{n+1} \\ \vdots & \vdots & \vdots & & \vdots \\ \widetilde \tau_{n-2} & \widetilde \tau_{n-1} & \widetilde \tau_{n} & \dots & \widetilde \tau_{2n-1} \\ 1 & x & x^2 & \dots & x^{n} \end{array} \right|\, .

Таким образом, теорема дает еще одно представление интерполяционного полинома — альтернативное формам Лагранжа и Ньютона. Практическая польза от такой альтернативы неочевидна, поскольку вычисление определителя, зависящего от параметра — та еще проблема! К счастью, имеются соображения, упрощающие процедуру вычисления ганкелевого полинома.

П

Пример 1. Построить интерполяционный полином по таблице

\begin{array}{c||c|c|c|c|c|c|c} x & -2 & -1 & 0 & 1 & 2 & 3 & 4 \\ \hline y & 30 & -7 & 8 & 9 & 11 & 35 & 60 \end{array}

Решение. Вычисляем последовательность

\widetilde \tau_0=\frac{\scriptstyle 48569}{\scriptstyle 19958400}, \ \widetilde \tau_1=-\frac{\scriptstyle 1501}{\scriptstyle 1247400},\ \widetilde \tau_2=\frac{\scriptstyle 1021}{\scriptstyle 249480},\ \widetilde \tau_3=\frac{\scriptstyle 1733}{\scriptstyle 311850}, \dots, \widetilde \tau_{12}=\frac{\scriptstyle 168257557}{\scriptstyle 623700}

и вычисляем первые два ганкелевых полинома

\mathcal H_1(x;\{\widetilde \tau\}) \equiv \frac{48569}{19958400}x+\frac{1501}{1247400}\, ,
\mathcal H_2(x;\{\widetilde \tau\}) \equiv \underbrace{\frac{79273}{9313920000}}_{\widetilde h_{2,0} }x^2 \underbrace{-\frac{128867}{6985440000}}_{\widetilde h_{2,1} }x\underbrace{-\frac{40927}{1746360000}}_{\widetilde h_{2,2} } \, .

Вычисление \mathcal H_3(x;\{\widetilde \tau\}) может быть организовано с помощью JJ-тождества:

\mathcal H_3(x;\{\widetilde \tau\}) \equiv - \left(\frac{\widetilde h_{3,0}}{\widetilde h_{2,0}}\right)^2 \mathcal H_1(x;\{\widetilde \tau\})+ \frac{h_{\widetilde 3,0}}{\widetilde h_{2,0}}\left(x-\frac{\widetilde h_{2,1}}{\widetilde h_{2,0}}+\frac{\widetilde h_{3,1}}{\widetilde h_{3,0}} \right)\mathcal H_2(x;\{\widetilde \tau\})

где все константы известны, кроме \widetilde h_{3,0}=H_3(\{\widetilde \tau\}) and \widetilde h_{3,1}. Для нахождения этих констант, используем разложение определителя H_3 по последней строке:

\widetilde h_{3,0}= \widetilde \tau_4 \widetilde h_{2,0}+ \widetilde \tau_3 h_{2,1}+ \widetilde \tau_2 \widetilde h_{2,2}=-\frac{7159}{111767040000} \, ,
\widetilde h_{3,1}=-(\widetilde \tau_5 \widetilde h_{2,0}+ \widetilde \tau_4 \widetilde h_{2,1}+\widetilde \tau_3 \widetilde h_{2,2})=\frac{277}{1128960000} \, .

Поэтому

\mathcal H_3(x;\{\widetilde \tau\}) \equiv \frac{1}{111767040000}(-7159\,x^3+27423\,x^2-21498\,x-40400)\, ,

Продолжая рекурсивное применение JJ-тождества, получаем ганкелевы полиномы следующих порядков:

\begin{array}{lll} \mathcal H_4(x;\{\widetilde \tau\}) & \equiv & -\frac{1}{\scriptstyle 1451520000}(4\,x^4-7\,x^3+3\,x^2-2\,x-16) \ , \\ \mathcal H_5(x;\{\widetilde \tau\}) & \equiv & \frac{1}{\scriptstyle 9313920000}\left(-\frac{77}{12}\,x^5+\frac{77}{2}\,x^4-\frac{61}{4}\,x^3-\frac{715}{6}\,x^2+124\,x+\frac{304}{3}\right) , \\ \mathcal H_6(x;\{\widetilde \tau\})& \equiv & \frac{1}{\scriptstyle 349272000} \bigg(\frac{3}{80}\,x^6-\frac{59}{80}\,x^5+\frac{51}{16}\,x^4-\frac{9}{16}\,x^3-\frac{409}{40}\,x^2+\frac{93}{10}\,x+8\bigg) \end{array}

и интерполяционный полином равен 349272000 \mathcal H_6(x;\{\widetilde \tau\}).

П

Пример 2. В предположении, что интерполяционная таблица предыдущего примера, изначально содержащая значения полинома f(x) степени \deg f(x) \le 2, впоследствии могла подвергнуться искажениям в количестве не более двух значений, определить места искажений и восстановить истинный полином f(x).

Решение. Если проверить полиномы из решения предыдущего примера на приводимость над \mathbb Z, то оказывается, что у полинома \mathcal H_4 выделяются линейные множители:

\mathcal H_4(x;\{\widetilde \tau\})\equiv -\frac{1}{1451520000}(x-2)(x+1)(4\,x^2-3\,x+8) \, .

Значения x=-1 и x=2 показывают места искажений значений полинома f(x)= 4\,x^2-3\,x+8:

f(-1) \ne - 7, f(2) \ne 11 \, .

Но значения f(x) во всех остальных узлах интерполяции совпадают с табличными!

Рассмотрим теперь другую последовательность ганкелевых полиномов — ту, что порождается последовательностью

\tau_k = \sum_{j=1}^{N} y_j \frac{x_j^{k}}{W^{\prime}(x_j)} \quad npu \ k \in \{0,\dots, 2\,n-2 \}

Поисследуем как она будет вести себя на интерполяционных таблицах — избыточных, но содержащих некоторое количество ошибочных значений.

П

Пример 3. Таблица

\begin{array}{c||c|c|c|c|c|c|c} x & -2 & \mathbf{-1} & 0 & 1 & 2 & 3 & 4 \\ \hline y & 30 & {\color{RubineRed} 12 } & 8 & 9 & 18 & 35 & 60 \end{array}

состоит из значений полинома f(x)=4\,x^2-3\, x+ 8 за исключением единственного ошибочного занчения при x_2=-1. Последовательность ганкелевых полиномов

\mathcal H_{1}(x; \{ \tau \})\equiv \frac{1}{40}(x+1), \ \mathcal H_{2}(x; \{ \tau \})\equiv 0, \ \mathcal H_{3}(x; \{ \tau \})\equiv -\frac{2}{5}(x+1), \dots

и можно наблюдать, что ошибочный узел оказывается корнем двух ганкелевых полиномов: \mathcal H_{1}(x;\{\tau\}) и \mathcal H_{3}(x;\{\tau\}).

Т

Теорема. Пусть полином f(x)=a_0x^n+\dots+a_n имеет степень n< N-2. Пусть таблица удовлетворяет условиям

1. y_j\ne 0 при j \in \{1,\dots, N\},

2. y_j=f(x_j) при j \in \{1,\dots, N\} \setminus \{ e \},

3. \widehat y_{e}:=f(x_{e}) \ne y_{e} и \widehat y_{e} \ne 0.

Тогда

\mathcal H_{1} (x) \equiv \frac{( y_{e} - \widehat y_{e})}{W^{\prime}(x_{e})} (x-x_{e}) \, .

Доказательство. Будем считать x_{e} = x_1 и положим \varepsilon := y_1 - \widehat y_1. Имеем:

\tau_{\ell}=\frac{x_1^{\ell}y_1}{W^{\prime}(x_1)}+\frac{x_2^{\ell}y_2}{W^{\prime}(x_2)}+\dots+\frac{x_N^{\ell}y_N}{W^{\prime}(x_N)}=
=\left(\frac{x_1^{\ell}\widehat y_1}{W^{\prime}(x_1)}+ \frac{\varepsilon x_1^{\ell}}{W^{\prime}(x_1)}\right) +\frac{x_2^{\ell}y_2}{W^{\prime}(x_2)}+\dots+\frac{x_N^{\ell}y_N}{W^{\prime}(x_N)}=
=\sum_{j=1}^N \frac{f(x_j)x_j^{\ell}}{W^{\prime}(x_j)}+\frac{\varepsilon x_1^{\ell}}{W^{\prime}(x_1)}=

и, на основании равенства Эйлера-Лагранжа:

=\frac{ \varepsilon x_1^{\ell}}{W^{\prime}(x_1)} \quad npu \quad \ell \in \{ 0 ,1\} \, .

Таким образом,

\mathcal H_{1} (x) \equiv \left|\begin{array}{cc} \tau_0 & \tau_1 \\ 1 & x \end{array} \right| \equiv \left|\begin{array}{cc} \varepsilon /W^{\prime}(x_1) & \varepsilon x_1 /W^{\prime}(x_1) \\ 1 & x \end{array} \right|=\frac{\varepsilon}{W^{\prime}(x_1)}(x-x_1) \, ,

что и доказывает утверждение теоремы.

П

Пример 4. Таблица

\begin{array}{c||c|c|c|c|c|c|c} x & -2 & \mathbf{-1} & 0 & 1 & \mathbf{2} & 3 & 4 \\ \hline y & 30 & {\color{RubineRed} -7 } & 8 & 9 & {\color{RubineRed} 11 } & 35 & 60 \end{array}

состоит из значений полинома f(x)=4\,x^2-3\, x+ 8 за исключением двух ошибочных значений в узлах x_2=-1 and x_5=2. Последовательность ганкелевых полиномов

\mathcal H_{1}(x; \{ \tau \})\equiv \frac{1}{80}(3\,x+38), \ \mathcal H_{2}(x; \{ \tau \} )\equiv -\frac{77}{320}(x+1)(x-2), \dots

в этот раз обнаруживает ошибочные узлы в корнях полинома \mathcal H_{2}(x; \{ \tau \}).

Безнаказанно последовательно портить интерполяционную таблицу не получится.

П

Пример 5. Таблица

\begin{array}{c||c|c|c|c|c|c|c} x & -2 & \mathbf{-1} & 0 & 1 & \mathbf{2} & \mathbf{3} & 4 \\ \hline y & 30 & {\color{RubineRed} -7 } & 8 & 9 & {\color{RubineRed} 11 } & {\color{RubineRed} -1 } & 60 \end{array}

— всё та же злополучная, порожденная полиномом f(x)=4\,x^2-3\, x+ 8, но испорченная теперь уже в трёх узлах. Неповрежденных значений достаточно (даже избыточно!) для восстановления этого полинома. Тем не менее, последовательность \{\mathcal H_{k}(x; \{ \tau \}) \} уже не локализует ошибочные узлы. И это — принципиальный порог. Ту же самую таблицу можно посчитать полученной в результате порчи таблицы

\begin{array}{c||c|c|c|c|c|c|c} x & -2 & -1 & 0 & 1 & 2 & 3 & 4 \\ \hline y & -31 & -7 & 8 & 14 & 11 & -1 & -22 \end{array}

сгенерированной полиномом f_1(x)= -9/2 \, x^2+21/2\,x+8 и впоследствии испорченной в значениях f_1(-2), f_1(1) и f_1(4).


2019/02/21 10:43 редактировал au