{  diff22.pas ver.1       }
{ (c) Kenrou Adachi, 1993 }

program diffnorm(input, output);

{ for TP4 crt --> crtansi }

uses crt, graph;

{$I data.pas}

procedure data_in(var a, xo, yo, ts, te : real);
begin
    writeln;
    writeln('input a for the range of x-axes for [-a,a]');
    write('    a = '); readln(a);
    writeln('input xo, dxo for initial value of x, dx/dt');
    write('    xo = '); readln(xo);
    write('   dxo = '); readln(yo);
    writeln('input ts and te for time range [ts,te]');
    write('    ts = '); readln(ts);
    write('    te = '); readln(te)
end;

procedure wind(a, ts, te : real; var dx, dt, ox, ot : real);
begin
    dx := 2.0 * a;
    dt := te - ts;
    ox := 400.0 * a / dx;
    ot := -640 * ts / dt
end;

procedure axes;
begin
    setcolor(BLUE);
    line(0, 0, 0, 400);
    line(0, 200, 640, 200)
end;

procedure draw(xo, yo, ts, te, dx, dt, ox, ot : real);
var t, h : real;
    i : integer;

    procedure pre_draw_graph(xo, ts, dx, dt, ox, ot : real);
    var oldx, oldt : real;
    begin
	oldx := -400.0 * xo / dx + ox;
	oldt := 640.0 * ts / dt + ot;
	moveto(round(oldt), round(oldx))
    end;

    procedure run_kut(var xo, yo, t : real; h : real);
    var x1, x2, x3, x4, y1, y2, y3, y4 : real;
    begin
	x1 := h * fx(xo, yo, t);
	y1 := h * fy(xo, yo, t);
	x2 := h * fx(xo + x1 / 2.0, yo + y1 / 2.0, t + h / 2.0);
	y2 := h * fy(xo + x1 / 2.0, yo + y1 / 2.0, t + h / 2.0);
	x3 := h * fx(xo + x2 / 2.0, yo + y2 / 2.0, t + h / 2.0);
	y3 := h * fy(xo + x2 / 2.0, yo + y2 / 2.0, t + h / 2.0);
	x4 := h * fx(xo + x3, yo + y3, t + h);
	y4 := h * fy(xo + x3, yo + y3, t + h);

        xo := xo + (x1 + 2.0 * x2 + 2.0 * x3 + x4) / 6.0;
	yo := yo + (y1 + 2.0 * y2 + 2.0 * y3 + y4) / 6.0;
	t := t + h
    end;

    procedure draw_graph(xo, dx, dt, ox, ot, t : real);
    begin
        setcolor(YELLOW);
	lineto(round(640.0 * t / dt + ot), round(-400.0 * xo / dx + ox))
    end;

begin
    h := (te - ts) / 640.0;
    t := ts;
    pre_draw_graph(xo, t, dx, dt, ox, ot);
    for i := 1 to 640 do
    begin
	run_kut(xo, yo, t, h);
	draw_graph(xo, dx, dt, ox, ot, t)
    end
end;

var a, xo, yo, ts, te : real;
    dx, dt, ox, ot : real;
    graphdriver, graphmode : integer;
begin
    clrscr;
    graphdriver := DETECT;
    initgraph(graphdriver, graphmode, '');
    cleardevice;

    repeat
	repeat
	    data_in(a, xo, yo, ts, te);
	    wind(a, ts, te, dx, dt, ox, ot)
        until (a > 0.0) and (a <= 128)
    until (abs(xo) <= abs(a));

    axes;
    draw(xo, yo, ts, te, dx, dt, ox, ot);

    gotoxy(0, 35);
    writeln('put return key !');
    readln
end.