Материал из eSyr's wiki.

Предыдущая лекция | Следующая лекция

[править] Глава 2. Интерполирование и приближение функций

[править] Параграф 3. Разделённые разности

Доказательство формулы, которая разделённая разность k-го порядка f.

Задача: выразить значение функции в k-м узле через значение функции в 0-м узле и разделённую разность k-го порядка. Решается она несложно.

• k = 1:
  • f(x₀, x₁) = f(x₀)/(x₀ − x₁) + f(x₁)/(x₁ − x₀) = (f(x₁) − f(x₀))/(x₁ − x₀)
  • (x₁ − x₀) × f(x₀, x₁) = −f(x₀) + f(x₁)
  • f(x₁) = f(x₀) + (x₁ − x₀) × f(x₀, x₁) (1)
• k = 2
  • f(x₀, x₁, x₂) = f(x₀)/((x₀ − x₁)(x₀ − x₂)) + f(x₁)/((x₁ − x₀)(x₁ − x₂)) + f(x₂)/((x₂ − x₀)(x₂ − x₁))
  • ((x₂ − x₀)(x₂ − x₁)) × f(x₀, x₁, x₂) = (−f(x₀)(x₂ − x₁))/(x₀ − x₁) { = γ} + (f(x₁)(x₂ − x₀))/(x₀ − x₁) { = α} + f(x₂)
  • α = (f(x₁)(x₂ − x₀))/(x₀ − x₁) = (x₂ − x₀)/(x₀ − x₁) × (f(x₀) + (x₁ − x₀)f(x₀, x₁)) = ((x₂ − x₀)f(x₀))/(x₀ − x₁) { = β} − (x₂ − x₀)f(x₀, x₁)
  • β − γ = f(x₀) × ((x₂ − x₀ − x₂ + x₁)/(x₀ − x₁)) = −f(x₀)
  • ((x₂ − x₀)(x₂ − x₁)) × f(x₀, x₁, x₂) = −f(x₀) − (x₂ − x₀)f(x₀, x₁) + f(x₂) (2)
• По аналогии, общая формула:
  • f(x_k) = f(x₀) + (x_k − x₀)f(x₀, x₁) + (x_k − x₀)(x_k − x₁)f(x₀, x₁, x₂) + … + (x_k − x₀)(x_k − x₁)…(x_k − x_{k − 1})f(x₀, x₁, …, x_k) (3)

[править] Параграф 4. Интерполяционная формула Ньютона

{x_i}₀^N — узлы интерполирования (n + 1 штук). Есть интерполируемая функция f(x), f(x_i) = f_i, i = 0…n.

Заменим в (3) x_k на x:

N_n(x) = f(x₀) + (x − x₀)f(x₀, x₁) + (x − x₀)(x − x₁)f(x₀, x₁, x₂) + … + (x − x₀)(x − x₁)…(x − x_{n − 1})f(x₀, x₁, …, x_n)

Чтобы эта формула была интерполяционной функцией, N_n(x_i) = f(x_i), i = 0…n

N_n(x_i) = f(x₀) + (x_i − x₀)f(x₀, x₁) + (x_i − x₀)(x_i − x₁)f(x₀, x₁, x₂) + … + (x_i − x₀)(x_i − x₁)…(x_i − x_{i − 1})f(x₀, x₁, …, x_i) = f(x_i), i = 0…n

Погрешность та же самая, но записана внешне по-другому:

|ψ_Nn(x)| = |f(x) − N_n(x)| ≤ M_{n + 1}/(n + 1)! × |ω(x)|, M_{n + 1} = sup_{a ≤ x≤ b} |f(x)^{(n + 1)}|

Замечание. Когда лучше использовать формулу Ньютона: очевидны два примера, которые показывают, что в одном случае удобно применять Ньютона, в другом — Лагранжа. У нас в Ньютоне всюду стоит функция и её разделённая разность, и она напоминает формулу Тейлора, некий её дискретный аналог, а узлы мы можем увеличивать сколько хотим, и если для специальной функции у нас она рассчитана грубая сетка, и нам надо её интерполировать. Лагранж удобен, когда у нас количество значений фиксировано, например, есть несколько датчиков.

[править] Параграф 5. Интерполирование с кратными узлами. Интерполяционная формула Эрмита

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

Пусть у нас есть узлы x₀, … x_m — m + 1 узел. Информация в узлах такая: в узле может быть задано значение производных:

• a_i ∈ N ( натуральное), i = 0…m
• a_i — кратность i-го узла

Цель: построить H_n(x)

Замечание. H_n⁽ⁱ⁾ = f⁽ⁱ⁾(x_k), k = 0…m, i = 0…a_k − 1 (1)

При выполнении условия a₀ + a₁ + … + a_m = n + 1 мы можем построить полином Эрмита, причём единственный.

Полином Эрмита строится аналогично полиному Лагранжа.

