[b]### Treść zadania ###[/b] [cytat]Na ile podprzedziałów należy podzielić przedział 0 do 5 aby funkcja f(x)=1/(1+x^2) była interpolowana z dokładnością co najmniej 1e-5.

Normę błędu liczyć jako L2. 

Dla uproszczenia dokonywać podziału na 1,2,4,8,16,...,2^n podprzedziałów.

Interpolacja liniowa[/cytat]

[b]UWAGA:[/b] Program kompilowany w Dev-C++

[code=cpp]#include #include #include #include #include

using namespace std;

//=================================================================================

struct wspolrzedne { double x; double y; };

// Nasza funkcja double f(double x) { return 1/(1+x*x); }

// Wzor funkcji liniowej przechodzącej przez dwa punkty [x0, f(x0)] i [x1, f(x1)] double p(double x, double x0, double x1) { return ( f(x0) + ( (f(x1)-f(x0))/(x1-x0) ) * (x-x0) ); }

/* ### UWAGA: Funkcja pomocnicza. Nie jest wykorzystywana w programie.

Całkowanie numeryczne - metoda Simpsona

Metoda ma zastosowanie do funkcji stablicowanych w nieparzystej liczbie równo odległych punktów 
(wliczając końce przedziału całkowania). Metoda opiera się na przybliżaniu funkcji całkowanej 
przez interpolację wielomianem drugiego stopnia.

xp - x początkowe
xk - x końcowe
n - liczba punktów podziałowych dla przedziału od xp do xk

*/ double calka_simpsona_z_f(double xp, double xk, int n) { double x; double s = 0; double st = 0; double h = (xk - xp) / n; // odleglosc miedzy kolejnym wezlem

for(int i=1; i<=n; i++)
{
	x = xp+i*h;
	st += f(x-h/2);
	if(i<n) {
		s += f(x);
	}
}
return h / 6 * (f(xp) + f(xk) + 2 * s + 4 * st);	

}

/* ### UWAGA: Funkcja pomocnicza. Nie jest wykorzystywana w programie.

Całkowanie metodą Simpsona funkcji liniowej p przechodzącej przez dwa punkty x0 i x1

*/ double calka_simpsona_z_p(double xp, double xk, int n, double x0, double x1) { double x; double s = 0; double st = 0; double h = (xk - xp) / n; // odleglosc miedzy kolejnym wezlem

for(int i=1; i<=n; i++)
{
	x = xp+i*h;
	st += p(x-h/2, x0, x1);
	if(i<n) {
		s += p(x, x0, x1);
	}
}
return h / 6 * (p(xp, x0, x1) + p(xk, x0, x1) + 2 * s + 4 * st);	

}

/* Funkcja interpolująca. Jest złożona z wielu funkcji liniowych. x - argument dla ktorego liczymy wartosc tab - tablica wspolrzednych punktow wezlowych

*/ double finterp(double x, vector &tab) { double y = 0; // przechowuje obliczoną wartosc funkcji interpolujacej

for(int i=0; i<tab.size(); i++)
{
	// szukamy między jakimi punktami węzłowymi znajduje się podany jako argument x
	// a następnie obliczamy wartosc funkcji liniowej przechodzącej przez te punkty węzłowe dla podanego x
	if(x <= tab[i+1].x)
	{
		y = p(x, tab[i].x, tab[i+1].x);
		break;	// po znalezieniu odpowiedniego przedziału przerywamy pętle
	}
}
return y;

}

/* Funkcja g(x) = (f(x) - finterp(x))^2 Jest to funkcja pod całkowa do wzoru na norme euklidesową - L2 */ double g(double x, vector &tab) { double tmp = f(x)-finterp(x, tab); return pow(tmp, 2); }

