Упражнение 2

страницата се нуждае от дописване/преглеждане


Пример за построяване на интерполационен полином на Лагранж за дадена функция

Следната програма на Mathematica построява интерполационния полином на Лагранж за функцията $f(x) = e^x$ , след което се изчертава графиката на грешката по модул заедно с графиката на функцията g. (интересна графика се получава :) )

n = 3;

Do[x[k] = -1 + 2 k/n, {k, 0, n}];
(* или това,или другото отдолу - дават различна точност - различни \
интерполационни възли с други думи *)

(* Do[x[k]=Cos[(2k+1) Pi/(2n+2)],{k,0,n}]; *)

Do[l[k_, t_] := (s = 1;
    Do[If[j != k, s *= (t - x[j])/(x[k] - x[j])], {j, 0, n}]; s), {k, 
   0, n}];
L[f_, t_] := Sum[l[k, t] f[x[k]], {k, 0, n}];

f[t_] := Exp[t];
(* Gamma е гама-функцията на Ойлер *)
g[t_] := E/Gamma[n + 2] Abs[Product[t - x[k], {k, 0, n}]];
Plot[{Abs[f[t] - L[f, t]], g[t]}, {t, -1, 1}]

Задача за границата на грешката

Нека f(x) притежава производни от всякакъв ред в интервала $[a,b]$, и $|f^{(n)}(x)| \le C.M^n , \forall x \in [a,b], n = 1, 2, ...$ където $C, M$ - положителни константи.
Да се докаже, че при произволен избор на интерполационни възли (за всяко n), $a \le x_0 \le x_1 \le \dots < x_n \le b$, е изпълнено, че $\max_{x \in [a, b]} |f(x) - L_n(f; x)| \rightarrow 0$, когато $n \rightarrow \infty$ .

Доказателство:

(1)
\begin{align} |f(x) - L_n(f;x)| = \cfrac {|f^{(n+1)}(\xi)|} {(n+1)!} . |\omega(x)| \le \cfrac {C.M^{n+1}} {(n+1)!} . \cfrac {|x - x_0|} {(\le b-a)} . \cfrac {|x - x_1|} {(\le b-a)} \ . \dots . \cfrac {|x - x_n|} {(\le b-a)} \le \cfrac {C(M.(b - a))^{n+1}} {(n+1)!} \rightarrow 0 \end{align}

, когато $n \rightarrow \infty$.

Защо клони към нула? Косвено доказателство:

Необходимо условие за един числов ред $\underset {k=1} {\overset \infty \sum} a_k$ да е сходящ е $lim_{k \rightarrow \infty} a_k = 0$.

Series[Exp[x], {x, 0, 10}] (* McLoren series *)

$\underset {k=0} {\overset n \sum} \cfrac {|x_0|^k} {k!}$ - числов ред, $a_k = \cfrac {|x_0|^k} {k!} , \cfrac {a_{k+1}} {a_k} = \cfrac {|x_0|} {k+1} \underset {k \rightarrow \infty} {\rightarrow} 0 < 1$, следователно редът е сходящ. $=> \cfrac {|x_0|^k} {k!} \underset {k \rightarrow \infty} {\rightarrow} 0 => lim_{k \le \infty} \cfrac {(M(b-a))^n} {n!} = 0$

Някои примери на функции, удовлетворяващи условието на задачата:

1.
$f(x) = e^{Mx} , M = const$ ; $f^{(n)} = M^n . e^{Mx}$;
$max_{x \in [a,b]} |f^{(n)}(x)| \le (max_{x \in [a,b]} e^{Mx}).|M|^n$ (първият множител е C - константа)

2.
$f(x) = sin(Mx)$ ($f'(x) = M. cos(Mx)$ )

Съществува цял клас примери. За "хубави" функции грешката клони към 0, за "нехубави" - не.
[TODO: EXPLAIN]

[TODO: WRITE]

Ще използваме свойството на разделените разлики.

Задачи: Да се докажат тъждествата:
а) $\underset {k=0} {\overset n \sum} \cfrac 1 {\omega'(x_k) . x_k^m} = 0$ за m = 0, …, n - 1
б) $\underset {k=0} {\overset n \sum} \cfrac 1 {\omega'(x_k) . x_k^n} = 1$
в) $\underset {k=0} {\overset n \sum} \cfrac {\omega''(x_k)} {\omega\(x_k)} = 0$
г) $\underset {k=0} {\overset n \sum} x_k . \cfrac {\omega''(x_k)} {\omega'(x_k)} = (n+1) * n$

Лявата страна е разделена разлика от ред n за следните функции:
а) $f(x) = x^m \in \Pi_{n-1}$
б) $f(x) = x^n = 1.x^n$
в) $f(x) = \omega''(x) \in \Pi_{n-1}$
г) $f(x) = x.\omega''(x)$ - полином от степен n с коефициент (n+1)n пред $x^n$

