Program ODE;

CONST
     h=0.5;

FUNCTION f(xn:real; yn:real):real;
BEGIN
	f:=(1+cos(xn))*yn;
END;

FUNCTION euler(xn:real; yn:real):real;
BEGIN
	euler:=yn+h*f(xn,yn);
END;

FUNCTION midpoint(xn:real; yn:real):real;
BEGIN
	midpoint:=yn+h*f(xn+h/2,yn+h/2*f(xn,yn));
END;

FUNCTION heun(xn:real; yn:real):real;
BEGIN
	heun:=yn+h/2*(f(xn,yn)+f(xn+h,yn+h*f(xn,yn)))
END;

FUNCTION rk4(xn:real; yn:real):real;
VAR
	k1,k2,k3,k4:real;
BEGIN
	k1:=h*f(xn,yn);
	k2:=h*f(xn+h/2,yn+k1/2);
	k3:=h*f(xn+h/2,yn+k2/2);
	k4:=h*f(xn+h,yn+k3);
	rk4:=yn+k1/6+k2/3+k3/3+k4/6;
END;

VAR
	i:integer;
	xn,ynE,ynM,ynH,ynR:real;
	fileout:text;

BEGIN
	{Pomoci Eulerovy metody}

	xn:=0;
	ynE:=1;
	ynM:=1;
	ynH:=1;
	ynR:=1;

	ASSIGN(fileout,'reseni.dat');
	REWRITE(fileout);

	FOR i:=1 TO 11 DO
	BEGIN
		WRITELN(fileout,xn:11:4,ynE:11:4,ynM:11:4,ynH:11:4,ynR:11:4);
		ynE:=euler(xn,ynE);
		ynM:=midpoint(xn,ynM);
		ynH:=heun(xn,ynH);
		ynR:=rk4(xn,ynR);
		xn:=xn+h;
	END;
	
	CLOSE(fileout);
        WRITELN('Hotovo');
	READLN;
END.