/* Całka Simpsona z funkcji g Całka ta zostanie wykorzystana do liczenia normy L2

Zalecana liczba podprzedzialow n to iloczyn liczby podprzedziałów całej funkcji 
i ilości podziałów w danym podprzedziale.

*/ double calka_simpsona_z_g(double xp, double xk, int n, vector &tab) { double x; double s = 0; double st = 0; double h = (xk - xp) / n; // odleglosc miedzy kolejnym wezlem

for(int i=1; i<=n; i++)
{
	x = xp+i*h;
	st += g(x-h/2, tab); // suma wartosci dla wczesniejszych węzłów
	if(i<n) { // zabezpieczamy sie przed wyjsciem poza przedzial
		s += g(x, tab); // suma wartosci dla następnych węzłów
	}
}
return h / 6 * (g(xp, tab) + g(xk, tab) + 2 * s + 4 * st);	

}

//=================================================================================

int main(int argc, char *argv[]) { cout«"### (c) Webook.pl (2010/2011) - Interpolacja liniowa ###nn"; //cout.setf(ios::showpoint); // uzupelnia liczby zerami po przecinku cout.precision(30); // wyswietla liczby zmienno przecinkowe do 30 miejsc po przecinku

double eps = pow(10.,-5); // wymagana dokładność
	
double a = 0.; // początek przedziału
double b = 5.; // koniec przedziału

int licznik = 2;
double n = pow(2., licznik); // liczba podprzedzialow od ktorej startujemy 

double h; // odleglosc miedzy kolejnym wezlem na osi x
double L2; // zmienna przechowujaca norme euklidesowa

wspolrzedne	tmp_wsp; // zmienna pomocnicza do tworzenia struktury zawierajacej wspolrzedne punktu [x,y]

vector<wspolrzedne> tab; // tworzymy vector złożony ze struktury dwóch zmiennych x i y

// powtarzamy pętlę aż zostanie wywołany break
while(1) 
{
	tab.clear(); // kasujemy stan vectora przy rozpoczeciu kazdej kolejnej iteracji
	
	h = (b-a)/n; // liczymy odleglosc miedzy kolejnym wezlem na osi x
	
	cout<<"=== Krok "<<licznik-1<<" ================n"
		<<"n = "<<n<<" (liczba podprzedzialow)n"
		<<"h = "<<h<<" (odleglosc miedzy kolejnym wezlem)nn";
		
	cout<<"Wspolrzedne punktow wezlowych:n";		
	
	// mamy n+1 punktow węzłowych dlatego zaczynamy od 0
	for(double i=0; i<=n; i++)
	{
		tmp_wsp.x = h*i;
		tmp_wsp.y = f(h*i);
		tab.push_back(tmp_wsp);	// wstawiamy na koniec vectora nasze wspolrzedne
		
		// podgląd wspolrzednych
		cout<<"Punkt "<<i<<" # x: "<<tab[i].x<<", y: "<<tab[i].y<<"n";		
	}
	
	// każdy podprzedział jest dzielony na kolejne 1000 podprzedzialow w czasie calkowania
	L2 = sqrt(calka_simpsona_z_g(a, b, n*1000, tab));
	
	cout<<"nNorma L2 = "<<L2<<"nnn";

	if(L2 <= eps)
	{
		break; // przerywamy pętle gdy norma jest mniejsza lub rowna zadanej dokładności	
	} 
	else
	{			
		n = pow(2., ++licznik); // zwiekszamy liczbe podprzedzialow
	}

}

// po wyjściu z pętli wypisujemy odpowiedź.
cout<<"Wymagana liczba podprzedzialow to: "<<n;

/*
	W ostatnim kroku zostanie wypisane:

	Norma L2 = 9.4489457156375221e-006

	Wymagana liczba podprzedzialow to: 512
*/




// testowe sprawdzenie jak działają funkcje:	
//cout<<"g: "<<g(1, tab);
//cout<<"finterp: "<<finterp(1.25, tab);

cout<<"nn";
system("PAUSE");
return EXIT_SUCCESS;

} [/code]