Предишните разсъждения малко не вървят, затова ще използваме други.

Задача: Да се намери $\underset {k=0} {\overset n \sum} \cfrac 1 {\omega'(x_k) . x_k^{n+1}} = f[x_0, ..., x_n], f(x) = x^{n+1}$

Решение 1:
$x^{n+1} - L_n(f; x) = \omega(x)$, защото лявата страна е полином от степен n + 1 със старши коефициент 1, който се анулира в $x = x_0, ..., x_n$.
От дефиницията на разделена разлика: Коефициентът пред $x^n$ отляво е $- f[x_0, ..., x_n]$, а отдясно е $\omega(x) = (x - x_0)...(x - x_n) = x^{n + 1} - (x_0 + ... + x_n).x^n + \dots$.

Така получаваме, че разделената разлика е сумата от $x_i$.

Решение 2: Ще използваме формулата за разделена разлика на произведение от две функции, от миналата лекция. Формулата се казваше Лема на Поповичу-Стефенсон.
$(f.g)[x_o, ..., x_n] = \underset {k=0} {overset n \sum} f[x_0, ..., x_k] g[x_k, ..., x_n]$.

Ще я приложим, като представим $x^{n+1}$ като $x^n.x (f.g)$ . Тогава разделената разлика от $x^{n+1}$ ще се представи като
$(f.g)[x_0, ..., x_n] = f[x_0, ..., x_{n-1}]g[x_{n-1}, x_n] + f[x_1, ..., x_n].g[x_n]$
Останалите се анулират, защото g е функция от ред 1 - разделените разлики от по-висок ред са нули, сещате се:
$g[x_k, ..., x_n] = 0$ при $k < n-1$

Получихме:
$(x^{n+1})[x_0, ..., x_n] = (x^n)[x_0, ..., x_{n-1}] + x_n$ - рекурентна връзка. Като итерираме разсъждението за по-малки n:

(2)
\begin{align} (x^{n+1})[x_0, ..., x_n] = (x^n)[x_0, ..., x_{n-1}] + x_n = (x^{n-1})[x_0, ..., x_{n-2}] + x_{n-1} + x_n = \dots = (x)[x_0] + x_1 + x_2 + ... + x_n = x_0 + x_1 + ... + x_n \end{align}

- което и искахме да докажем.

Пример

Да се намери полином $p(x) \in \Pi_n$, такъв, че $p(-1) = 0, p(0) = 3, p(1) = 2, p(2) = -1$

$L_n(f;x) = f(x_0) + \underset {k = 1} {\overset n \sum} f[x_0, ..., x_n] (x - x_0)(x - x_1)...(x - x_{k - 1})$

Правим си таблица като на лекции:

$x_i$ $p[x_i]$ $p[x_i, x[{i+1}]$ $p[x_i, x_{i+1}, x_{i+1}, x_{i+2}]$ $p[x_i, x_{i+1}, x_{i+1}, x_{i+2}, x_{i+3}]$
-1 0 3
0 3 -1 -2
1 2 -3 -1 $\frac 1 3$
2 -1

Така $p(x) = 0 + 3(x+1) - 2(x+1)x + \frac 1 3 (x + 1)x(x - 1)$

Тази схема, освен че е лесна за програмиране, има и предимството, че, ако одобавим още един възел - примерно p(3) = 6, можем да получим новия полином от степен, ненадминаваща 4, като използваме направените до момента сметки и просто добавяме по една разделена разлика във всеки стълб. В нашия случай става:

Правим си таблица като на лекции:

$x_i$ $p[x_i]$ $p[x_i, x[{i+1}]$ $p[x_i, x_{i+1}, x_{i+1}, x_{i+2}]$ $p[x_i, x_{i+1}, x_{i+1}, x_{i+2}, x_{i+3}]$ $p[x_i, x_{i+1}, x_{i+1}, x_{i+2}, x_{i+3}, x_{i+4}]$
-1 0 3
0 3 -1 -2
1 2 -3 -1 $\frac 1 3$ $\frac 5 {12}$
2 -1 7 5 2
3 6

и полиномът се получава като $p_1(x) = p(x) + \frac 5 {12} (x + 1)x(x - 1)(x - 2)$

При това няма значение къде ще вмъкнем възлите - наредбата няма значение. Важното е само възлите да са еднакви. Когато, обаче, имаме интерполация по Ермит - със съвпадащи възли - тогава трябва да се грижим схемата да е в нарастващ ред на възлите.

Unless otherwise stated, the content of this page is licensed under Creative Commons Attribution-NonCommercial-ShareAlike 3.0 License