Общий вид полинома Эрмита: H_n(x) = ∑_{k = 0}^m ∑_{i = 0}^{a_k − 1} c_{k, i}(x)f⁽ⁱ⁾(x_k)

• c_{k, i}(x) — полином n-й степени
• x_i — различные

В общем виде мы строить не будем, построим H₃(x), у которого будет кратным только средний узел. Он потребуется для оценки в квадратурной формуле Гаусса. В чём фишка (тонкость): позволяет получать точную оценку квадратурной формулы.

[править] Построение H₃(x)

Существует несколько способ построения H₃(x).

H₃(x):

• H₃(x₀) = f(x₀)
• H₃(x₁) = f(x₁)
• H₃(x₂) = f(x₂)
• H₃'(x₁) = f'(x₁) (1)

Это получается из Лагранжа путём предельного перехода.

H₃(x) будем искать в виде (2):

• H₃(x) = c₀(x)f(x₀) + c₁(x)f(x₁) + c₂(x)f(x₂) + b₁(x)f'(x₁) (2)
• c_i(x), b_i(x) — многочлены 3-й степени

• c₀(x) = k(x − x₁)²(x − x₂)
• c₀(x₀) = k(x₀ − x₁)²(x₀ − x₂)

Отсюда

• c₀(x) = ((x − x₁)²(x − x₂))/((x₀ − x₁)²(x₀ − x₂))

Аналогично

• c₂(x) = ((x − x₁)²(x − x₀))/((x₂ − x₁)²(x₂ − x₀))
• b₁(x) = b₀(x − x₀)(x − x₁)(x − x₂) = b₀[(x − x₀)(x − x₂)](x − x₁)
• b₁' = b₀([ ]'(x − x₁) + [])
• b₁'(x₁) = b₀(x₁ − x₀)(x₁ − x₂)
• b₁(x) = ((x − x₀)(x − x₁)(x − x₂))/((x₁ − x₀)(x₁ − x₂))
• c₁(x) = (x − x₀)(x − x₂)(ax + b)
• c₁(x₁) = 1 = (x₁ − x₀)(x₁ − x₂)(ax₁ + b)
• c₁(x) = [(x − x₀)(x − x₂)](ax + b)
• c₁'(x) = [ ]'(ax + b) + a[ ]
• c₁'(x₁) = (2x₁ − x₀ − x₂)(ax + b) + a(x₁ − x₀)(x₁ − x₂) = 0
• (2x₁ − x₀ − x₂)/((x₁ − x₀)(x₁ − x₂)) + a(x₁ − x₀)(x₁ − x₂) = 0
• a = −(2x₁ − x₀ − x₂)/((x₁ − x₀)²(x₁ − x₂)²)
• b(x₁ − x₀)(x₁ − x₂) − (x₁(2x₁ − x₀ − x₂))/((x₁ − x₀)(x₁ − x₂)) = 1
• b = 1/((x₁ − x₀)(x₁ − x₂)) × (1 + (x₁(2x₁ − x₀ − x₂))/((x₁ − x₀)²(x₁ − x₂)²))
• c₁(x) = (x − x₀)(x − x₂)/((x₁ − x₀)²(x₁ − x₂)²)) & times; [1 − ((x − x₁)(2x₁ − x₀ − x₂))/((x₁ − x₀)(x₁ − x₂))]

[править] Погрешность интерполяционной формулы H₃(x)

ω(x) = (x − x₀)(x − x₁)²(x − x₂)

Вводим g(s) = f(s) − H₃(s) − kω(s), x₀ ≤ s ≤ x₂ (3)

Выберем константу k из условия, что в некоторой несовпадающей с узлами точке g(x) = 0: f(x) − H₃(x) − kω(x) = 0, то есть k = (f(x) − H₃(x))/ω(x).

У g(s) видим 4 нуля, то есть по теореме Ролля у первой производной существует не менее 3 нулей, у второй — не менее 2, у третьей — не менее 1, а нам нужна четвёртая. Откуда возмётся ещё один ноль: у производной по теореме Ролля ноль будет не в x₁, но он есть из условия, то есть всего у первой производной не мкенее 4 нулей. Тогда существует ξ ∈ (x₀, x₂): g⁽⁴⁾(ξ) = 0.

• g⁽⁴⁾(ξ) = 0 = f⁽⁴⁾(ξ) − k × 4!
• k = f(x) − H₃(x))/ω(x) = f⁽⁴⁾(ξ)/4!
• ψ_H3(x) = f(x) − H₃(x) = f⁽⁴⁾(ξ)/4!
• M₄ = sup_{x0 ≤ x ≤ x2} |f⁽⁴⁾(x)|
• |f(x) − H₃(x)| ≤ M₄/4! × |ω(x)|

• |f(x) − H_n(x)| ≤ M_{n + 1}/(n + 1)! × |ω(x)|
• ω(x) = (x − x₀)^a₀(x − x₁)^a₁…(x − x_m)^a_m
• M_{n + 1} = sup_{a ≤ x ≤ b} |f^{(n + 1)}(x)|