Algoritmy interpolácie
Polynomická interpolácia je metóda[1], pri ktorej definujeme polynóm, ktorý obsahuje všetky vstupné uzlové body.
Majme množinu n+1 vstupných bodov (xi,yi) kde žiadne dve hodnoty xi nie sú rovnaké. Hľadáme polynóm p stupňa najviac n, pre ktorý platí:
- [math]p(x_i) = y_i,\; i=0,\ldots,n.[/math]
Obsah
Lagrangeov interpolačný polynóm
Lagrangeov polynóm[2], pomenovaný podľa Josepha Louisa Lagrangea, je v numerickej matematike interpolujúci polynóm pre danú množinu bodov v Lagrangeovom tvare.
Je povšimnutia hodné, že pre danú množinu bodov existuje len jeden polynóm (najmenšieho možného stupňa), ktorý interpoluje dané body. Preto je správnejšie o Lagrangeovom polynóme hovoriť ako o Lagrangeovom tvare interpolujúceho polynómu, než o Lagrangeovom interpolujúcom polynóme.
Definícia
Nech je daná množina n + 1 bodov
- [math](x_0, y_0),\ldots,(x_j, y_j),\ldots,(x_n, y_n)[/math]
kde žiadne dve hodnoty [math]x_j[/math] nie sú rovnaké. Potom interpolujúci polynóm v Lagrangeovom tvare pre túto množinu bodov je lineárna kombinácia
- [math]L(x) := \sum_{j=0}^{n} y_j \ell_j(x)[/math]
Lagrangeových bázických polynómov
- [math]\ell_j(x) := \prod_{i=0,\, i\neq j}^{n} \frac{x-x_i}{x_j-x_i} = \frac{(x-x_0)}{(x_j-x_0)} \cdots \frac{(x-x_{j-1})}{(x_j-x_{j-1})} \frac{(x-x_{j+1})}{(x_j-x_{j+1})} \cdots \frac{(x-x_{n})}{(x_j-x_{n})}.[/math]
Je povšimnutia hodné, že za predpokladu, že žiadne dve hodnoty [math]x_i[/math] nie sú rovnaké (a to ani nemôžu byť, keďže by daná úloha nedávala zmysel), platí [math]x_j - x_i \neq 0[/math], čiže daný výraz je vždy dobre definovaný.
Grafická reprezentácia
Na obrázku je kubický Lagrangeov interpolujúci polynóm L(x) (zobrazený čiernou farbou) pre body (−9, 5), (−4, 2), (−1, −2) a (7, 9). Tento polynóm je súčtom konštantných faktorov bázických polynómov y0ℓ0(x), y1ℓ1(x), y2ℓ2(x) a y3ℓ3(x). Interpolujúci polynóm prechádza cez všetky štyri body, každý z bázických polynómov prechádza jedným z daných bodov a na x-ových súradniciach daných ostatnými bodmi má hodnotu 0.
Implementácia v jazyku C
Návrh vhodnej dátovej štruktúry
Pre reprezentáciu vstupných uzlových bodov si vytvorme štruktúru Bod:
struct Bod{
double x,y;
};
kde x a y predstavujú x-ovú a y-ovú súradnicu bodu v rovine. Na reprezentáciu všetkých vstupných bodov si vytvorme pole štruktúr typu Bod o veľkosti n (n je počet vstupných uzlových bodov).
Bod *body=new Bod[n];
1 struct Bod{
2 double x,y;
3 };
4
5 double LagrangeInterpol(Bod *data, int n, double x)
6 { double lag,lag_pol=0;
7 int j,i;
8 for(j=0;j<n;j++)
9 { lag=1.0;
10 for(i=0;i<n;i++)
11 { if(i==j) continue;
12 lag*=(x-data[i].x) / (data[j].x-data[i].x) ;
13 }
14 lag_pol+= data[j].y*lag;
15 }
16 return lag_pol;
17 }
18
19 int main()
20 { int n,i;
21 Bod *body,A;
22 cin>>n;
23 body=new Bod[n];
24 for(i=0;i<n;i++)
25 { cin>>body[i].x; // x-ova suradnica
26 cin>>body[i].y; // y-ova suradnica
27 }
28 cout<<"zadaj bod zaujmu";
29 cin>>A.x;
30 A.y=LagrangeInterpol(body,n,A.x);
31 cout<<"Hodnota funkcie v bode"<< A.x<<" = "<< A.y<<endl;
32 }
Newtonov interpolačný polynóm
Definícia
Nech je daná množina n + 1 bodov
- [math](x_0, y_0),\ldots,(x_j, y_j),\ldots,(x_n, y_n)[/math]
kde žiadne dve hodnoty [math]x_j[/math] nie sú rovnaké. Potom interpolujúci polynóm v Newtonovom tvare pre túto množinu bodov je lineárna kombinácia
- [math]N(x) := \sum_{j=0}^{n} a_{j} n_{j}(x)[/math]
Newtonových bázických polynómov definovaných ako
- [math]n_j(x) := \prod_{i=0}^{j-1} (x - x_i)[/math]
a koeficienty aj sú definované ako:
- [math]a_j := [y_0,\ldots,y_j][/math]
kde
- [math][y_0,\ldots,y_j][/math]
je zápis rozdielových diferencií.
Newtonov polynóm môžeme zapísať ako
- [math]N(x) = [y_0] + [y_0,y_1](x-x_0) + \cdots + [y_0,\ldots,y_k](x-x_0)(x-x_1)\cdots(x-x_{k-1}).[/math]
alebo
- [math]\begin{align} & N\left( x \right)={{f}_{0}}+\sum\limits_{j=1}^{n}{{{n}_{j-1}}\left( x \right)\left[ {{x}_{0}},{{x}_{1}},...,{{x}_{j}} \right]} \\ & {{n}_{j}}\left( x \right)=\prod\limits_{i=0}^{j}{\left( x-{{x}_{i}} \right)} \\ \end{align}[/math]
Rozdielové diferencie môžeme vyjadriť nasledovne:
Rozdielová diferencie [math]f\left[ {{x}_{0}},{{x}_{1}},{{x}_{2}},...,{{x}_{n}} \right][/math], označované aj ako [math]\left[ {{x}_{0}},{{x}_{1}},{{x}_{2}},...,{{x}_{n}} \right][/math] vypočítaná z n+1 bodov [math]x_0, x_1, ..., x_n[/math] funkcie f(x) je definovaná nasledovne:
- [math]f\left[ {{x}_{0}} \right]\equiv f\left( {{x}_{0}} \right)[/math]
- [math]f\left[ {{x}_{0}},{{x}_{1}},...,{{x}_{n}} \right]=\frac{f\left[ {{x}_{0}},...,{{x}_{n-1}} \right]-f\left[ {{x}_{1}},{{x}_{2}},...,{{x}_{n}} \right]}{{{x}_{0}}-{{x}_{n}}}[/math]
pre n>=1.
Prvé dve diferencie sú nasledovné:
- [math]f\left[ {{x}_{0}},{{x}_{1}} \right]=\frac{f\left( {{x}_{0}} \right)-f\left( {{x}_{1}} \right)}{{{x}_{0}}-{{x}_{1}}}[/math]
- [math]f\left[ {{x}_{0}},{{x}_{1}},{{x}_{2}} \right]=\frac{f\left[ {{x}_{0}},{{x}_{1}} \right]-f\left[ {{x}_{1}},{{x}_{2}} \right]}{{{x}_{0}}-{{x}_{2}}}[/math]
Implementácia v jazyku C
1 struct Bod{
2 double x,y;
3 };
4
5 double n_k(Bod *data,int j, double x)
6 {
7 double n_j=1;
8 for(int i=0;i<j;i++)
9 n_j*=(x-data[i].x);
10 return n_j;
11 }
12
13 double NewtonPol(Bod *data, int n,double x)
14 {
15 double val=0;
16 double *koef=new double[n];
17 for(int i=0;i<n;i++)
18 koef[i]=data[i].y;
19
20 for(int j=0;j<n-1;j++)
21 { for(int i=n-1;i>j;i--)
22 koef[i]=(koef[i]-koef[i-1])/(data[i].x-data[i-j-1].x);
23 }
24 for(int i=0;i<n;i++)
25 val+=koef[i]*n_k(data,i,x);
26
27 delete koef;
28 return val;
29 }
30 int main()
31 { int n,i;
32 Bod *body,A;
33 ifstream in;
34 in.open("body.txt");
35 in>>n;
36 body=new Bod[n];
37
38 for(i=0;i<n;i++)
39 { in>>body[i].x; // x-ova suradnica
40 in>>body[i].y; // y-ova suradnica
41 }
42 A.x=1.5;
43 A.y=NewtonPol(body,n,A.x);
44 cout<<"V bode "<<A.x<<" ma Netonov iterpolacny polynom hodnotu "<<A.y;
45 }