Algoritmy numerickej aproximácie (riešené príklady): Rozdiel medzi revíziami

Z Kiwiki
Skočit na navigaci Skočit na vyhledávání
 
(7 medziľahlých úprav od 2 ďalších používateľov nie je zobrazených)
Riadok 1: Riadok 1:
[[Kategória:Študijné materiály]]
 
[[Kategória:Programovanie]]
 
[[Kategória:jazyk C]]
 
 
{{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 = 6.9628 e^(1.054x)
+
#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

Imbox draft.png
Toto je projekt, na ktorom sa ešte stále pracuje!!

Aj keď sú v tomto dokumente použiteľné informácie, ešte nie je dokončený. Svoje návrhy môžete vyjadriť v diskusii o tejto stránke.

Algoritmy a programovanie - zbierka úloh


Štruktúry

Rekurzia

Dynamická alokácia pamäti

Vyhľadávanie

Triedenie

Lineárny zoznam

Binárny strom

Numerické algoritmy

>Algoritmy numerickej interpolácie
>Algoritmy numerickej aproximácie
>Numerické integrovanie
>Numerické derivovanie

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:

  1. Lineárnu aproximáciu v tvare [math]y=ax+b[/math]
  2. Logaritmickú aproximáciu v tvare [math]y=a \ln{x} + b[/math]
  3. Exponenciálnu aproximáciu v tvare [math]y=b e^{ax}[/math]
  4. 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

  1. y = 213.634 x - 467.195
  2. y = 784.982 ln x - 572.83
  3. y = 0.9628 e^(1.054x)
  4. 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
Vzor graf.jpg