Program p59r2008;
USES
	crt;
CONST
	krok=0.001;
	presnost=0.00001;
	x1=0;	{ levy okraj intervalu }
	x2=2;	{ pravy okraj intervalu }
TYPE
	vek2=array[1..2] of real;	{ pro fce }
	vek3=array[1..3] of real; 	{ pro strelbu }
PROCEDURE f(x:real; p:vek2; var dp:vek2);	{ zadani fci, [1] -> y, [2] -> z (resp. y') }
	BEGIN
		dp[1]:=p[2];
		dp[2]:=4*x*p[2]+(2-4*sqr(x))*p[1];
	END;
PROCEDURE RKstep(x:real; var p:vek2);
	VAR
		k1,k2,k3,k4,dp,pomprom:vek2;
		i:integer;
	BEGIN
		f(x,p,dp);
		FOR i:=1 TO 2 DO k1[i]:=krok*dp[i];
		FOR i:=1 TO 2 DO pomprom[i]:=p[i]+k1[i]/2;
		f(x+krok/2,pomprom,dp);
		FOR i:=1 TO 2 DO k2[i]:=krok*dp[i];
		FOR i:=1 TO 2 DO pomprom[i]:=p[i]+k2[i]/2;
		f(x+krok/2,pomprom,dp);
		FOR i:=1 TO 2 DO k3[i]:=krok*dp[i];
		FOR i:=1 TO 2 DO pomprom[i]:=p[i]+k3[i];
		f(x+krok,pomprom,dp);
		FOR i:=1 TO 2 DO k4[i]:=krok*dp[i];
		FOR i:=1 TO 2 DO p[i]:=p[i]+(k1[i]+2*k2[i]+2*k3[i]+k4[i])/6;
	END;
VAR
	x,z,alpha,beta:real;
	p:vek2;
	str,stropt:vek3;
	ans:string[1];
	fileout:text;
BEGIN
	REPEAT
		ClrScr;
		WRITELN('Rovnice zadana ve tvaru:');
		WRITELN('   ddy-4*x*dy+(4*x^2-2)*y=0');
		WRITELN('kde x je z (0;2).');
		WRITELN;
		WRITELN('Podminky:');
		WRITELN('   y(0)=alpha');
		WRITELN('   y(2)=beta');
		WRITELN;
		WRITELN('----------------------------');
		WRITELN;
		WRITE('Zadej alpha= ');
		READLN(alpha);
		WRITE('Zadej beta= ');
		READLN(beta);
		WRITELN;
		WRITELN('----------------------------');
		WRITELN;

		stropt[1]:=1;
		REPEAT
			x:=x1;
			p[1]:=alpha;
			p[2]:=stropt[1];
			WHILE (x<x2) DO BEGIN
				RKstep(x,p);
				x:=x+krok;
			END;
			stropt[2]:=p[1];
			x:=x1;
			p[1]:=alpha;
			p[2]:=-stropt[1];
			WHILE (x<x2) DO BEGIN
				RKstep(x,p);
				x:=x+krok;
			END;
			stropt[3]:=p[1];
			stropt[1]:=10*stropt[1];
		UNTIL ((stropt[2]-beta)*(stropt[3]-beta) < 0);

		str[1]:=stropt[1];
		str[2]:=-stropt[1];
		str[3]:=(str[1]+str[2])/2;
		REPEAT
			x:=x1;
			p[1]:=alpha;
			p[2]:=str[3];
			WHILE (x<x2) DO BEGIN
				RKstep(x,p);
				x:=x+krok;
			END;
			z:=str[3];
			{ zjisteni, v ktere "pulce" jsme }
			IF(p[1]-beta<0) THEN BEGIN
				str[2]:=str[3];
				str[3]:=(str[1]+str[2])/2;
			END ELSE BEGIN
				str[1]:=str[3];
				str[3]:=(str[1]+str[2])/2;
			END;
		UNTIL (abs(p[1]-beta)<presnost);

		ASSIGN(fileout,'vysl.dat');
		REWRITE(fileout);
			x:=x1;
			p[1]:=alpha;
			p[2]:=z;
			WRITELN(fileout,x:1:3,' ',p[1]:1:5);
			WHILE (x<x2) DO BEGIN
				RKstep(x,p);
				x:=x+krok;
				WRITELN(fileout,x:1:3,' ',p[1]:1:5);
			END;
		CLOSE(fileout);
		WRITELN('Reseni zapsano do souboru VYSL.DAT');
		WRITELN;
		WRITELN('Pocatecni podminka: ',str[3]:1:5);
		WRITELN;
		WRITELN('----------------------------');
		WRITELN;
		REPEAT {Jedem znovu?}
			WRITE('Znovu [y/n]? ');
			READLN(ans);
		UNTIL ((ans='y') or (ans='n'));
	UNTIL (ans='n');
END.