Algoritmy numerickej aproximácie (riešené príklady): Rozdiel medzi revíziami
(7 medziľahlých úprav od 2 ďalších používateľov nie je zobrazených) | |||
Riadok 1: | Riadok 1: | ||
− | |||
− | |||
− | |||
{{Draft}} | {{Draft}} | ||
{{Skripta programovanie (zbierka úloh)}} | {{Skripta programovanie (zbierka úloh)}} | ||
− | ==Zadanie== | + | ==Zadanie 1== |
Riešte problém aproximácie dát. K dispozícii máme ''n'' bodov v rovine (ich súradnice ''x'' a ''y''). Úlohou bude vypočítať rovnicu aproximujúcej krivky metódou najmenších štvorcov pre: | Riešte problém aproximácie dát. K dispozícii máme ''n'' bodov v rovine (ich súradnice ''x'' a ''y''). Úlohou bude vypočítať rovnicu aproximujúcej krivky metódou najmenších štvorcov pre: | ||
#Lineárnu aproximáciu v tvare <math>y=ax+b</math> | #Lineárnu aproximáciu v tvare <math>y=ax+b</math> | ||
Riadok 14: | Riadok 11: | ||
Porovnajte výsledok funkcie s výsledkom v tabuľkovom procesore, kde do (bodového) grafu zo vstupných bodov pridáte spojnicu trendu typu lineárna, resp. logaritmická, mocninová exponenciálna. | Porovnajte výsledok funkcie s výsledkom v tabuľkovom procesore, kde do (bodového) grafu zo vstupných bodov pridáte spojnicu trendu typu lineárna, resp. logaritmická, mocninová exponenciálna. | ||
− | ==Vstupné údaje== | + | ===Vstupné údaje=== |
Ako vzorku vstupných bodov budeme uvažovať body z intervalu <1,8>. | Ako vzorku vstupných bodov budeme uvažovať body z intervalu <1,8>. | ||
''Pre generovanie vstupných hodnôt použijem funkciu <math>f(x)=x^4-exp(x)</math>'' | ''Pre generovanie vstupných hodnôt použijem funkciu <math>f(x)=x^4-exp(x)</math>'' | ||
− | ===Vzorový vstup=== | + | ====Vzorový vstup==== |
{| class=datatable | {| class=datatable | ||
|- | |- | ||
Riadok 47: | Riadok 44: | ||
|1115.04 | |1115.04 | ||
|} | |} | ||
− | ===Vzorový výstup=== | + | ====Vzorový výstup==== |
#y = 213.634 x - 467.195 | #y = 213.634 x - 467.195 | ||
#y = 784.982 ln x - 572.83 | #y = 784.982 ln x - 572.83 | ||
− | #y = | + | #y = 0.9628 e^(1.054x) |
#y = 0.2609 x^(4.443) | #y = 0.2609 x^(4.443) | ||
− | ==Analýza matematických vzťahov== | + | ===Analýza matematických vzťahov=== |
Poznámka: v tomto texte nie sú uvádzané celé postupy výpočtu, nakoľko v tomto predmete nie sú primárnym cieľom matematické dôkazy. Študent si postup môže ľahko overiť sám. | Poznámka: v tomto texte nie sú uvádzané celé postupy výpočtu, nakoľko v tomto predmete nie sú primárnym cieľom matematické dôkazy. Študent si postup môže ľahko overiť sám. | ||
− | ===Lineárna aproximácia: y=ax+b=== | + | ====Lineárna aproximácia: y=ax+b==== |
Pre lineárnu aproximáciu platia vzťahy: | Pre lineárnu aproximáciu platia vzťahy: | ||
:<math>\begin{align} | :<math>\begin{align} | ||
Riadok 81: | Riadok 78: | ||
\end{align}</math> | \end{align}</math> | ||
− | ===Logaritmická aproximácia: <math>y=a*\ln x + b</math>=== | + | ====Logaritmická aproximácia: <math>y=a*\ln x + b</math>==== |
Hľadáme najlepšiu aproximáciu pre logaritmickú rovnicu: <math>y=a*\ln x + b</math>. Postup výpočtu je rovnaký ako v prvom prípade. | Hľadáme najlepšiu aproximáciu pre logaritmickú rovnicu: <math>y=a*\ln x + b</math>. Postup výpočtu je rovnaký ako v prvom prípade. | ||
Jediná zmena v riešení je tá, že os x je v logaritmickej mierke. Dostávame teda sústavu 2 rovníc: | Jediná zmena v riešení je tá, že os x je v logaritmickej mierke. Dostávame teda sústavu 2 rovníc: | ||
Riadok 102: | Riadok 99: | ||
\end{align}</math> | \end{align}</math> | ||
− | ===Exponenciálna aproximácia: <math>y=b{{e}^{ax}}</math>=== | + | ====Exponenciálna aproximácia: <math>y=b{{e}^{ax}}</math>==== |
Rovnicu <math>y=b{{e}^{ax}}</math> si upravíme tak, že ju zapíšeme v logaritmickom tvare: | Rovnicu <math>y=b{{e}^{ax}}</math> si upravíme tak, že ju zapíšeme v logaritmickom tvare: | ||
:<math>\begin{align} | :<math>\begin{align} | ||
Riadok 132: | Riadok 129: | ||
kde <math>B=\ln b</math>. | kde <math>B=\ln b</math>. | ||
− | ===Mocnimová aproximácia: <math>y=b{{x}^{a}}</math>=== | + | ====Mocnimová aproximácia: <math>y=b{{x}^{a}}</math>==== |
Rovnicu <math>y=b{{x}^{a}}</math> si upravíme tak, že ju zapíšeme v logaritmickom tvare: | Rovnicu <math>y=b{{x}^{a}}</math> si upravíme tak, že ju zapíšeme v logaritmickom tvare: | ||
:<math>\begin{align} | :<math>\begin{align} | ||
Riadok 161: | Riadok 158: | ||
kde <math>B=\ln b</math>. | kde <math>B=\ln b</math>. | ||
− | ==Analýza programátorského riešenia== | + | ===Analýza programátorského riešenia=== |
− | ===Návrh dátových štruktúr=== | + | ====Návrh dátových štruktúr==== |
V úlohe budeme pracovať s bodmi v rovine, ktoré sú charakterizované x a y súradnicou. Preto si vytvoríme štruktúru Bod, ktorá bude takýto bod reprezentovať: | V úlohe budeme pracovať s bodmi v rovine, ktoré sú charakterizované x a y súradnicou. Preto si vytvoríme štruktúru Bod, ktorá bude takýto bod reprezentovať: | ||
<source lang="c"> | <source lang="c"> | ||
Riadok 188: | Riadok 185: | ||
</source> | </source> | ||
− | ==Riešenie v jazyku C== | + | ===Riešenie v jazyku C=== |
Vytvoríme funkciu ''Krivka minSq(Bod *body, int n)'' kde v poli body sú vstupné údaje, n je počet bodov. Funkcia bude vracať dátovú struktúru typu Krivka, kde budú informácie o tvare aproximujúcej krivky. | Vytvoríme funkciu ''Krivka minSq(Bod *body, int n)'' kde v poli body sú vstupné údaje, n je počet bodov. Funkcia bude vracať dátovú struktúru typu Krivka, kde budú informácie o tvare aproximujúcej krivky. | ||
Riadok 203: | Riadok 200: | ||
Krivka k; | Krivka k; | ||
k.a=(n*sxy-sy*sx)/(n*sxx-sx*sx); | k.a=(n*sxy-sy*sx)/(n*sxx-sx*sx); | ||
− | k.b=(sy-a*sx)/n; | + | k.b=(sy-k.a*sx)/n; |
k.typ=LIN; | k.typ=LIN; | ||
return k; | return k; | ||
Riadok 254: | Riadok 251: | ||
} | } | ||
</source> | </source> | ||
− | ===Použitie v programe=== | + | ====Použitie v programe==== |
Pre názornosť nebudeme v programe dáta načítavať z klávesnice, ale vstupné údaje budú zadané priamo v zdrojovom kóde. Modifikácia zdrojového kódu s načítaním vstupných hodnôt je jednoduchá a určite ju každý zvládne sám. | Pre názornosť nebudeme v programe dáta načítavať z klávesnice, ale vstupné údaje budú zadané priamo v zdrojovom kóde. Modifikácia zdrojového kódu s načítaním vstupných hodnôt je jednoduchá a určite ju každý zvládne sám. | ||
Riadok 316: | Riadok 313: | ||
</source> | </source> | ||
− | ==Overenie vypočítaných hodnôt== | + | ===Overenie vypočítaných hodnôt=== |
Na nasledujúcich obrázkoch sú vstupné body z príkladu znázornené v tabuľkovom procesore a je vypočítaná najlepšia aproximujúca krivka (MS Excel ju nazýva spojnica trendu). | Na nasledujúcich obrázkoch sú vstupné body z príkladu znázornené v tabuľkovom procesore a je vypočítaná najlepšia aproximujúca krivka (MS Excel ju nazýva spojnica trendu). | ||
<gallery heights=220px widths=400px perrow=2> | <gallery heights=220px widths=400px perrow=2> | ||
Riadok 324: | Riadok 321: | ||
Súbor:min SQ moc.png | Súbor:min SQ moc.png | ||
</gallery> | </gallery> | ||
+ | |||
+ | ==Zadanie 2== | ||
+ | Vytvorte funkciu, ktorá pre zadaná vstupy vypočíta všetky typy aproximácií (lineárnu, logaritmickú, exponenciálnu a mocninovú) a určí ktorá aproximácia je optimálna pre zadané vstupné body. | ||
+ | |||
+ | ===Analýza=== | ||
+ | Pri vyhodnocovaní optimálnej aproximácie musíme vypočítať všetky možné aproximačné krivky a zvoliť kritérium voľby najlepšej aproximácie. Kritériom posudzovania môže byť hľadanie minimálych hodnôt: | ||
+ | *vzdialenosť uzlových (vstupných) bodov od aproximujúcej krivky | ||
+ | *štvorce vzdialeností uzlových (vstupných) bodov od aproximujúcej krivky | ||
+ | *plocha medziproximujúcou krivkou a interpolačným polynómom, ktorých prechádza cez uzlové body | ||
+ | |||
+ | ===Riešenie v jazyku C=== | ||
+ | Kriteriom na vyhodnotenie najlepšej aproximovanej krivky bude hľadanie minimálych hodnôt a vzdialeností uzlových (vstupných) bodov od aproximujúcej krivky. | ||
+ | |||
+ | Definovanie štruktúr potrebných pre spracovanie súradníc: | ||
+ | <source lang="c" line> | ||
+ | struct Bod{ | ||
+ | double x,y; | ||
+ | }; | ||
+ | |||
+ | enum TYP_KRIVKY{ | ||
+ | LIN, LOG, EXP, MOC | ||
+ | }; | ||
+ | |||
+ | struct Krivka{ | ||
+ | TYP_KRIVKY typ; | ||
+ | double a,b; | ||
+ | }; | ||
+ | </source> | ||
+ | |||
+ | Táto štruktúra pre zjednodušenie poslúži na uloženie všetkých typov aproximácie do jednej premennej typu Aprox: | ||
+ | <source lang="c" line> | ||
+ | struct Aprox{ | ||
+ | Krivka Lin, Log, Exp, Moc; | ||
+ | }; | ||
+ | </source> | ||
+ | |||
+ | Funkcia aproximacie vypočíta všetky 4 typy aproximácií. Návratová hodnota funkcie je typu Aprox: | ||
+ | <source lang="c" line> | ||
+ | Aprox aproximacie(Bod *body, int n, Aprox B) //vypocita vsetky | ||
+ | //aproximacie | ||
+ | { | ||
+ | double sx=0,sxx=0,sy=0,sxy=0,slxly; | ||
+ | double slx=0,slxx=0,sly=0,slxy=0,sxly=0; | ||
+ | for(int i=0;i<n;i++) | ||
+ | { sx+=body[i].x; | ||
+ | slx+=log(body[i].x); | ||
+ | sy+=body[i].y; | ||
+ | sly+=log(body[i].y); | ||
+ | sxy+=body[i].x * body[i].y; | ||
+ | slxy+=log(body[i].x) * body[i].y; | ||
+ | sxly+=body[i].x * log(body[i].y); | ||
+ | slxly+=log(body[i].x) *log(body[i].y); | ||
+ | sxx+=body[i].x * body[i].x; | ||
+ | slxx+=log(body[i].x) * log(body[i].x); | ||
+ | } | ||
+ | |||
+ | B.Lin.a=(n*sxy-sy*sx)/(n*sxx-sx*sx); | ||
+ | B.Lin.b=(sy-B.Lin.a*sx)/n; | ||
+ | |||
+ | B.Log.a=(n*slxy-sy*slx)/(n*slxx-slx*slx); | ||
+ | B.Log.b=(sy-B.Log.a*slx)/n; | ||
+ | |||
+ | B.Exp.a=(n*sxly-sly*sx)/(n*sxx-sx*sx); | ||
+ | B.Exp.b=(sly-B.Exp.a*sx)/n; | ||
+ | B.Exp.b=exp(B.Exp.b); | ||
+ | |||
+ | B.Moc.a=(n*slxly-sly*slx)/(n*slxx-slx*slx); | ||
+ | B.Moc.b=(sly-B.Moc.a*slx)/n; | ||
+ | B.Moc.b=exp(B.Moc.b); | ||
+ | |||
+ | return B; | ||
+ | } | ||
+ | </source> | ||
+ | |||
+ | Na vypísanie najvhodnejšieho typu rovnice je použitá funkcia vypisRovnicuKrivky: | ||
+ | <source lang="c" line> | ||
+ | void vypisRovnicuKrivky(Krivka k) //vypise rovnicu urciteho typu | ||
+ | { | ||
+ | cout<<"y = "; | ||
+ | switch(k.typ) | ||
+ | { | ||
+ | case LIN: | ||
+ | cout<<k.a<<" x"; | ||
+ | if(k.b>0) cout<<"+"; | ||
+ | cout<<k.b<<" -> LINEARNA"<<endl; | ||
+ | break; | ||
+ | case LOG: | ||
+ | cout<<k.a<<" ln x"; | ||
+ | if(k.b>0) cout<<"+"; | ||
+ | cout<<k.b<<" -> LOGARITMICKA"<<endl; | ||
+ | break; | ||
+ | case EXP: | ||
+ | cout<<k.b<<" e^("; | ||
+ | cout<<k.a<<" x)"<<" -> EXPONENCIALNA"<<endl; | ||
+ | break; | ||
+ | case MOC: | ||
+ | cout<<k.b<<" x^("; | ||
+ | cout<<k.a<<" )"<<" -> MOCNINOVA"<<endl; | ||
+ | break; | ||
+ | } | ||
+ | } | ||
+ | </source> | ||
+ | |||
+ | Pre zistenie optimálnej aproximovanej krivky pre dané vstupné body je použitý jednoduchý algoritmus daný vo funkcii zistiOptimum. | ||
+ | *Postup riešenia algoritmu spočíva v zistení vzdialenosti y-onovej súradnice vstupného bodu od aproximovanej y-onovej súradnice, ktorú musíme najskôr vypočítať, a porovnaním ako sú tieto body vzdialené od jednotlivých bodov kriviky pre každý typ, zistíme optimum. | ||
+ | *Dosadením x-ových súradníc vstupných bodov do príslušnej rovnice krivky sa vypočíta aproximovaná y-onová súradnica. Absolútnu hodnotu rozdielu vzdialenosti medzi y-onovými súradnicami(vypočítanými a vstupnými) sčíta a uloží túto sumu do premennej sumLin (pre každý ďalší jeden typ aproximácie do - sumLog, sumExp, sumMoc ). Toto sa opakuje pre každý jeden vstupný bod a pre každú rovnicu krivky. | ||
+ | *Na záver sa sumy všetkých typov aproximácií porovnajú. Najmenšia suma predstavuje najmenšie odchýlky pre danú krivku v porovnaní so vstupnými bodmi. Funkcia vráti Optimálnu krivku aproximácie s najmenšou sumou: | ||
+ | <source lang="c" line> | ||
+ | Krivka zistiOptimum(Bod *body, int n, Aprox B) | ||
+ | { | ||
+ | double sumLin=0, sumLog=0, sumExp=0, sumMoc=0, y1=0,y2=0,y3=0,y4=0, Naj=0; | ||
+ | Krivka Optim; | ||
+ | for(int i=0; i < n; i++) | ||
+ | { | ||
+ | y1=((B.Lin.a*body[i].x)+B.Lin.b)-body[i].y; | ||
+ | if(y1<0) y1=-y1; | ||
+ | sumLin+=y1; | ||
+ | |||
+ | y2=((B.Log.a*(log(body[i].x)))+B.Log.b)-body[i].y; | ||
+ | if(y2<0) y2=-y2; | ||
+ | sumLog+=y2; | ||
+ | |||
+ | y3=(B.Exp.b*(exp(B.Exp.a*body[i].x))); | ||
+ | y3=log(y3); | ||
+ | if(y3<0) y3=-y3; | ||
+ | sumExp+=y3; | ||
+ | |||
+ | y4=(B.Moc.b*(pow(body[i].x, B.Moc.a)))-body[i].y; | ||
+ | if(y4<0) y4=-y4; | ||
+ | sumMoc+=y4; | ||
+ | |||
+ | } | ||
+ | Naj=min(sumLin,sumLog); | ||
+ | Naj=min(Naj,sumExp); | ||
+ | Naj=min(Naj,sumMoc); | ||
+ | |||
+ | if(Naj==sumLin) {B.Lin.typ=LIN; return B.Lin;} | ||
+ | if(Naj==sumLog) {B.Log.typ=LOG; return B.Log;} | ||
+ | if(Naj==sumExp) {B.Exp.typ=EXP; return B.Exp;} | ||
+ | |||
+ | B.Moc.typ=MOC; | ||
+ | return B.Moc; | ||
+ | } | ||
+ | </source> | ||
+ | |||
+ | ===Modifikácia zdrojového kódu=== | ||
+ | Pre názornosť nebudeme v programe dáta načítavať z klávesnice, ale vstupné údaje budú zadané priamo v zdrojovom kóde | ||
+ | <source lang="c" line> | ||
+ | #include <iostream.h> | ||
+ | #include <conio.h> | ||
+ | #include <math.h> | ||
+ | #include <cmath> | ||
+ | |||
+ | struct Bod{ | ||
+ | double x,y; | ||
+ | }; | ||
+ | enum TYP_KRIVKY{ | ||
+ | LIN, LOG, EXP, MOC | ||
+ | }; | ||
+ | struct Krivka{ | ||
+ | TYP_KRIVKY typ; | ||
+ | double a,b; | ||
+ | }; | ||
+ | struct Aprox{ | ||
+ | Krivka Lin, Log, Exp, Moc; | ||
+ | }; | ||
+ | |||
+ | //Aprox aproximacie(Bod *body, int n, Aprox B) | ||
+ | //void vypisRovnicuKrivky(Krivka k) | ||
+ | //Krivka zistiOptimum(Bod *body, int n, Aprox B) | ||
+ | |||
+ | int main(int argc, char* argv[]) | ||
+ | { | ||
+ | const int n=8; | ||
+ | Bod body[n]={ | ||
+ | {1.50,0.58}, | ||
+ | {2.00,8.61}, | ||
+ | {3.00,60.91}, | ||
+ | {4.00,201.40}, | ||
+ | {5.00,476.59}, | ||
+ | {6.00,892.57}, | ||
+ | {7.00,1304.37}, | ||
+ | {8.00,2715.04} | ||
+ | }; | ||
+ | Aprox B; | ||
+ | B = aproximacie(body, n, B); | ||
+ | cout<<"Najvhodnejsia: "; | ||
+ | vypisRovnicuKrivky(zistiOptimum(body, n, B)); | ||
+ | |||
+ | getch(); | ||
+ | return 0; | ||
+ | } | ||
+ | </source> | ||
+ | |||
+ | ===Vzorový Príklad=== | ||
+ | Program vypočítal, že pre uvedené vstupné body je navhodnejšia Mocninová aproximácia. | ||
+ | {| class=datatable | ||
+ | |- | ||
+ | !<math>x_i</math> | ||
+ | !<math>y_i</math> | ||
+ | |- | ||
+ | |1.50 | ||
+ | |0.58 | ||
+ | |- | ||
+ | |2.00 | ||
+ | |8.61 | ||
+ | |- | ||
+ | |3.00 | ||
+ | |60.91 | ||
+ | |- | ||
+ | |4.00 | ||
+ | |201.40 | ||
+ | |- | ||
+ | |5.00 | ||
+ | |476.59 | ||
+ | |- | ||
+ | |6.00 | ||
+ | |892.57 | ||
+ | |- | ||
+ | |7.00 | ||
+ | |1304.37 | ||
+ | |- | ||
+ | |8.00 | ||
+ | |2715.04 | ||
+ | |} | ||
+ | |||
+ | [[Súbor:Vzor_graf.jpg|center|framed]] |
Aktuálna revízia z 22:31, 16. august 2010
Obsah
Zadanie 1
Riešte problém aproximácie dát. K dispozícii máme n bodov v rovine (ich súradnice x a y). Úlohou bude vypočítať rovnicu aproximujúcej krivky metódou najmenších štvorcov pre:
- Lineárnu aproximáciu v tvare [math]y=ax+b[/math]
- Logaritmickú aproximáciu v tvare [math]y=a \ln{x} + b[/math]
- Exponenciálnu aproximáciu v tvare [math]y=b e^{ax}[/math]
- Mocninovú aproximáciu v tvare [math]y=bx^a[/math]
Porovnajte výsledok funkcie s výsledkom v tabuľkovom procesore, kde do (bodového) grafu zo vstupných bodov pridáte spojnicu trendu typu lineárna, resp. logaritmická, mocninová exponenciálna.
Vstupné údaje
Ako vzorku vstupných bodov budeme uvažovať body z intervalu <1,8>. Pre generovanie vstupných hodnôt použijem funkciu [math]f(x)=x^4-exp(x)[/math]
Vzorový vstup
[math]x_i[/math] | [math]y_i[/math] |
---|---|
1.50 | 0.58 |
2.00 | 8.61 |
3.00 | 60.91 |
4.00 | 201.40 |
5.00 | 476.59 |
6.00 | 892.57 |
7.00 | 1304.37 |
8.00 | 1115.04 |
Vzorový výstup
- y = 213.634 x - 467.195
- y = 784.982 ln x - 572.83
- y = 0.9628 e^(1.054x)
- y = 0.2609 x^(4.443)
Analýza matematických vzťahov
Poznámka: v tomto texte nie sú uvádzané celé postupy výpočtu, nakoľko v tomto predmete nie sú primárnym cieľom matematické dôkazy. Študent si postup môže ľahko overiť sám.
Lineárna aproximácia: y=ax+b
Pre lineárnu aproximáciu platia vzťahy:
- [math]\begin{align} & b\sum\limits_{i=1}^{n}{{{x}_{i}}}+a\sum\limits_{i=1}^{n}{{{x}^{2}}_{i}}=\sum\limits_{i=1}^{n}{{{x}_{i}}{{y}_{i}}} \\ & nb+a\sum\limits_{i=1}^{n}{{{x}_{i}}}=\sum\limits_{i=1}^{n}{{{y}_{i}}} \\ \end{align}[/math]
Čo môžeme pomocou maticového zápisu zapísať nasledovne:
- [math]\left[ \begin{matrix} \sum\limits_{i=1}^{n}{{{x}^{2}}_{i}} & \sum\limits_{i=1}^{n}{{{x}_{i}}} \\ \sum\limits_{i=1}^{n}{{{x}_{i}}} & n \\ \end{matrix} \right]\cdot \left[ \begin{matrix} a \\ b \\ \end{matrix} \right]=\left[ \begin{matrix} \sum\limits_{i=1}^{n}{{{x}_{i}}{{y}_{i}}} \\ \sum\limits_{i=1}^{n}{{{y}_{i}}} \\ \end{matrix} \right][/math]
Riešenie tejto sústavy rovníc je:
- [math]\begin{align} & a=\frac{n\sum\limits_{i=1}^{n}{{{x}_{i}}{{y}_{i}}}-\sum\limits_{i=1}^{n}{{{x}_{i}}}\sum\limits_{i=1}^{n}{{{y}_{i}}}}{n\sum\limits_{i=1}^{n}{x_{i}^{2}}-{{\left( \sum\limits_{i=1}^{n}{{{x}_{i}}} \right)}^{2}}} \\ & b=\frac{\sum\limits_{i=1}^{n}{{{y}_{i}}}-a\sum\limits_{i=1}^{n}{{{x}_{i}}}}{n} \\ \end{align}[/math]
Logaritmická aproximácia: [math]y=a*\ln x + b[/math]
Hľadáme najlepšiu aproximáciu pre logaritmickú rovnicu: [math]y=a*\ln x + b[/math]. Postup výpočtu je rovnaký ako v prvom prípade. Jediná zmena v riešení je tá, že os x je v logaritmickej mierke. Dostávame teda sústavu 2 rovníc:
- [math]\left[ \begin{matrix} \sum\limits_{i=1}^{n}{{{\left( \ln {{x}_{i}} \right)}^{2}}} & \sum\limits_{i=1}^{n}{\left( \ln {{x}_{i}} \right)} \\ \sum\limits_{i=1}^{n}{\left( \ln {{x}_{i}} \right)} & n \\ \end{matrix} \right]\cdot \left[ \begin{matrix} a \\ b \\ \end{matrix} \right]=\left[ \begin{matrix} \sum\limits_{i=1}^{n}{\left( \ln {{x}_{i}} \right){{y}_{i}}} \\ \sum\limits_{i=1}^{n}{{{y}_{i}}} \\ \end{matrix} \right][/math]
Riešenie tejto sústavy rovníc je:
- [math]\begin{align} & a=\frac{n\sum\limits_{i=1}^{n}{\left( \ln {{x}_{i}} \right){{y}_{i}}}-\sum\limits_{i=1}^{n}{\ln {{x}_{i}}}\sum\limits_{i=1}^{n}{{{y}_{i}}}}{n\sum\limits_{i=1}^{n}{{{\left( \ln {{x}_{i}} \right)}^{2}}}-{{\left( \sum\limits_{i=1}^{n}{\ln {{x}_{i}}} \right)}^{2}}} \\ & b=\frac{\sum\limits_{i=1}^{n}{{{y}_{i}}}-a\sum\limits_{i=1}^{n}{\ln {{x}_{i}}}}{n} \\ \end{align}[/math]
Exponenciálna aproximácia: [math]y=b{{e}^{ax}}[/math]
Rovnicu [math]y=b{{e}^{ax}}[/math] si upravíme tak, že ju zapíšeme v logaritmickom tvare:
- [math]\begin{align} & \ln y=\ln \left( b{{e}^{ax}} \right) \\ & \ln y=\ln b+\ln {{e}^{ax}} \\ & \ln y=ax+\ln b\,\,\left[ B=\ln b \right] \\ & \ln y=ax+B \\ \end{align}[/math]
Opäť dostávame sústavu 2 rovníc o 2 neznámych. V rovniciach je použitá substitúcia [math]B=\ln b[/math].
- [math]\left[ \begin{matrix} \sum\limits_{i=1}^{n}{{{x}^{2}}_{i}} & \sum\limits_{i=1}^{n}{{{x}_{i}}} \\ \sum\limits_{i=1}^{n}{{{x}_{i}}} & n \\ \end{matrix} \right]\cdot \left[ \begin{matrix} a \\ B \\ \end{matrix} \right]=\left[ \begin{matrix} \sum\limits_{i=1}^{n}{{{x}_{i}}\ln {{y}_{i}}} \\ \sum\limits_{i=1}^{n}{\ln {{y}_{i}}} \\ \end{matrix} \right][/math]
Riešenie tejto sústavy rovníc je:
- [math]\begin{align} & a=\frac{n\sum\limits_{i=1}^{n}{{{x}_{i}}\ln {{y}_{i}}}-\sum\limits_{i=1}^{n}{{{x}_{i}}}\sum\limits_{i=1}^{n}{\ln {{y}_{i}}}}{n\sum\limits_{i=1}^{n}{x_{i}^{2}}-{{\left( \sum\limits_{i=1}^{n}{{{x}_{i}}} \right)}^{2}}} \\ & B=\frac{\sum\limits_{i=1}^{n}{\ln {{y}_{i}}}-a\sum\limits_{i=1}^{n}{{{x}_{i}}}}{n} \\ \end{align}[/math]
kde [math]B=\ln b[/math].
Mocnimová aproximácia: [math]y=b{{x}^{a}}[/math]
Rovnicu [math]y=b{{x}^{a}}[/math] si upravíme tak, že ju zapíšeme v logaritmickom tvare:
- [math]\begin{align} & \ln y=\ln \left( b{{x}^{a}} \right) \\ & \ln y=\ln b+\ln {{x}^{a}} \\ & \ln y=\ln b+a\ln x\,\,\left[ B=\ln b \right] \\ & \ln y=B+a\ln x \\ \end{align}[/math]
Dostávame sústavu 2 rovníc o 2 neznámych. V rovniciach je použitá substitúcia [math]B=\ln b[/math].
- [math]\left[ \begin{matrix} \sum\limits_{i=1}^{n}{\ln {{x}_{i}}^{2}} & \sum\limits_{i=1}^{n}{\ln {{x}_{i}}} \\ \sum\limits_{i=1}^{n}{\ln {{x}_{i}}} & n \\ \end{matrix} \right]\cdot \left[ \begin{matrix} a \\ B \\ \end{matrix} \right]=\left[ \begin{matrix} \sum\limits_{i=1}^{n}{\ln {{x}_{i}}\ln {{y}_{i}}} \\ \sum\limits_{i=1}^{n}{\ln {{y}_{i}}} \\ \end{matrix} \right][/math]
Riešením je:
- [math]\begin{align} & a=\frac{n\sum\limits_{i=1}^{n}{\ln {{x}_{i}}\ln {{y}_{i}}}-\sum\limits_{i=1}^{n}{\ln {{x}_{i}}}\sum\limits_{i=1}^{n}{\ln {{y}_{i}}}}{n\sum\limits_{i=1}^{n}{{{\left( \ln {{x}_{i}} \right)}^{2}}}-{{\left( \sum\limits_{i=1}^{n}{\ln {{x}_{i}}} \right)}^{2}}} \\ & B=\frac{\sum\limits_{i=1}^{n}{\ln {{y}_{i}}}-a\sum\limits_{i=1}^{n}{\ln {{x}_{i}}}}{n} \\ \end{align}[/math]
kde [math]B=\ln b[/math].
Analýza programátorského riešenia
Návrh dátových štruktúr
V úlohe budeme pracovať s bodmi v rovine, ktoré sú charakterizované x a y súradnicou. Preto si vytvoríme štruktúru Bod, ktorá bude takýto bod reprezentovať:
struct Bod{
double x,y;
};
V zadaní je úlohou vypočítať rovnicu najlepšie aproximujúcej krivky (lineárnej, logaritmickej, exponenciálnej alebo mocninovej). Tieto krivky sú charakterizované (okrem svojej matematickej funkcie, ktorú vyjadríme neskôr) koeficientami a a b. Na to, aby sme rozlíšili o akú krivku ide, vytvorme si vymenovaný zoznam TYP_KRIVKY, ktorý bude ma 4 symbolické hodnoty:
enum TYP_KRIVKY{
LIN, LOG, EXP, MOC
};
Význam týchto skratiek si zadefinujme nasledovne:
- LIN - lineárna krivka (polynóm 1 stupňa) s rovnicou [math]y=ax+b[/math]
- LOG - logaritmická krivka s rovnicou [math]y=a \ln{x} + b[/math]
- EXP - exponenciálna krivka s rovnicou [math]y=b e^{ax}[/math]
- MOC - mocninová krivka s rovnicou [math]y=bx^a[/math]
Teraz si môžeme zadefinovať štruktúru Krivka ako trojicu: typ (TYP_KRIVKY), koeficient a a koeficient b.
struct Krivka{
TYP_KRIVKY typ;
double a,b;
};
Riešenie v jazyku C
Vytvoríme funkciu Krivka minSq(Bod *body, int n) kde v poli body sú vstupné údaje, n je počet bodov. Funkcia bude vracať dátovú struktúru typu Krivka, kde budú informácie o tvare aproximujúcej krivky.
1 Krivka minSq_lin(Bod *body, int n)
2 { // y=a+bx
3 double sx=0,sxx=0,sy=0,sxy=0;
4 for(int i=0;i<n;i++)
5 { sx+=body[i].x;
6 sy+=body[i].y;
7 sxy+=body[i].x * body[i].y;
8 sxx+=body[i].x * body[i].x;
9 }
10 Krivka k;
11 k.a=(n*sxy-sy*sx)/(n*sxx-sx*sx);
12 k.b=(sy-k.a*sx)/n;
13 k.typ=LIN;
14 return k;
15 }
Podobným spôsobom môžeme vytvoriť ďalšie 3 funkcie minSq_log, minSq_exp minSq_moc. Ukážeme si ale všeobecnejší tvar funkcie minSq, kde pridáme ďalší parameter, ktorý bude hovoriť aký typ aproximácie sa má vypočítať.
1 Krivka minSq(Bod *body, int n, TYP_KRIVKY typ_proximacie)
2 {
3 double sx=0,sxx=0,sy=0,sxy=0,slxly;
4 double slx=0,slxx=0,sly=0,slxy=0,sxly=0;
5 for(int i=0;i<n;i++)
6 { sx+=body[i].x;
7 slx+=log(body[i].x);
8 sy+=body[i].y;
9 sly+=log(body[i].y);
10 sxy+=body[i].x * body[i].y;
11 slxy+=log(body[i].x) * body[i].y;
12 sxly+=body[i].x * log(body[i].y);
13 slxly+=log(body[i].x) *log(body[i].y);
14 sxx+=body[i].x * body[i].x;
15 slxx+=log(body[i].x) * log(body[i].x);
16 }
17 Krivka k;
18 switch(typ_proximacie)
19 {
20 case LIN:
21 k.a=(n*sxy-sy*sx)/(n*sxx-sx*sx);
22 k.b=(sy-k.a*sx)/n;
23 break;
24 case LOG:
25 k.a=(n*slxy-sy*slx)/(n*slxx-slx*slx);
26 k.b=(sy-k.a*slx)/n;
27 break;
28 case EXP:
29 k.a=(n*sxly-sly*sx)/(n*sxx-sx*sx);
30 k.b=(sly-k.a*sx)/n;
31 k.b=exp(k.b);
32 break;
33 case MOC:
34 k.a=(n*slxly-sly*slx)/(n*slxx-slx*slx);
35 k.b=(sly-k.a*slx)/n;
36 k.b=exp(k.b);
37 break;
38 }
39 k.typ=typ_proximacie;
40 return k;
41 }
Použitie v programe
Pre názornosť nebudeme v programe dáta načítavať z klávesnice, ale vstupné údaje budú zadané priamo v zdrojovom kóde. Modifikácia zdrojového kódu s načítaním vstupných hodnôt je jednoduchá a určite ju každý zvládne sám.
V zdrojovom kóde pribudla ešte funkcia vypisRovnicuKrivky, ktorá vypíše rovnicu krivky v čitateľnom tvare pre človeka. Vo výpise sú vynechané všetky funkcie, ktoré sú boli v tomto príklade vysvetlené a ich zdrojový kód bol napísaný skôr.
1 #include <iostream.h>
2 #include <conio.h>
3 #include <math.h>
4
5 struct Bod {double x,y;};
6 enum TYP_KRIVKY{ LIN, LOG, EXP, MOC};
7 struct Krivka{
8 TYP_KRIVKY typ;
9 double a,b;
10 };
11
12 //Krivka minSq(Bod *body, int n, TYP_KRIVKY typ_proximacie)
13
14 void vypisRovnicuKrivky(Krivka k)
15 {
16 cout<<"y=";
17 switch(k.typ)
18 {
19 case LIN:
20 cout<<k.a<<" x";
21 if(k.b>0) cout<<"+";
22 cout<<k.b;
23 break;
24 case LOG:
25 cout<<k.a<<" ln x";
26 if(k.b>0) cout<<"+";
27 cout<<k.b;
28 break;
29 case EXP:
30 cout<<k.b<<" e^(";
31 cout<<k.a<<" x)";
32 break;
33 case MOC:
34 cout<<k.b<<" x^(";
35 cout<<k.a<<" )";
36 break;
37 }
38 }
39 int main()
40 { const int n=8;
41 Bod body[n]={
42 {1.50,0.58},
43 {2.00,8.61},
44 {3.00,60.91},
45 {4.00,201.40},
46 {5.00,476.59},
47 {6.00,892.57},
48 {7.00,1304.37},
49 {8.00,1115.04}
50 };
51 Krivka Q;
52 Q=minSq(body,n,EXP); // LOG, LIN, EXP
53 vypisRovnicuKrivky(Q);
54 getch();
55 }
Overenie vypočítaných hodnôt
Na nasledujúcich obrázkoch sú vstupné body z príkladu znázornené v tabuľkovom procesore a je vypočítaná najlepšia aproximujúca krivka (MS Excel ju nazýva spojnica trendu).
Zadanie 2
Vytvorte funkciu, ktorá pre zadaná vstupy vypočíta všetky typy aproximácií (lineárnu, logaritmickú, exponenciálnu a mocninovú) a určí ktorá aproximácia je optimálna pre zadané vstupné body.
Analýza
Pri vyhodnocovaní optimálnej aproximácie musíme vypočítať všetky možné aproximačné krivky a zvoliť kritérium voľby najlepšej aproximácie. Kritériom posudzovania môže byť hľadanie minimálych hodnôt:
- vzdialenosť uzlových (vstupných) bodov od aproximujúcej krivky
- štvorce vzdialeností uzlových (vstupných) bodov od aproximujúcej krivky
- plocha medziproximujúcou krivkou a interpolačným polynómom, ktorých prechádza cez uzlové body
Riešenie v jazyku C
Kriteriom na vyhodnotenie najlepšej aproximovanej krivky bude hľadanie minimálych hodnôt a vzdialeností uzlových (vstupných) bodov od aproximujúcej krivky.
Definovanie štruktúr potrebných pre spracovanie súradníc:
1 struct Bod{
2 double x,y;
3 };
4
5 enum TYP_KRIVKY{
6 LIN, LOG, EXP, MOC
7 };
8
9 struct Krivka{
10 TYP_KRIVKY typ;
11 double a,b;
12 };
Táto štruktúra pre zjednodušenie poslúži na uloženie všetkých typov aproximácie do jednej premennej typu Aprox:
1 struct Aprox{
2 Krivka Lin, Log, Exp, Moc;
3 };
Funkcia aproximacie vypočíta všetky 4 typy aproximácií. Návratová hodnota funkcie je typu Aprox:
1 Aprox aproximacie(Bod *body, int n, Aprox B) //vypocita vsetky
2 //aproximacie
3 {
4 double sx=0,sxx=0,sy=0,sxy=0,slxly;
5 double slx=0,slxx=0,sly=0,slxy=0,sxly=0;
6 for(int i=0;i<n;i++)
7 { sx+=body[i].x;
8 slx+=log(body[i].x);
9 sy+=body[i].y;
10 sly+=log(body[i].y);
11 sxy+=body[i].x * body[i].y;
12 slxy+=log(body[i].x) * body[i].y;
13 sxly+=body[i].x * log(body[i].y);
14 slxly+=log(body[i].x) *log(body[i].y);
15 sxx+=body[i].x * body[i].x;
16 slxx+=log(body[i].x) * log(body[i].x);
17 }
18
19 B.Lin.a=(n*sxy-sy*sx)/(n*sxx-sx*sx);
20 B.Lin.b=(sy-B.Lin.a*sx)/n;
21
22 B.Log.a=(n*slxy-sy*slx)/(n*slxx-slx*slx);
23 B.Log.b=(sy-B.Log.a*slx)/n;
24
25 B.Exp.a=(n*sxly-sly*sx)/(n*sxx-sx*sx);
26 B.Exp.b=(sly-B.Exp.a*sx)/n;
27 B.Exp.b=exp(B.Exp.b);
28
29 B.Moc.a=(n*slxly-sly*slx)/(n*slxx-slx*slx);
30 B.Moc.b=(sly-B.Moc.a*slx)/n;
31 B.Moc.b=exp(B.Moc.b);
32
33 return B;
34 }
Na vypísanie najvhodnejšieho typu rovnice je použitá funkcia vypisRovnicuKrivky:
1 void vypisRovnicuKrivky(Krivka k) //vypise rovnicu urciteho typu
2 {
3 cout<<"y = ";
4 switch(k.typ)
5 {
6 case LIN:
7 cout<<k.a<<" x";
8 if(k.b>0) cout<<"+";
9 cout<<k.b<<" -> LINEARNA"<<endl;
10 break;
11 case LOG:
12 cout<<k.a<<" ln x";
13 if(k.b>0) cout<<"+";
14 cout<<k.b<<" -> LOGARITMICKA"<<endl;
15 break;
16 case EXP:
17 cout<<k.b<<" e^(";
18 cout<<k.a<<" x)"<<" -> EXPONENCIALNA"<<endl;
19 break;
20 case MOC:
21 cout<<k.b<<" x^(";
22 cout<<k.a<<" )"<<" -> MOCNINOVA"<<endl;
23 break;
24 }
25 }
Pre zistenie optimálnej aproximovanej krivky pre dané vstupné body je použitý jednoduchý algoritmus daný vo funkcii zistiOptimum.
- Postup riešenia algoritmu spočíva v zistení vzdialenosti y-onovej súradnice vstupného bodu od aproximovanej y-onovej súradnice, ktorú musíme najskôr vypočítať, a porovnaním ako sú tieto body vzdialené od jednotlivých bodov kriviky pre každý typ, zistíme optimum.
- Dosadením x-ových súradníc vstupných bodov do príslušnej rovnice krivky sa vypočíta aproximovaná y-onová súradnica. Absolútnu hodnotu rozdielu vzdialenosti medzi y-onovými súradnicami(vypočítanými a vstupnými) sčíta a uloží túto sumu do premennej sumLin (pre každý ďalší jeden typ aproximácie do - sumLog, sumExp, sumMoc ). Toto sa opakuje pre každý jeden vstupný bod a pre každú rovnicu krivky.
- Na záver sa sumy všetkých typov aproximácií porovnajú. Najmenšia suma predstavuje najmenšie odchýlky pre danú krivku v porovnaní so vstupnými bodmi. Funkcia vráti Optimálnu krivku aproximácie s najmenšou sumou:
1 Krivka zistiOptimum(Bod *body, int n, Aprox B)
2 {
3 double sumLin=0, sumLog=0, sumExp=0, sumMoc=0, y1=0,y2=0,y3=0,y4=0, Naj=0;
4 Krivka Optim;
5 for(int i=0; i < n; i++)
6 {
7 y1=((B.Lin.a*body[i].x)+B.Lin.b)-body[i].y;
8 if(y1<0) y1=-y1;
9 sumLin+=y1;
10
11 y2=((B.Log.a*(log(body[i].x)))+B.Log.b)-body[i].y;
12 if(y2<0) y2=-y2;
13 sumLog+=y2;
14
15 y3=(B.Exp.b*(exp(B.Exp.a*body[i].x)));
16 y3=log(y3);
17 if(y3<0) y3=-y3;
18 sumExp+=y3;
19
20 y4=(B.Moc.b*(pow(body[i].x, B.Moc.a)))-body[i].y;
21 if(y4<0) y4=-y4;
22 sumMoc+=y4;
23
24 }
25 Naj=min(sumLin,sumLog);
26 Naj=min(Naj,sumExp);
27 Naj=min(Naj,sumMoc);
28
29 if(Naj==sumLin) {B.Lin.typ=LIN; return B.Lin;}
30 if(Naj==sumLog) {B.Log.typ=LOG; return B.Log;}
31 if(Naj==sumExp) {B.Exp.typ=EXP; return B.Exp;}
32
33 B.Moc.typ=MOC;
34 return B.Moc;
35 }
Modifikácia zdrojového kódu
Pre názornosť nebudeme v programe dáta načítavať z klávesnice, ale vstupné údaje budú zadané priamo v zdrojovom kóde
1 #include <iostream.h>
2 #include <conio.h>
3 #include <math.h>
4 #include <cmath>
5
6 struct Bod{
7 double x,y;
8 };
9 enum TYP_KRIVKY{
10 LIN, LOG, EXP, MOC
11 };
12 struct Krivka{
13 TYP_KRIVKY typ;
14 double a,b;
15 };
16 struct Aprox{
17 Krivka Lin, Log, Exp, Moc;
18 };
19
20 //Aprox aproximacie(Bod *body, int n, Aprox B)
21 //void vypisRovnicuKrivky(Krivka k)
22 //Krivka zistiOptimum(Bod *body, int n, Aprox B)
23
24 int main(int argc, char* argv[])
25 {
26 const int n=8;
27 Bod body[n]={
28 {1.50,0.58},
29 {2.00,8.61},
30 {3.00,60.91},
31 {4.00,201.40},
32 {5.00,476.59},
33 {6.00,892.57},
34 {7.00,1304.37},
35 {8.00,2715.04}
36 };
37 Aprox B;
38 B = aproximacie(body, n, B);
39 cout<<"Najvhodnejsia: ";
40 vypisRovnicuKrivky(zistiOptimum(body, n, B));
41
42 getch();
43 return 0;
44 }
Vzorový Príklad
Program vypočítal, že pre uvedené vstupné body je navhodnejšia Mocninová aproximácia.
[math]x_i[/math] | [math]y_i[/math] |
---|---|
1.50 | 0.58 |
2.00 | 8.61 |
3.00 | 60.91 |
4.00 | 201.40 |
5.00 | 476.59 |
6.00 | 892.57 |
7.00 | 1304.37 |
8.00 | 2715.04 |