{$N+}PROGRAM INT;

CONST
	pi=3.141592654;
	m=5;

TYPE
        vektor=ARRAY[1..m] OF double;
	pole=ARRAY[1..m,1..m] OF double;

FUNCTION f(x:double):double;
BEGIN
	f:=sin(x)*exp(cos(x));
        {f:=x*x*x*x*ln(10+abs(x))*cos(x)*exp(-abs(x));}
END;

FUNCTION power(x,i:integer):integer;
VAR
	j:integer;
	y:integer;
BEGIN
	y:=1;
	FOR j:=1 TO i DO y:=y*x;
	power:=y;
END;

FUNCTION licho(a,b:double; n:integer):double;
VAR
	h,t,y:double;
	i:integer;

BEGIN
	h:=(b-a)/n;
	y:=0;
	FOR i:=1 to n DO BEGIN
                t:=a+i*h;
                y:=y+f(t);
	END;
	licho:=(y-0.5*f(a)-0.5*f(b))*h;
END;


VAR
	a,b,koef:double;
	romberg1,romberg2,romberg3:double;
	integral:pole;
        int1:vektor;
	i,j,kroku:integer;

BEGIN

	a:=0;
	b:=pi;

	FOR i:=1 TO m DO BEGIN
		kroku:=power(2,i-1);
                int1[i]:=licho(a,b,kroku);
		integral[i,1]:=int1[i];
		WRITELN('Pocet kroku ',kroku);
		WRITELN('Integral z funkce sin(x)*exp(cos(x)) na intervalu <',a:2:1,',',b:2:1,'> je ',integral[i,1]:8:7);
                READLN;
	END;

        {Romberg s pouzitim explicitnich vztahu - viz. Maple}

	romberg1:=4/3*int1[m]-1/3*int1[m-1];
	romberg2:=64/45*int1[m]-20/45*int1[m-1]+1/45*int1[m-2];
        romberg3:=4096/2835*int1[m]-64/135*int1[m-1]+4/135*int1[m-2]-1/2835*int1[m-3];

        WRITELN('S pouzitim explicitnich vzorcu (viz. Maple)');
	WRITELN('Integral Rombergovou metodou z integralu s kroky 16 a 8 je ',romberg1:8:7);
	WRITELN('Integral Rombergovou metodou z integralu s kroky 16, 8 a 4 je ',romberg2:8:7);
	WRITELN('Integral Rombergovou metodou z integralu s kroky 16, 8, 4 a 2 je ',romberg3:8:7);

        READLN;

        {Romberg pomoci rekurentniho vztahu}

        FOR j:=1 TO (m-1) DO BEGIN
                koef:=power(4,j);
                FOR i:=1 TO (m-j) DO integral[i,j+1]:=(koef*integral[i+1,j]-integral[i,j])/(koef-1);
        END;

	{Romberg}
	WRITELN;

        WRITELN('S pouzitim rekurentniho vzorce');
        WRITE('Romberg z integralu s pocety kroku 16 a 8');
        WRITELN(' je ',integral[m-1,2]:8:7);

        WRITE('Romberg z integralu s pocety kroku 16, 8 a 4');
        WRITELN(' je ',integral[m-2,3]:8:7);

        WRITE('Romberg z integralu s pocety kroku 16, 8, 4 a 2');
        WRITELN(' je ',integral[m-3,4]:8:7);

        WRITE('Romberg z integralu s pocety kroku 16, 8, 4, 2 a 1');
        WRITELN(' je ',integral[m-4,5]:8:7);

	READLN;

	WRITELN('Presna hodnota ma byt ',(exp(1)-exp(-1)):8:7);
        WRITELN('Primitivni funkce je -exp(cos(x))');

	READLN;

END.