Общая информация

Интерполирование по Эрмиту есть специальный вид многоинтервальной интерполяции, при котором интерполирующий многочлен $H_n(х)$ обеспечивает не только равенство $H_n(х)$ значениям $f(x)$, но и совпадение некоторого количества производных в узлах интерполяции.

Многочлен Эрмита единственный, это следует из однозначности его построения.

Формула для остаточного члена при интерполяции с кратными узлами для функции $f ∈ C$ $m[a, b]$ приведена в следующей теореме:

Теорема. Если функция $f(x)$ является $m = a_1+a_2+\cdots +a_n$ раз непрерывно дифференцируемой, то существует точка $ξ$ такая, что

$$r(x)= f(x) - H(x) = \frac{f^{(m)}(\xi )}{m!}\Omega (x),$$

где $\Omega (x) = (x-x_1)^{a_1}(x-x_2)^{a_2}\cdots (x-x_n)^{a_n}$.

Алгоритм построения

Представим себе, что на плоскости в координатах x,y нарисовано несколько точек. Задача соединить их плавной кривой так чтобы в итоге получился график гладкой функции y=ƒ(x).

На каждом отдельном участке между соседними точками это должен быть многочлен третьей степени и при этом соседние многочлены должны каким-то образом гладко сопрягаться.

Для того, чтобы решить эту задачу, сначала решим частную задачу: представим себе что на отрезке [0; 1] у нас задана функция - многочлен третьей степени:

ƒ(x)=а₀+a₁x+a₂x²+a₃x³

Нам неизвестны коэффициенты a₀,a₁,a₂,a₃ и мы хотим определить их по следующим данным: ƒ(0)=y₀, ƒ(1)=y₁, ƒ' (0)=z₀, ƒ'(1)=z₁

Производная функции ƒ определяется следующим образом:

ƒ'(x)=a₁+2a₂x+3a₃x²

Запишем наши уравнения:

f(0)=a₀=y₀

f'(0)=a₁=z₀

f(1)=a₀+a₁+a₂+a₃=y₁

f'(1)=a₁+2a₂+3a₃=z₁

Было бы удобно записать эти уравнения в матричном виде:
$$
P = \begin{bmatrix}
1 & 0 & 0 & 0 \\
0 & 1 & 0 & 0 \\
1 & 1 & 1 & 1 \\
0 & 1 & 2 & 3 \\
\end{bmatrix}
$$

Вектор-столбец неизвестных переменных обозначим
$$
A = \begin{bmatrix}
a_0 \\
a_1 \\
a_2 \\
a_3 \\
\end{bmatrix}
$$

Правая часть системы линейных уравнений, это столбец с координатами y₀,z₀,y₁,z₁
$$
B = \begin{bmatrix}
y_0 \\
z_0 \\
y_1 \\
z_1 \\
\end{bmatrix}
$$

Тогда нам нужно решить систему линейных уравнений P*A=B

В общем случае решением будет:

Воспользуемся алгоритмом Гаусса: для этого в левой части таблицы запишем матрицу P, a справа приписывается единичная матрица

$$
\begin{bmatrix}
\left.\begin{matrix}
1& 0 & 0 & 0 \\
0& 1& 0& 0 \\
1& 1& 1& 1\\
0& 1& 2& 3\\
\end{matrix}\right| &
\begin{matrix}
1& 0& 0& 0 \\
0& 1& 0& 0\\
0& 0& 1& 0\\
0& 0& 0& 1\\
\end{matrix}
\end{bmatrix}
$$

Задача привести матрицу некоторыми преобразованиями к такому виду, чтобы слева оказалась единичная матрица, а то, что получим в правой половине и будет обратной матрицей

Итак,

По значениям в концах отрезка [0,1] и их производных, мы восстановили всю функцию:

a₀=y₀
a₁=z₀
a₂=-3y₀-2z₀+3y₁-z₁
a₃=2y₀+z₀-2y₁+z₁
ƒ(x)=y₀+z₀x+(-3y₀-2z₀+3y₁-z₁)x²+(2y₀+z₀-2y₁+z₁)x³

Пусть нужно произвести интерполяцию на отрезке [x₀;x₁], тогда график функции необходимо сдвинуть на x₀ и растянуть в x₁-x₀ раз, а поскольку z₀,z₁ в x₁-x₀ раз меньше, чем требуется, то их необходимо домножить.

Для удобства введем k=(x-x₀)/(x₁-x₀) и h=x₁-x₀.

Функция примет вид:

ƒ(k)=y₀+z₀hk+(-3y₀-2z₀h+3y₁-z₁h)k²+(2y₀+z₀h-2y₁+z₁h)k³

Аналогичным образом определяется многочлен любой степени.

При добавлении значений производных (до некоторого порядка m) в точках добавляются дополнительные слагаемые.
Точность растет с увеличением количества известных значений производных в точках.