1. Wstęp

Celem niniejszego ćwiczenia było napisanie algorytmu aproksymacji funkcji y = sin(2x). Program powinien czytać dane wejściowe z pliku tekstowego, jak również eksportować do takiego pliku dane wynikowe. Na koniec obliczeń program ma określić błąd aproksymacji.

2. Wprowadzenie teoretyczne

Aproksymacja polega na przybliżeniu funkcji f(x) (funkcji aproksymowanej) inną funkcją Q(x) (funkcją aproksymującą). Funkcja aproksymująca jest kombinacją liniową tzw. funkcji bazowych. Jako funkcje bazowe stosuje się: jednomiany, funkcje trygonometryczne, oraz tzw wielomiany ortogonalne. W naszym programie będziemy używać jednomianów.

Do aproksymacji średniokwadratowej, punktowej wielomianowej używamy następującej postaci funkcji aproksymującej:

Wielomian aproksymujący musi być tak dobrany aby błąd średniokwadratowy określony wzorem:

był jak najmniejszy.

W celu uzyskania współczynników
należy rozwiązać następujący układ równań liniowych:

Gdzie parametry
oraz
definiowane są w następujący sposób:

Parametry
są elementami macierzy głównej rozpatrywanego układu równań liniowych, a parametry

3. Opis programu.

Program wymaga podania przez użytkownika następujących danych:

