Aproksymacja Funkcji



Program aproksymuje funkcję x^3+3*x+3 metodą najmniejszych kwadratów w 100 losowo wybranych punktach z przedziału (0;10)
Po wybraniu punktów i obliczeniu wartości funkcji w tych punktach, program dodaje jeszcze do tych wartości losowe zakłócenie z przedziału (-0.1;0.1) (-0.5;0.5) (-1;1) (-2;2) (-10;10). Własnie te wartości z zakłóceniami aproksymuje i sprawdza jak różne zakłócenia wpływają na aproksymację.


#include <iostream>
#include <cstdlib>
#include <cmath>

using namespace std;

double fodx(double X);
double mojrand(double zakr);
double Ptga(double lcba, int dopot);
void aprox(double zakbl);

#define n 99
#define m 2

// m - stopien wielomianu
// n - ilosc pkt od zera

double S[2*m],T[m],x[n+1],p[n+1],y[n+1][n],yx[n+1],yapr[n+1],yg[n+1];
int rozm;
double rzad=1.0;
int main(int argc, char *argv[])
{
  // losowanie x-ow i obliczanie y-ow poprawnych
  long int a;
  for (int ia=0;ia<=n;ia++)
    {
    a = rand();
    x[ia] = double(a)/2147483647*10;
    yg[ia] = fodx(x[ia]);
    }
  aprox(0.1);
  aprox(0.5);
  aprox(1.0);
  aprox(2.0);
  aprox(10.0);
  return EXIT_SUCCESS;
}

void aprox(double zakbl)
{
rzad = 1.0;
  // losowanie bledu i odejmownie od y-ow
  for (int ih=0;ih<=n;ih++)
    yx[ih]=yg[ih]-mojrand(zakbl);

// obliczanie S[k]
for (int ib=0; ib<=2*m; ib++)
   {
   S[ib]=0;
   for (int ic=0; ic<=n; ic++)
        S[ib]+=Ptga(x[ic],ib);
   }

// obliczanie T[k]
for (int id=0; id<=m; id++)
   {
   T[id]=0;
   for (int ie=0; ie<=n; ie++)
        T[id]+=Ptga(x[ie],id)*yx[ie];
   }

// przepisanie y-ow i x-ow na macierz do obliczenia ukladu rownan
rozm=m+1;
for (int jj=0; jj<=m; jj++)
   for (int gg=jj,hh=0; gg<=jj+m; gg++,hh++)
      y[hh][jj]=S[gg];
for (int dd=0;dd<=m;dd++)
   y[rozm][dd]=T[dd];

// rozwiazywanie ukladu
 for (int b=0;b<rozm;b++)
  {
  //porzadkowanie macierzy
  for (int na=b;na<rozm;na++)
     if (y[b][na]!=0)
        for (int p=rozm;p>=b;p--)
           y[p][na]/=y[b][na];
  //odejmowanie wierszy (jak trzeba)
  for (int v=b+1;v<rozm;v++)
     if (y[b][v]!=0)
        for (int vv=b;vv<=rozm;vv++)
           y[vv][v]-=y[vv][b];
  }
for (int e=0;e<rozm;e++)
    rzad*=y[e][e];
 if (rzad==0)
   printf(" Rzad ukladu nie jest rowny liczbie niewiadomych, rownanie ma nieskonczenie wiele rozwiazan. Podaj inny uklad rownan.");
 else
   {
   //wszystko ok
   for (int ij=0;ij<=n;ij++)
     {
    yapr[ij]=0;
    // obliczanie y aproksymowanych
    for (int iv=0,ptg=rozm-1;iv<rozm;iv++,ptg--)
       {
       p[iv]=y[rozm][rozm-iv-1];
       for (int jv=0;jv<iv;jv++)
         p[iv]-=y[rozm-jv-1][rozm-iv-1]*p[jv];
       //printf("(%lf)*x^%d\n",p[iv],ptg);
       yapr[ij]+=p[iv]*Ptga(x[ij],ptg);
       }
     }
   }

  double S1=0;

  for (int iz=0;iz<=n;iz++)
    {
    S1+=Ptga((yg[iz]-yapr[iz]),2);
    //printf("x[i] = %lf y[x] = %lf y[z] = %lf yap = %lf \n",x[iz],yg[iz],yx[iz],yapr[iz]);
    //printf("%lf\n",yapr[iz]);
    }
  S1=sqrt(S1);
  //printf("\n");
  //printf("%lf\n",yapr[iz]);
  printf("dla bledu z zakresu %lf S = %lf\n",zakbl,S1);
}

double fodx(double X)
{
return 3*X*X*X+3*X+3;
}

double mojrand(double zakr)
{
double wyn;
long int a = random();
wyn = double(a)/2147483647;
wyn *= zakr*2;
wyn -= zakr;
return wyn;
}

double Ptga(double lcba, int dopot)
{
double wynik=1;
if (dopot != 0)
  {
  if (lcba == 0)
    wynik = 0;
  else
    for(int ki=1;ki<=abs(dopot);ki++)
      wynik*=lcba;
  }
if (dopot<0)
  wynik=1/wynik;
return wynik;
}

Powrot do Programów
Powrót na stronę główną