(*$B- Aufgabe 7.22*)
    PROGRAM systemheun;
    CONST nn=10;
    TYPE vektor = ARRAY[1..nn] of real;
    VAR x,x0,h,xend,xmax, xmin, ymin, ymax : real ;
        m,n,i,k,schritte, fall, xg, yg, xt, yt, oldx, oldy, farbe, keine: integer;
        zeich : boolean;
        y,y0,ys,k1,k2 : vektor; tf:text; line : string[120];

    (*$I plotproz *)

    PROCEDURE f(x:real; y:vektor; VAR ys:vektor);
    BEGIN
       CASE fall OF
       (*$I fdgl.pas*)
       END;
    END;

    BEGIN
       assign(tf,'fdgl.pas'); reset(tf);
       while not eof(tf) do
       begin readln(tf,line); writeln(line)  end;
       writeln('Fall eingeben:'); read(fall);
       writeln('Anzahl Gleichungen ?'); read(n);
       writeln('Anfangsbedingungen eingeben x='); read(x0);
       FOR i := 1 TO n DO BEGIN writeln('y[',i,']='); read(y0[i]) END;
       writeln('bis wohin integrieren? xend='); read(xend);
       writeln('wieviele Integrationsschritte ?'); read(schritte);
       writeln('Wertetabelle (1) oder Zeichnung (2) ?');
       read(i); zeich := i=2;
       IF zeich THEN
       BEGIN
         writeln('welche Funktion plotten ?'); read(m);
         writeln('ymin, ymax eingeben'); readln( ymin, ymax);
         xmin:=x0; xmax :=xend;
         plotbegin( xg, yg, xt, yt, farbe, keine);
         achsen(xmin, xmax, ymin, ymax, xg, yg, xt, yt, farbe, keine);
       END;
       h:=(xend-x0)/schritte; x:=x0; y:=y0;
       FOR k := 1 TO schritte DO
       BEGIN
          f(x,y,k1); FOR i:=1 TO n DO ys[i]:=y[i]+h*k1[i];
          x:=x0 + k*h;
          f(x,ys,k2);
          FOR i:=1 TO n DO y[i] := y[i] + h/2*(k1[i]+k2[i]);
          IF zeich THEN pl(x,y[m],farbe)
          ELSE
          BEGIN
             writeln('x=',x:7:5);
             FOR i:=1 TO n DO writeln('  y[',i:2,']=',y[i]);
          END
       END
    END.