• parametr n (ilość punktów używanych w aproksymacji

• parametr m (stopień wielomianu aproksymującego).

Program odczytuje również z pliku tekstowego, współrzędne x punktów aproksymacji. Ścieżka tego pliku określona jest na sztywno jako „C:\x.txt”. Plik ten musi zostać przygotowany przed rozpoczęciem programu.

Po odczycie danych wejściowych, program inicjalizuje wszystkie tablice, które są wykorzystywane w programie.

Następnie konstruujemy w programie główny układ równań liniowych który służy nam do znalezienia współczynników wielomianu aproksymującego Q(x).

Układ ten rozwiązywany jest metodą eliminacji Gaussa, zaimplementowaną w dwóch funkcjach:

ZrobMacierzTrojkatna() - transformuje macierz główna układu równań do postaci górnotrójkątnej, dokonując jednocześnie zmian w kolejności elementów wektorów B oraz X zgodnie ze zmianaw kolejności wierszy i kolumn macierzy głównej układu.

RozwiazMacierzTrojkatna() - dokonuje rozwiązania, układu z macierzą górnotrójkątną i odtwarza pierwotną kolejność elementów wektora rozwiązań.

Po otrzymaniu jako wynik współczynników wielomianu aproksymującego, dokonujemy obliczeń wartości wielomianu aproksymującego w punktach aproksymacji i zapisujemy je wraz z oryginalnymi wartościami funkcji w pliku tekstowym „C:\wyniki.txt”.

Na koniec program oblicza błąd aproksymacji ze wzoru:

wzór z instrukcji

4. Testowanie i wnioski

Program został przetestowany na dwóch zestawach punktów liczących odpowiednio 5 oraz 31 elementów. Każdy zestaw z osobna był testowany dla wielomianów aproksymujących stopnia 5, 10 oraz 15. Wyniki testów przedstawione są poniżej.

Dla 5 punktów:

m = 5

Δ = 0.000000 Błąd jest większy od zera, tylko precyzja wyświetlania wyniku jest zbyt niska.

(zielonym kolorem przedstawiono przebieg funkcji aproksymującej)

m = 10

Δ = 0.000000 Błąd jest większy od zera, tylko precyzja wyświetlania wyniku jest zbyt niska.

(zielonym kolorem przedstawiono przebieg funkcji aproksymującej)

m = 15

Δ = 0.000000 Błąd jest większy od zera, tylko precyzja wyświetlania wyniku jest zbyt niska.

(zielonym kolorem przedstawiono przebieg funkcji aproksymującej)

Dla 31 punktów:

m = 5

Δ =
0.004327

(zielonym kolorem przedstawiono przebieg funkcji aproksymującej)

m = 10

Δ = 0.000002

(zielonym kolorem przedstawiono przebieg funkcji aproksymującej)

m = 15

Δ = 0.000000 Błąd jest większy od zera, tylko precyzja wyświetlania wyniku jest zbyt niska.

(zielonym kolorem przedstawiono przebieg funkcji aproksymującej)

Jak widać z powyżej zamieszczonych wyników algorytm aproksymacji został poprawnie zaimplementowany. Bardzo dobre wyniki dla małej ilości punktów funkcji aproksymowanej można tłumaczyć tym że były one rozłożone bardzo blisko siebie tak że przebieg funkcji był nieomal liniowy, stad wielomian piątego rzędu był w zupełności wystarczający dla odtworzenia przebiegu tego fragmentu funkcji.

Przy dłuższym i bardziej zróżnicowanym odcinku funkcji, można obserwować zależność między stopniem wielomianu aproksymującego a dokładnością aproksymacji. Im wyższy stopień wielomianu aproksymującego tym dokładniejsza reprodukcja oryginalnego przebiegu funkcji.

5. Listing programu

#include "stdafx.h"

#include "iostream"

#include "stdio.h"

#include "math.h";

double blad;

double *a;

double *x;

double *wynik;

double *f;

int n,m;

double **A,*B,*X;

double pi;

int rozmiar;

int *kolejnosc;

double *tymczasowy_X;

void ReadXWector()

{

FILE *file;

file = fopen("c:\\x.txt","r");

for(int i = 0 ; i <=n; i++)

fscanf(file,"%lf",&x[i]);

}

int ZrobMacierzTrojkatna()

{

double max,element_redukujacy;

int max_x,max_y,index_tymczasowy;

double *wiersz_tymczasowy,*kolumna_tymczasowa,element_tymczasowy;

wiersz_tymczasowy = new double[rozmiar];

kolumna_tymczasowa = new double[rozmiar];

for(int i = 0; i <rozmiar; i++)

kolejnosc[i]=i;

for(int i = 0; i <rozmiar; i++)

{

//Szuka w podmacierzy elelementu o najwiekszym module

//oraz jego lokalizacji

max = fabs(A[i][i]);

max_x = i;

max_y = i;

for(int y = i;y<rozmiar;y++)

for(int x = i; x<rozmiar;x++)

{

if(fabs(A[y][x]) > max)

{

max = A[y][x];

max_x = x;

max_y = y;

}

}

//zamienia wiersz i-ty macierzy A z tym w którym jest element o najwiekszym module

for(int j =0;j < rozmiar;j++)

{

wiersz_tymczasowy[j] = A[i][j];

A[i][j] = A[max_y][j];

A[max_y][j] = wiersz_tymczasowy[j];

}

//zamienia kolumne i-ta macierzy A z tym w którym jest element o najwiekszym module

for(int j =0;j < rozmiar;j++)

{

kolumna_tymczasowa[j] = A[j][i];

A[j][i] = A[j][max_x];

A[j][max_x] = kolumna_tymczasowa[j];

}

//analogiczna zmiane jak zmiane kolejnosci wierszy wykonujemy w wektorze B

element_tymczasowy = B[i];

B[i] = B[max_y];

B[max_y] = element_tymczasowy;

//Zapisujemy ta zmiane kolejnosci kolumn by odtworzyc wlasciwa kolejnosc X

index_tymczasowy = kolejnosc[i];

kolejnosc[i] = kolejnosc[max_x];

kolejnosc[max_x] = index_tymczasowy;

//Sprawdzamy czy element aii = 0 ? jeżeli tak sygnalizujemy blad i konczymy

if(A[i][i] == 0)

return 0;

//Dokonujemy redukcji elementow i-tej kolumny w wierszach od i+1 do n wz 7

for(int y = i+1;y<rozmiar;y++)

{

element_redukujacy = A[y][i];

for(int x = i; x<rozmiar;x++)

{

A[y][x]=A[y][x] - element_redukujacy*A[i][x]/A[i][i];

}

B[y] = B[y] - element_redukujacy*B[i]/A[i][i];

}

}

return 1;

}

void RozwiazMacierzTrojkatna()

{

//wyliczanie poszczególnych X

for(int i = rozmiar - 1; i >= 0; i--)

{

X[i] = B[i]/A[i][i];

for(int j = i+1; j < rozmiar; j++)

{

X[i] = X[i]-A[i][j]*X[j]/A[i][i];

}

}

//porzadkowanie wektora X w celu przywrocenia pierwotnej kolejnosci

//najpierw zapis do wektora tymczasowego

for(int i = 0; i < rozmiar; i++)

tymczasowy_X[i]=X[i];

//nastepnie wpisywanie x we wlasciwej kolejnosci

for(int i = 0; i < rozmiar; i++)

X[kolejnosc[i]]=tymczasowy_X[i];

}

double Function(double x)

{

return sin(2*x);

}

int _tmain(int argc, _TCHAR* argv[])

{

double pi = 3.14159265358979323846;

double sum;

printf("Podaj n ");

scanf("%d",&n);

printf("Podaj m ");

scanf("%d",&m);

rozmiar = m+1;

kolejnosc = new int[rozmiar];

tymczasowy_X = new double[rozmiar];

A = new double*[rozmiar];//tworzenie tablicy wskaznikow

A[0] = new double[rozmiar*rozmiar];//tablica elementow double

for(int i=1;i<rozmiar;i++)

A[i] = A[i-1] + rozmiar;

B = new double[m+1];

X = new double[m+1];

wynik = new double[n+1];

a = new double[n+1];

x = new double[n+1];

f = new double[n+1];

ReadXWector();

for(int i = 0; i <=n;i++)

f[i] = Function(x[i]);

for(int i=0;i<=m;i++)

{

sum = 0;

for(int j = 0 ; j <=n; j++)

sum += f[j]*pow(x[j],i);

B[i] = sum;

}

for(int i = 0;i<=m;i++)

{

for(int j = 0; j <= m ;j++)

{

sum = 0;

for(int k = 0;k<=n;k++)

sum += pow(x[k],i+j);

A[i][j]= sum;

}

}

ZrobMacierzTrojkatna();

RozwiazMacierzTrojkatna();

for(int i = 0 ; i <= n;i++)

{

sum = 0;

for ( int j = 0 ; j <= m; j++)

sum += X[j]*pow(x[i],j);

wynik[i] = sum;

}

blad = 0 ;

for(int i = 0; i <=n;i++)

{

blad += pow(wynik[i]-f[i],2);

/*for(int j = 0 ; j <= m ; j++)

{

temp += wynik[i] - f[i];

}

blad += temp*temp;*/

}

blad = sqrt(blad/(n+1));

printf(" Blad aproksymacji delta = %lf",blad);

FILE *file;

file = fopen("c:\\wyniki.txt","w");

for(int i = 0 ; i <=n; i++)

fprintf(file,"%lf\t%lf\t%lf\n",x[i],f[i],wynik[i]);

fclose(file);

system("pause");

return 0;

}

Zabrze, dnia 03.03.2008 r.

Politechnika Śląska

Wydział Automatyki Elektroniki i Informatyki

WSZInf.

Laboratorium Metod Numerycznych.

Temat ćwiczenia: Aproksymacja funkcji y = sin(2x).

Skład sekcji:

Dyrek Adam

Buczkowski Arkadiusz

Kentnowski Marek

---