program agmath3b; {produziert Juliamenge f�r f(z)=z^2-l und Orbits}
                 {mit Grafik-Schrift f�r windows xp; G.Helmberg, 27.11.2006}
uses crt,graph,graph0,dialog2,xpzeich1;

type vector = array[1..2] of real;
     bit = 0..1;
     vectorfile = file of vector;
     zarray = array[1..4] of vector;
     boolearray = array[1..4] of boolean;
     triple = array[1..3] of word; {(tr[1],tr[2])=Pixel, tr[3]=Farbe}
     screenfile = file of triple;  {nimmt Pixelpunkt mit Farbe auf}

var z,l,zero:vector;
    ch,chm,chc:char;
    am,n,nenner,numits,minxpaint,maxxpaint,minypaint,maxypaint,i,j,
      graphdriver,graphmode,x0,y0,x1,y1,k:integer;
    driverpath:string20;
    links,rechts,unten,oben,r,r0,r1,r2:real;
    h:bit;
    out,conn,call,save,col:boolean;
    larray:array[1..16] of vector;
    fix:zarray;
    doit,attfix:boolean;
    att:boolearray;
    filename:string;
    periodfile,juliaold,julianew:vectorfile;
    screen:screenfile;

procedure readrsafe(var r:real);
 var code:integer;
     st:string;
     x,y:integer;
begin
x:=wherex;
y:=wherey;
repeat
  gotoxy(x,y);
  reset(input);
  read(st);
  val(st,r,code);
  until code=0;
end;  { readrsafe }

procedure metrik(var x0,y0,x1,y1:integer);
begin
x0:=bildx(links,rechts,0);
y0:=bildy(unten,oben,0);
x1:=bildx(links,rechts,1);
y1:=bildy(unten,oben,1);
end;  {metrik}

procedure bild(ch:char;var links,rechts,unten,oben:real;
               var x0,y0,x1,y1:integer);
begin
case ch of
  'b':begin
    links:=-2.4;
    rechts:=4.8;
    unten:=-2.7;
    oben:=2.7;
    end;
  'w':begin
    links:=-0.4;
    rechts:=0.8;
    unten:=-0.05;
    oben:=0.85;
    end;  {'w'}
  'm':begin
    links:=-0.1;
    rechts:=1.1;
    unten:=-0.05;
    oben:=0.85;
    end;  {'m'}
  'i':begin
    links:=-2.4;
    rechts:=2.4;
    unten:=-1.8;
    oben:=1.8;
    end;  {'i'}
  'z':begin
    links:=-4;
    rechts:=4;
    unten:=-3;
    oben:=3;
    end;  {'z'}
  'j':begin
    links:=-2.4;
    rechts:=2.4;
    unten:=-1.8;
    oben:=1.8;
    end;  {'j'}
  'a':begin
    links:=-1.5;
    rechts:=2.5;
    unten:=-1.5;
    oben:=1.5;
    end;  {'a'}
  end;  {case ch}
metrik(x0,y0,x1,y1);
end;  {bild}

procedure xpnumberxy(z:real;a,d,x,y:integer;var xarr:xarray;color:word);
  var i,j,k,xz,yz,n,code:integer;  {Ausgabe der Zahl z im Grafikmodus f�r}
      st:string;                   {das Betriebssystem XP; Anfang im linken}
      ch,chi:char;                 {oberen Pixelpunkt (x,y)}
      oldcolor:word;
begin
oldcolor:=getcolor;
setcolor(color);
xar(color,xarr);
str(z:a:d,st);
for i:=1 to length(st) do
  begin
  xz:=x+(i-1)*width;
  yz:=y;
  chi:=st[i];
  n:=ord(chi);
  for j:=0 to 7 do
    for k:=0 to 7 do
      begin
      xz:=x+(i-1)*8+k;
      yz:=y+j;
      putpixel(xz,yz,xarr[n,j,k]);
      end;   {k-Schleife}
  end;  {i-Schleife}
setcolor(oldcolor);
end;  {xpnumberxy}

procedure xptextxy(st:string;x,y:integer;var xarr:xarray;color:word);
  var i,j,k,xz,yz,n:integer;  {Ausgabe des Strings st im Grafikmodus f�r das}
      ch,chi:char;            {Betriebssystem XP; Anfang im linken}
      oldcolor:word;                        {oberen Pixelpunkt (x,y)}
begin
oldcolor:=getcolor;
setcolor(color);
xar(color,xarr);
for i:=1 to length(st) do
  begin
  xz:=x+(i-1)*width;
  yz:=y;
  chi:=st[i];
  n:=ord(chi);
  if n>63 then outtextxy(xz,yz,chi)    {chi = gew”hnlicher Buchstabe}
    else                               {chi = Zeichen}
    for j:=0 to 7 do
      for k:=0 to 7 do
        begin
        xz:=x+(i-1)*8+k;
        yz:=y+j;
        putpixel(xz,yz,xarr[n,j,k]);
        end;   {k-Schleife}
  end;  {i-Schleife}
setcolor(oldcolor);
end;  {xptextxy}

procedure axes(color:word);
begin
setcolor(color);
line(width,y0,maxx-width,y0);
line(x0,height,x0,maxy-height);
horimarkexy(x0,y1);
vertmarkexy(x1,y0);
vertmarkexy(2*x1-x0,y0);
vertmarkexy(2*x0-x1,y0);
xpnumberxy(0,1,0,x0-width-3,y0+height,xarr,color);
xpnumberxy(1,1,0,x0-width-2,y1-height,xarr,color);
xpnumberxy(1,1,0,x1-width,y0+height,xarr,color);
xpnumberxy(2,1,0,2*x1-x0-width-3,y0+height,xarr,color);
xpnumberxy(-1,1,0,2*x0-x1-2*width,y0+height,xarr,color);
end;  {axes}



function nmod16(n:integer):integer;
begin
nmod16:=n mod(16);
end;

procedure abfrage(str:string;var doit:boolean);
  var ch:char;
      len:integer;
begin
len:=72-length(str);
clearwindow(len*width,maxy-2*height,maxx,maxy);
xptextxy(str+' (j/n) ?',len*width,maxy-2*height,xarr,7);
repeat ch:=readkey until ch in ['j','n'];
doit:=ch='j';
clearwindow(len*width,maxy-2*height,maxx,maxy);
end;  {abfrage}

procedure fmap(var z,l:vector); {bildet z=(x+iy) ab in f(z)=z^2-l}
  var xein,yein:real;
begin
xein:=z[1];
yein:=z[2];
z[1]:=sqr(xein)-sqr(yein)-l[1];      { Re(f(z)) }
z[2]:=2*xein*yein-l[2];              { Im(f(z)) }
end;  {fmap}

procedure zeroescape(r:real;var n,numits:integer;l:vector;var out:boolean);
  var x,y,r2,r2xy:real;
      z:vector;
begin                {gibt Zeitpunkt n an, zu dem |fmap^n(x,y,lambda)|>r}
z:=zero;
r2:=sqr(r);
n:=0;
out:=false;
repeat
  n:=n+1;              { Fluchtzeit }
  fmap(z,l);           { z=(x+iy)->z^2-l }
  r2xy:=sqr(x)+sqr(y);
  until ((r2xy>=r2) or (n>=numits));
out:=n<numits;         { l geh”rt zur Menge, f�r die Fluchtzeit(0)=n }
end;  {zeroescape}

procedure coord(xpix,ypix:integer;var z:vector);
begin     { Koordinaten eines Punktes mit Pixel-Koordinaten (xpix,ypix) }
z[1]:=links+xpix*(rechts-links)/maxx;
z[2]:=unten+(maxy-ypix)*(oben-unten)/maxy;
end;  { coord}

procedure pixel(z:vector;var xpix,ypix:integer);
begin                              {Pixelkoordinaten eines Punktes (x,y)}
xpix:=bildx(links,rechts,z[1]);
ypix:=bildy(unten,oben,z[2]);
end;  {pixel}

procedure circles(z:vector;r,r0:real;attfix:boolean);
  var x,y:integer;
      rad,rad0,col:word;
begin
col:=getcolor;
setcolor(4);
rad:=bildx(links,rechts,links+r);
circle(x0,y0,rad);
if attfix then
  begin
  pixel(z,x,y);
  rad0:=bildx(links,rechts,links+r0);
  circle(x,y,rad0);
  end;
setcolor(col);
end;  {circles}

function absz(z:vector):real;           {Absolutbetrag |z| von z}
  var r2:real;
begin
r2:=sqr(z[1])+sqr(z[2]);
if r2=0 then absz:=0
  else absz:=sqrt(r2);
end;  { absz }

procedure complexroot(var z,w:vector);  { berechtnet w=sqrt(z) }
  var xl,yl,radx,rady:real;             { mit positiver w[1]=Re(w) }
begin
xl:=z[1];
yl:=z[2];
radx:=xl+sqrt(sqr(xl)+sqr(yl));
rady:=-xl+sqrt(sqr(xl)+sqr(yl));
if radx>0 then w[1]:=sqrt(radx/2)     { Verhinderung von Laufzeitfehlern }
  else w[1]:=0;                       { infolge Wurzel aus negativer Zahl }
if rady>0 then w[2]:=sqrt(rady/2)
  else w[2]:=0;
if yl<0 then w[2]:=-w[2];
end;    {complexroot}

procedure inversmap(h:bit;var l,z,w:vector);
  var z0:vector;    { ausgegeben wird w=w(z)=+-sqrt(z+l) }
      xl,yl,radx,rady:real;
begin
z0[1]:=z[1]+l[1];
z0[2]:=z[2]+l[2];
complexroot(z0,w);
if h=1 then
  begin
  w[1]:=-w[1];
  w[2]:=-w[2];
  end;  {i=1}
end;    {inversmap}

procedure applyinversmap(l:vector;var screen:screenfile);
  var tr:triple;   {auf die Punkte z, die den Bildschirmpixeln entsprechen,}
      i,j,p:word;    {werden die Abbildungen w=+-û(z+l) ausge�bt; die den}
      dx,dy:real;  {Bildpunkten w entsprechenden Pixel werden}
      z,w:vector;  {zusammen mit dem Farbwert als Tripel tr in}
      ibit:bit;    {der Datei screen0 abgelegt}
begin
dx:=(rechts-links)/maxx;    {Pixel-Inkrement in x-Richtung}
dy:=(oben-unten)/maxy;      {Pixel-Inkrement in y-Richtung}
assign(screen,'screen0');
rewrite(screen);            {Tripel-Datei 'screen0' nimmt Bild-Tripel auf}
for i:=0 to maxx do         {Bildschirm wird Spaltenweise bearbeitet}
  begin
  z[1]:=dx*(i-x0);
  for j:=height to maxy-2*height do
    begin
    p:=getpixel(i,j);
    if getpixel(i,j)<>0 then   {zur Urbildmenge geh”rige Pixel werden erfaát}
      begin
      z[2]:=dy*(y0-j);
      for ibit:=0 to 1 do
        begin
        tr[3]:=ibit+1;          {ihre Bilder mit der Farbe ibit+1 versehen}
        inversmap(ibit,l,z,w);            {die zugeh”rigen Punkte abgebildet}
        tr[1]:=bildx(links,rechts,w[1]);  {und die zugeh”rigen Pixel}
        tr[2]:=bildy(unten,oben,w[2]);
        if ((tr[1]>maxx) or (tr[2]>maxy)) then
          ch:=readkey;
        write(screen,tr);                 {in 'screen0' gespeichert}
        end;  {ibit-Schleife}
      end;  {getpixel(i,j)<>0}
    end;  {j-Schleife}
  end;  {i-Schleife}
close(screen);
cleardevice;                         {in den geleerten Bildschirm}
reset(screen);
repeat
  read(screen,tr);                   {werden die aus 'screen0' eingelesenen}
  putpixel(tr[1],tr[2],tr[3]);       {Bildpixel eingesetzt}
  until eof(screen);                 {und so das Bild ausgegeben}
close(screen);
end;  {applyinversmap}

procedure initiateset(fix:zarray;r,r0:real;var screen:screenfile;
                      var attfix:boolean;att:boolearray;var k:integer);
  var paint:boolean;    { Initialisierung der Datei 'screen0' f�r die }
        i,j:word;       { Kontraktionsdarstellung der Julia-Menge }
      rad0,rad:real;
      z,z0:vector;
      tr:triple;
begin
k:=0;
attfix:=false;
assign(screen,'screen0');
rewrite(screen);
for i:=0 to maxx do
  for j:=0 to maxy do
    begin
    coord(i,j,z);
    rad:=sqrt(z[1]*z[1]+z[2]*z[2]);   { Keisscheibe um 0 mit Radius rad }
    paint:=rad<=r;
    if (att[1] or att[2]) then
      begin
      attfix:=true;
      if att[1] then
        begin
        z0:=fix[1];
        k:=1;
        end  {att[1]}
        else
          begin
          z0:=fix[2];
          k:=2;
          end;  {att[2]}
      rad0:=sqrt((z[1]-z0[1])*(z[1]-z0[1])+(z[2]-z0[2])*(z[2]-z0[2]));
      paint:=paint and (rad0>=r0);  { Auáengebiet der Kreisscheibe }
      end;  {att}                    { mit Radius rad0 um Fixpunkt z0 }
    if paint then
      begin
      putpixel(i,j,1);
      tr[1]:=i;
      tr[2]:=j;
      tr[3]:=1;
      write(screen,tr);
      end;  {paint}
    end;   {j-Schleife}
  close(screen);
end;   { initiateset }

procedure fixpoints(var l:vector;var fix:zarray);
  var rad,root:vector;      {Fixpunkte und Periode 2-Punkte von f}
begin
rad[1]:=1+4*l[1];           {rad = Radikand f�r Fixpunkte}
rad[2]:=4*l[2];
complexroot(rad,root);
fix[1,1]:=(1+root[1])/2;    {Fixpunkte}
fix[1,2]:=root[2]/2;
fix[2,1]:=(1-root[1])/2;
fix[2,2]:=-root[2]/2;
rad[1]:=4*l[1]-3;           {rad = Radikand f�r Periode 2 - Punkte}
rad[2]:=4*l[2];
complexroot(rad,root);
fix[3,1]:=(-1+root[1])/2;   {Periode 2 Orbitpunkte}
fix[3,2]:=root[2]/2;
fix[4,1]:=(-1-root[1])/2;
fix[4,2]:=-root[2]/2;
end;  {fixpoints}

procedure showpoint(w:vector;color:word);
  var xpix,ypix:integer;
begin
pixel(w,xpix,ypix);
xpix:=bildx(links,rechts,w[1]);
ypix:=bildy(unten,oben,w[2]);
putpixel(xpix,ypix,color);
end;  { showpoint }

procedure showkreuz(w:vector;color:word);
  var xpix,ypix:integer;
begin
pixel(w,xpix,ypix);
kreuz(xpix,ypix,color);
end;  { showkreuz }

procedure blinkkreuz(z:vector;color1,color2:word);
  var xpix,ypix:integer;
begin
reset(input);
pixel(z,xpix,ypix);
repeat
  kreuz(xpix,ypix,color1);
  delay(200);
  kreuz(xpix,ypix,color2);
  delay(200);
  until keypressed;
end;  {blinkkreuz}

procedure galerieinfo;
  var l,z1,z2,z3,z4:vector;
      fix:zarray;
      g:integer;
begin
g:=1;
repeat
  clrscr;
  writeln;
  write(' Daten der Juliamenge nr. ');
  repeat readgwritexy(wherex,wherey,2,g)
    until ((g>0) and (g<17));
  writeln;
  writeln;
  l:=larray[g];
  fixpoints(l,fix);
  z1:=fix[1];
  z2:=fix[2];
  z3:=fix[3];
  z4:=fix[4];
  z:=z1;
  if 2*absz(z1)<=1 then z:=z2;   { Wahl eines abstoáenden Fixpunktes }
  writeln('c = (',l[1]:8:5,',',l[2]:8:5,')');
  writeln;
  writeln('Fixpunkte:');
  write('z1 = (',fix[1,1]:8:5,',',fix[1,2]:8:5,'); |f''(z1)| = ',
    2*absz(fix[1]):4:2);
  if 2*absz(z1)<1 then writeln('; attraktiv')
    else writeln;
  write('z2 = (',fix[2,1]:8:5,',',fix[2,2]:8:5,'); |f''(z2)| = ',
    2*absz(fix[2]):4:2);
  if 2*absz(z2)<1 then writeln('; attraktiv')
    else writeln;
  writeln;
  writeln('Orbit mit Periode 2:');
  writeln('z3 = (',fix[3,1]:8:5,',',fix[3,2]:8:5,
            '); |fý''(z3)| = |fý''(z4)| =', 4*absz(fix[3])*absz(fix[4]):4:2);
  write('z4 = (',fix[4,1]:8:5,',',fix[4,2]:8:5,')');
  if 4*absz(z3)*absz(z4)<1 then writeln('; attraktiv')
    else writeln;
  writeln;
  writeln(' Ausgabe weiterer Daten (j/n)?');
  repeat ch:=readkey until ch in ['j','n'];
  until ch='n';
end;  {galerieinfo}

procedure showfix(l:vector;var fix:zarray);
begin
fixpoints(l,fix);
showkreuz(l,7);                { Ausgabe der Funktionskonstante in grau}
showkreuz(fix[1],4);           { Graphikausgabe der Fixpunkte in rot}
showkreuz(fix[2],4);
showkreuz(fix[3],2);           { Ausgabe der 2-er Periode in gr�n}
showkreuz(fix[4],2);
end;  {showfix}

procedure juliainfo(l:vector;var fix:zarray;var att:boolearray;
                    var k:integer);
  var z,z1,z2,z3,z4:vector;
  strl1,strl2,strz11,strz12,strz21,strz22,strz31,strz32,strz41,strz42,
  strz1p,strz2p,strz3p:string;
  i:integer;
begin
for i:=1 to 4 do att[i]:=false;
k:=0;
showfix(l,fix);
z1:=fix[1];
z2:=fix[2];
z3:=fix[3];
z4:=fix[4];
str(l[1]:7:4,strl1);
str(l[2]:7:4,strl2);
str(z1[1]:7:4,strz11);
str(z1[2]:7:4,strz12);
str(z2[1]:7:4,strz21);
str(z2[2]:7:4,strz22);
str(z3[1]:7:4,strz31);
str(z3[2]:7:4,strz32);
str(z4[1]:7:4,strz41);
str(z4[2]:7:4,strz42);
str(2*absz(z1):3:1,strz1p);
str(2*absz(z2):3:1,strz2p);
str(4*absz(z3)*absz(z4):3:1,strz3p);
z:=z1;
if 2*absz(z1)<=1 then z:=z2;   { Wahl eines abstoáenden Fixpunktes }
clearwindow(0,0,50,8*height);
setcolor(7);
xptextxy('c = ('+strl1+','+strl2+')',10*width,height,xarr,7);
xptextxy('Fixpunkte:',10*width,2*height,xarr,4);
xptextxy('z1=('+strz11+','+strz12+'); abs f''(z1) = '+strz1p,10*width,
  3*height,xarr,4);
if 2*absz(z1)<1 then xptextxy(' attraktiv',55*width,3*height,xarr,4);
xptextxy('z2=('+strz21+','+strz22+'); abs f''(z2) = '+strz2p,10*width,
  4*height,xarr,4);
if 2*absz(z2)<1 then
  xptextxy('attrakiv',55*width,4*height,xarr,4);
att[1]:=2*absz(z1)<1;
att[2]:=2*absz(z2)<1;
xptextxy('Orbit mit Periode 2:',10*width,5*height,xarr,2);
xptextxy(
  'z3=('+strz31+','+strz32+'); abs fý''(z3)=abs fý''(z4)='+strz3p,
  10*width,6*height,xarr,2);
xptextxy('z4=('+strz41+','+strz42+')',10*width,7*height,xarr,2);
if 4*absz(z3)*absz(z4)<1 then
  begin
  xptextxy('attraktiv',45*width,7*height, xarr,2);
  att[3]:=true;
  att[4]:=true;
  end;  {4*abbsz(z3)*absz(z4)<1}
for i:=1 to 4 do
  if att[i]=true then k:=i;
end;  {juliainfo}

procedure iterationsinfo(k:integer;r,r0:real);
begin
xptextxy('Radius r = ',33*width,height,xarr,7);
xpnumberxy(r,3,1,45*width,height,xarr,7);
if k<>0 then
  begin
  xptextxy(', Innen-Radius r0 = ',48*width,height,xarr,7);
  xpnumberxy(r0,3,1,67*width,height,xarr,7);
  xptextxy('um Fixpunkt z',50*width,2*height,xarr,7);
  xpnumberxy(k,1,0,63*width,2*height,xarr,7);
  end;  {k<>0}
end;  {iterationsinfo}

procedure fixinfo(k:integer);
begin
xptextxy('Fixpunkt z',33*width,height,xarr,7);
xpnumberxy(k,1,0,43*width,height,xarr,7);
end;  {fixinfo}

procedure fluchtinfo(r:real;n:integer);
begin
xptextxy('Fluchtradius r =',33*width,height,xarr,7);
xpnumberxy(r,3,1,51*width,height,xarr,7);
xptextxy(', Fluchzeit n =',54*width,height,xarr,7);
xpnumberxy(n,4,0,70*width,height,xarr,7);
end;  {fluchtinfo}

procedure fixpointjulia(l:vector;var juliaold,julianew:vectorfile);
  var z1,z2,z,w:vector;
      fix:zarray;
      strl1,strl2,strz11,strz12,strz21,strz22,strz1p,strz2p:string;
      ch:char;    { Grafik-Ausgabe der Julia-Menge durch Urbildpunkte}
      n,k:integer;  { eines abstoáenden Fixpunktes von f }
begin
clearviewport;
axes(7);
n:=0;
juliainfo(l,fix,att,k);
xptextxy('Fixpunkt Nummer k =',width,maxy-2*height,xarr,7);
readln(k);
z:=fix[k];
{if 2*absz(z)<1 then z:=fix[2];  Wahl eines abstoáenden Fixpunktes }
fixinfo(k);
inversmap(0,l,z,w);
assign(julianew,'newjulia');   { Aufnahmedatei f�r neue Urbildpunkte }
assign(juliaold,'oldjulia');   { Ablesedatei f�r Punkte, deren Urbilder }
rewrite(juliaold);             { noch zu bestimmen sind }
write(juliaold,w);             { dessen zwei Urbilder in juliaold }
inversmap(1,l,z,w);            { aufgenommen werden }
write(juliaold,w);
close(juliaold);
repeat
  clearwindow(0,maxy-2*height,maxx,maxy);
  n:=n+1;
  reset(juliaold);             { juliaold wird zum lesen von z ge”ffnet, }
  rewrite(julianew);           { julianew zum schreiben }
  repeat
    read(juliaold,z);
    inversmap(0,l,z,w);        { die zwei Urbilder w von z werden }
    showpoint(w,1);            { auf dem Bildschirm ausgegeben und }
    write(julianew,w);         { in julianew gespeichert }
    inversmap(1,l,z,w);
    showpoint(w,1);
    write(julianew,w);
    until eof(juliaold);
  close(juliaold);
  seek(julianew,0);           { Zeiger auf Beginn von julianew }
  rewrite(juliaold);          { šbertragung des Inhaltes auf juliaold }
  repeat
    read(julianew,w);
    write(juliaold,w);
    until eof(julianew);
  close(julianew);
(*  juliainfo(l,fix,att,k);*)
  xptextxy('Stufe:',width,maxy-height,xarr,7);
  xpnumberxy(n,2,0,9*width,maxy-height,xarr,7);
  xptextxy('weiter (j/n) ?',maxx-15*width,maxy-height,xarr,7);
  repeat ch:=readkey until ((ch='j') or (ch='n'));
  until ch='n';               { Abfrage Prozedur-Wiederholung }
end;  { fixpointjulia }

procedure pointescape(var r:real;var n,numits:integer;z,l:vector;
                      var out:boolean);
  var r2,r2xy:real;
begin                { gibt Zeitpunkt n an, zu dem |fmap^n(z,l)|>r }
r2:=sqr(r);          { wenn n=numits, wird fmap^n(z,l) beschr„nkt }
n:=0;                { angenommen }
out:=false;
repeat
  n:=n+1;              { Fluchtzeit }
  fmap(z,l);           { z=(x+iy)->z^2-l }
  r2xy:=sqr(z[1])+sqr(z[2]);
  until ((r2xy>=r2) or (n>=numits));
out:=n<numits;         { Fluchtzeit(z)=n }
end;  {pointescape}

procedure zykloid(color:word); {Graph von c = (e^(2i phi)-2e^(i phi))/4}
  var phi:real;
      i,j,m:integer;
      c:vector;
begin
for  m:=0 to 1080 do
  begin
  phi:=m*pi/540;
  c[1]:=(cos(2*phi)-2*cos(phi))/4;
  c[2]:=(sin(2*phi)-2*sin(phi))/4;
  i:=bildx(links,rechts,c[1]);
  j:=bildy(unten,oben,c[2]);
  putpixel(i,j,color);
  end;
end;  {zykloid}

procedure grid;  {�berdeckt das Einheitsquadrat mit Koordinatenlinien}
  var i:integer;
begin            {im Graphmodus}
setcolor(3);  {t�rkis}
for i:=-10 to 20 do
  begin
  line(bildx(links,rechts,i/10),y0,bildx(links,rechts,i/10),y1);
  if i mod 5 = 0 then
    xpnumberxy(i/10,3,1,bildx(links,rechts,i/10)-12,y1-2*height,xarr,7);
  end; {i-Schleife}
for i:=0 to 10 do
  begin
  line(bildx(links,rechts,-1),bildy(unten,oben,i/10),
       bildx(links,rechts,2),bildy(unten,oben,i/10));
  xpnumberxy(i/10,3,1,2*x1-x0+width,bildy(unten,oben,i/10),xarr,7);
  end;  {i-Schleife}
setcolor(1);  {blau}
end;  {grid}

procedure verhulst;
  var n,m,xx,y:integer;
      x,l:real;
begin
clearwindow(0,0,maxx,y0);
axes(7);
m:=bildy(unten,oben,1/2);
setcolor(4);
line(width,m,maxx-width,m);
xpnumberxy(1,1,0,x0-width,m-height,xarr,4);
xpnumberxy(0,1,0,x0-width,y0-height,xarr,4);
for xx:=bildx(links,rechts,-0.25) to bildx(links,rechts,2) do
  begin
  n:=0;
  x:=0;
  l:=links+xx*(rechts-links)/maxx;
  repeat
    n:=n+1;
    x:=x*x-l;
    if n>5000 then
      begin
      y:=bildy(unten,oben,-x/2);
      putpixel(xx,y,4);
      end;  {n>5000}
    until n=5120;
  end;  {xx-Schleife}
end;  {verhulst}

procedure juliamenge(r:real;numits:integer;l:vector;fix:zarray);
  var i,j,n:integer;
      anfang,ende,x,y:real;
      z:vector;
begin
for i:=bildx(links,rechts,-2) to bildx(links,rechts,2) do
  for j:=bildy(unten,oben,1.5) to bildy(unten,oben,-1.5) do
    begin
    coord(i,j,z);
    pointescape(r,n,numits,z,l,out);
    if not out then
      begin
      putpixel(i,j,1);
      end;
    end;  {j-Schleife}
axes(7);
showkreuz(l,7);
for i:=1 to 2 do
  showkreuz(fix[i],4);                {Fixpunkte}
for i:=3 to 4 do
  showkreuz(fix[i],2);               {Orbit mit Periode 2}
xptextxy('grid (j/n) ?',width,maxy-2*height,xarr,7);
repeat ch:=readkey until ch in ['j','n'];
if ch='j' then  grid;
clearwindow(0,maxy-2*height,maxx,maxy);
end;  {juliamenge}

procedure lwahl(var l:vector;chm:char);   {Wahl des Vektors l in f(z)=z^2-l}
  var abs:real;
      nr:integer;
begin                            {in Graphmode als Punkt der Zahlenebene}
nr:=1;
bild('a',links,rechts,unten,oben,x0,y0,x1,y1);
setcolor(7);
axes(7);
line(2*x0-x1,0,2*x0-x1,maxy);
circle(x0,y0,2*(x1-x0));
zykloid(7);
circle(x1,y0,(x1-x0) div 4);
repeat
  setcolor(1);
  xptextxy('f(z) = zý-l',3*width,4*height,xarr,1);
  if ((chm='k') or (chm='i')) then
    begin
    xptextxy('c[1] = ',3*width,5*height,xarr,1);
    readrsafe(l[1]);        {Eingabe von l[1]}
    xpnumberxy(l[1],8,6,10*width,5*height,xarr,1);
    xptextxy('c[2] = ',3*width,6*height,xarr,1);
    readrsafe(l[2]);        {Eingabe von l[2]}
    xpnumberxy(l[2],8,6,10*width,6*height,xarr,1);
    end  {chm='k','i'}
    else begin
      xptextxy('Nummer nr:',3*width,5*height,xarr,1);
      readln(nr);
      xpnumberxy(nr,2,0,13*width,5*height,xarr,1);
      l:=larray[nr];
      end;  {chm='g'}
  showkreuz(l,7);
  xptextxy('akzeptabel (j/n) ?',3*width,8*height,xarr,1);
  repeat ch:=readkey until ((ch='j') or (ch='n'));
  clearwindow(0,0,24*width,9*height);
  if ch='n' then
     kreuz(bildx(links,rechts,l[1]),bildy(unten,oben,l[2]),15);
  until ch='j';
end;  {lwahl}

procedure zwahl(var l,z:vector;var doit:boolean);
  var abs:real;                  {Wahl des Anfangs-Vektors z in f(z)=z^2-l}
      l1,l2,z1,z2:string;        {in Graphmode als Punkt der Zahlenebene}
begin
setcolor(7);
if doit then
  begin
  clearwindow(0,maxy-2*height,maxx,maxy);
  clearwindow(0,11*height,26*width,16*height);
  repeat
    setcolor(9);
    xptextxy('Orbit-Anfang z:',10*width,11*height,xarr,9);
    xptextxy('z[1] = ',10*width,12*height,xarr,9);
    readrsafe(z[1]);        {Eingabe von z[1]}
    xptextxy('z[1] = ',10*width,12*height,xarr,9);
    xpnumberxy(z[1],4,2,18*width,12*height,xarr,9);
    xptextxy('z[2] = ',10*width,13*height,xarr,9);
    readrsafe(z[2]);        {Eingabe von z[2]}
    xpnumberxy(z[2],4,2,18*width,13*height,xarr,9);
    xptextxy('akzeptabel?',10*width,15*height,xarr,9);
    blinkkreuz(z,14,9);
    repeat ch:=readkey until ((ch='j') or (ch='n'));
    clearwindow(0,12*height,23*width,16*height);
    if ch='j' then
      begin
      str(z[1]:4:2,z1);
      str(z[2]:4:2,z2);
      xptextxy('z=('+z1+';'+z2+')',10*width,12*height,xarr,9);
      end  {'j'}
        else showkreuz(z,1);
    until ch='j';
  end {doit}
  else clearwindow(0,maxy-2*height,maxx,maxy);
end;  {zwahl}

procedure attcycle(l,z:vector);
{Bildschirmausgabe eines attraktiven Zykels f�r f(z)=z^2-l mit Anfangswert z}
  var zold,zdif:vector;
      x,y,i:integer;
      n:longint;
      ch:char;
      rand,flucht:boolean;
      nstring,textstring:string;
      color:word;
      r:real;
begin
n:=0;
textstring:='';
clearwindow(0,maxy-3*height,maxx,maxy);
xptextxy('Einzelschritt oder Orbit (e/o) ?',width,maxy-2*height,xarr,7);
repeat ch:=readkey until ch in ['e','o'];
repeat
  n:=n+1;
  fmap(z,l);
  r:=absz(z);
  x:=bildx(links,rechts,z[1]);
  y:=bildy(unten,oben,z[2]);
  if n>5000 then color:=2 else color:=9;
  if ((n>5000) or (ch='o') or (ch='e')) then
    kreuz(x,y,color);
  if ch='e' then
    repeat ch:=readkey until ch in ['e','o'];
  flucht:=r>2;
  rand:=((x<0) or (x>maxx) or (y<0) or (y>maxy));
  until ((n=5120) or rand or flucht);
if flucht then textstring:='; Fluchtkreis�berschreitung';
if rand then textstring:='; Rand�berschreitung';
if (rand or flucht) then
  begin
  clearwindow(0,maxy-3*height,maxx,maxy);
  str(n,nstring);
  xptextxy('n='+nstring+textstring,width,maxy-2*height,xarr,7);
  ch:=readkey;
  end {rand}
  else begin
    clearwindow(0,maxy-3*height,maxx,maxy);
    xptextxy('Ausgabe der approximativen Perioden-Punkte (j/n) ?',
              30*width,maxy-2*height,xarr,7);
    repeat ch:=readkey until ch in ['j','n'];
    clearwindow(0,maxy-3*height,maxx,maxy);
    if ch='j' then
      begin
      n:=0;
      zold:=z;
      repeat
        n:=n+1;
        fmap(z,l);
        zdif[1]:=z[1]-zold[1];
        zdif[2]:=z[2]-zold[2];
        until ((absz(zdif)<1E-5) or (n>100));
      str(n,nstring);
      if n<=100 then
        begin
        xptextxy('Periodenl„nge: '+nstring+'; <Taste>',width,maxy-3*height
                      ,xarr,7);
        repeat
          z:=zold;
          for i:=1 to n do
            begin
            fmap(z,l);
            clearwindow(0,maxy-2*height,maxx,maxy);
            xptextxy('i = ',width,maxy-2*height,xarr,7);
            xpnumberxy(i,3,0,4*width,maxy-2*height,xarr,7);
            xptextxy('z = (',8*width,maxy-2*height,xarr,7);
            xpnumberxy(z[1],8,5,13*width,maxy-2*height,xarr,7);
            xptextxy(',',21*width,maxy-2*height,xarr,7);
            xpnumberxy(z[2],8,5,22*width,maxy-2*height,xarr,7);
            xptextxy(')',30*width,maxy-2*height,xarr,7);
            blinkkreuz(z,14,2);
            ch:=readkey;
            end;  {i-Schleife}
            xptextxy('Nochmals Ausgabe der Perioden-Punkte (j/n) ?',
              35*width,maxy-2*height,xarr,7);
          repeat ch:=readkey until ch in ['j','n'];
          until ch='n';
        end {n<=100}
        else begin
          xptextxy('unzureichende Konvergenz - Taste',width,maxy-3*height,
                        xarr,7);
          ch:=readkey;
          end;  {n>100}
      clearwindow(0,maxy-3*height,maxx,maxy);
      end;  {ch='j'}
    end;  {not rand}
end;  {attcycle}

procedure orbit(r:real;numits:integer);
  var l:vector;
      l1,l2:string;
      doit:boolean;
      ch:char;
begin
bild('a',links,rechts,unten,oben,x0,y0,x1,y1);
setgraphmode(graphmode);
setbkcolor(15);
repeat
  clearviewport;
  axes(7);
  line(2*x0-x1,0,2*x0-x1,maxy);
  circle(x0,y0,2*(x1-x0));
  zykloid(7);
  circle(x1,y0,(x1-x0) div 4);
  grid;
  repeat
    gotoxy(0,0);
    clearwindow(0,maxy-3*height,maxx,maxy);
    xptextxy('w=zý-c;',width,maxy-2*height,xarr,7);
    xptextxy('c[1] = ',10*width,maxy-2*height,xarr,7);
    readrsafe(l[1]);        {Eingabe von l[1]}
    xpnumberxy(l[1],4,2,16*width,maxy-2*height,xarr,7);
    xptextxy('; c[2] = ',20*width,maxy-2*height,xarr,7);
    readrsafe(l[2]);        {Eingabe von l[2]}
    xpnumberxy(l[2],4,2,29*width,maxy-2*height,xarr,7);
    xptextxy('akzeptabel (j/n) ?',60*width,maxy-2*height,xarr,7);
    blinkkreuz(l,14,7);
    repeat ch:=readkey until ch in ['j','n'];
    clearwindow(0,maxy-3*height,maxx,maxy);
    if ch='j' then
      begin
      str(l[1]:8:6,l1);
      str(l[2]:8:6,l2);
      outtextxy(width,maxy-2*height,'l=('+l1+';'+l2+')');
      end;
    until ch='j';
  attcycle(l,zero);
  abfrage('neuer Orbit',doit);
  until not doit;
end;  {orbit}

procedure zyk;
begin
bild('a',links,rechts,unten,oben,x0,y0,x1,y1);
setcolor(7);
axes(7);
line(2*x0-x1,0,2*x0-x1,maxy);
circle(x0,y0,2*(x1-x0));
zykloid(7);
circle(x1,y0,(x1-x0) div 4);
end;  {zyk}

procedure fluchtwahl(var r:real;var numits:integer);
  var ch:char;
begin
writeln;
repeat
  write(' Fluchtradius r = ');
  readrwritexy(wherex,wherey,3,1,r);
  writeln;
  write(' Fluchtzeit n = ');
  readgwritexy(wherex,wherey,4,numits);
  writeln;
  writeln(' akzeptabel (j/n)?');
  repeat ch:=readkey until ch in ['j','n'];
  until ch = 'j';
end;  {fluchtwahl}

procedure menu(var chm:char;var r:real;var n:integer);
begin
clrscr;
writeln;
writeln(' Auswahlmenu:');
writeln(' Zykloide .................. z');
writeln(' Galeriemenge .............. g');
writeln(' Konstanteneingabe ......... k');
writeln(' Fixpunkturbilder .......... u');
writeln(' Iteration ................. i');
writeln(' zero-Orbit ................ o');
writeln(' Galeriedaten .............. d');
writeln(' Programm-Ende: ............ e');
writeln;
repeat chm:=readkey until chm in ['d','e','g','i','k','o','u','z'];
if chm='k' then fluchtwahl(r,n);
if chm='u' then
  begin
  writeln(' f�r Galeriemengen oder Konstanteneingabe (g/k)?');
  repeat ch:=readkey until ch in ['g','k'];
  if ch='k' then chm:='u'
    else chm:='v';
  end;  {chm='u'}
end; {menu}

begin{main}
graphdriver:=9;        {bgi-Datei egavga.bgi muá im aktuellen Ordner sein}
graphmode:=2;
initgraph(graphdriver,graphmode,'');
graphvariables(maxx,maxy,width,height);
restorecrtmode;
bild('b',links,rechts,unten,oben,x0,y0,x1,y1);
r:=2;                                 { escape-time Radius }
k:=0;                                 { Nummer des attraktiven Fixpunktes }
numits:=100;                          { escape-time }
zero[1]:=0;
zero[2]:=0;
z:=zero;
larray[1,1]:=-0.31;                   { Koordinaten der Konstanten c }
larray[1,2]:=0.04;                    { f�r die Juliamengen der Galerie }
larray[2,1]:=0.11;
larray[2,2]:=0.6557;
larray[3,1]:=0.12;
larray[3,2]:=0.74;
larray[4,1]:=-0.3;
larray[4,2]:=0.55;
larray[5,1]:=0.194;
larray[5,2]:=0.6557;
larray[6,1]:=0.74543;
larray[6,2]:=0.11301;
larray[7,1]:=1.25;
larray[7,2]:=0;
larray[8,1]:=0.481762;
larray[8,2]:=-0.531657;
larray[9,1]:=0.39054;
larray[9,2]:=-0.58679;
larray[10,1]:=0.15652;
larray[10,2]:=-1.03225;
larray[11,1]:=-0.11031;
larray[11,2]:=-0.67037;
larray[12,1]:=-0.27334;
larray[12,2]:=0.00742;
larray[13,1]:=0;
larray[13,2]:=0.5;
larray[14,1]:=1;
larray[14,2]:=0.1;
larray[15,1]:=1.15;
larray[15,2]:=0.25;
larray[16,1]:=1.3;
larray[16,2]:=0.05;
repeat
  menu(chm,r,numits);
  if chm in ['g','i','k','o','u','v','z'] then
    begin
    setgraphmode(graphmode);
    setbkcolor(15);
    setcolor(1);
    end;  {ch<>'e','d'}
  case chm of
  'd':galerieinfo;   { Ausgabe der Daten f�r die Galeriemengen }
  'g','k':begin      { Galeriemenge, Konstanteneingabe }
    lwahl(l,chm);    { Fluchtzeitgrafik der Juliamenge }
    bild('j',links,rechts,unten,oben,x0,y0,x1,y1);
    repeat
      cleardevice;
      juliainfo(l,fix,att,k);
      fluchtinfo(r,numits);
      juliamenge(r,numits,l,fix);
      setcolor(7);
      abfrage('Ausgabe eines Orbits',doit);
      if doit then
        begin
        zwahl(l,z,doit);
        attcycle(l,z);
        abfrage(' Neuer Orbit',doit);
        end
      until not doit;
    end;  {'j'}
  'i':begin              { Kontraktionsiteration f�r die Juliamenge }
    n:=-1;
    r0:=0;
    k:=0;
    lwahl(l,chm);
    clearviewport;
    juliainfo(l,fix,att,k);
    bild('j',links,rechts,unten,oben,x0,y0,x1,y1);
    xptextxy('Auáenradius r = ',width,maxy-2*height,xarr,7);
    readln(r);
    xpnumberxy(r,3,1,18*width,maxy-2*height,xarr,7);
    if (att[1] or att[2]) then
      begin
      xptextxy(', Innen-Radius r0 = ',21*width,maxy-2*height,xarr,7);
      readln(r0);
      xpnumberxy(r0,3,1,40*width,maxy-2*height,xarr,7);
      end;  {att}
    clearwindow(60*width,maxy-2*height,maxx,maxy);
    xptextxy('Taste',70*width,maxy-2*height,xarr,7);
    clearwindow(0,maxy-2*height,maxx,maxy);
    iterationsinfo(k,r,r0);
    abfrage('weiter',doit);
    clearviewport;
    initiateset(fix,r,r0,screen,attfix,att,k);
    repeat
      n:=n+1;
      xptextxy('Stufe',width,maxy-2*height,xarr,7);
      xpnumberxy(n,3,0,7*width,maxy-2*height,xarr,7);
      abfrage('weiter',doit);
      if doit then
        begin
        clearwindow(0,maxy-2*height,maxx,maxy);
        applyinversmap(l,screen);
        {ch:=readkey;}
        end;  {doit}
      until not doit;
    clearwindow(0,maxy-2*height,maxx,maxy);
    axes(7);
    clearwindow(60*width,maxy-2*height,maxx,maxy);
    xptextxy('Taste',70*width,maxy-2*height,xarr,7);
    ch:=readkey;
    iterationsinfo(k,r,r0);
    clearwindow(40*width,maxy-2*height,maxx,maxy);
    xptextxy('Taste',70*width,maxy-2*height,xarr,7);
    ch:=readkey;
    if n>0 then
      begin
      juliainfo(l,fix,att,k);
      circles(fix[k],r,r0,attfix);
      end
      else showfix(l,fix);
    clearwindow(40*width,maxy-2*height,maxx,maxy);
    xptextxy('Periodenpunkte (j/n) ?',55*width,maxy-2*height,xarr,7);
    repeat ch:=readkey until ch in ['j','n'];
    if ch='j' then
      begin
      attcycle(l,zero);
      end;  {ch='j', attcycle}
    clearwindow(0,maxy-2*height,maxx,maxy);
    xptextxy('Taste',70*width,maxy-2*height,xarr,7);
    ch:=readkey;
    end;  {'i'}
  'o':orbit(r,numits);
  'u','v':begin                  { Fixpunkturbilder }
    if chm='u' then lwahl(l,'k') { Eingabe der c-Koordinaten }
    else lwahl(l,'v');           { šbernahme der Galerie-c-Koordinaten }
    bild('j',links,rechts,unten,oben,x0,y0,x1,y1);
    fixpointjulia(l,juliaold,julianew);
    end;  {'u'}
  'z':begin                      { Grafik der Zykloide }
    zyk;
    ch:=readkey;
    end;  {'z'}
  end;  {case chm}
  restorecrtmode;
  until chm='e';
ch:=readkey;
end. {agmath3b}

  'i':begin              { Kontraktionsiteration f�r die Juliamenge }
    n:=-1;
    r0:=0;
    k:=0;
    lwahl(l,chm);
    clearviewport;
    juliainfo(l,fix,att,k);
    bild('j',links,rechts,unten,oben,x0,y0,x1,y1);
    xptextxy('Auáenradius r = ',width,maxy-2*height,xarr,7);
    readln(r);
    xpnumberxy(r,3,1,18*width,maxy-2*height,xarr,7);
    if (att[1] or att[2]) then
      begin
      xptextxy(', Innen-Radius r0 = ',21*width,maxy-2*height,xarr,7);
      readln(r0);
      xpnumberxy(r0,3,1,40*width,maxy-2*height,xarr,7);
      end;  {att}
    clearwindow(60*width,maxy-2*height,maxx,maxy);
    xptextxy('Taste',70*width,maxy-2*height,xarr,7);
    clearwindow(0,maxy-2*height,maxx,maxy);
    iterationsinfo(k,r,r0);
    abfrage('weiter',doit);
    clearviewport;
    initiateset(fix,r,r0,screen,attfix,att,k);
    repeat
      n:=n+1;
      xptextxy('Stufe',width,maxy-2*height,xarr,7);
      xpnumberxy(n,3,0,7*width,maxy-2*height,xarr,7);
      abfrage('weiter',doit);
      if doit then
        begin
        clearwindow(0,maxy-2*height,maxx,maxy);
        applyinversmap(l,screen);
        {ch:=readkey;}
        end;  {doit}
      until not doit;
    clearwindow(0,maxy-2*height,maxx,maxy);
    axes(7);
    clearwindow(60*width,maxy-2*height,maxx,maxy);
    xptextxy('Taste',70*width,maxy-2*height,xarr,7);
    ch:=readkey;
    iterationsinfo(k,r,r0);
    clearwindow(40*width,maxy-2*height,maxx,maxy);
    xptextxy('Taste',70*width,maxy-2*height,xarr,7);
    ch:=readkey;
    if n>0 then
      begin
      juliainfo(l,fix,att,k);
      circles(fix[k],r,r0,attfix);
      end
      else showfix(l,fix);
    clearwindow(40*width,maxy-2*height,maxx,maxy);
    xptextxy('Periodenpunkte (j/n) ?');
    repeat ch:=readkey until ch in ['j','n'];
    if ch='j' then
      begin
      attcycle(l,zero);
      end;  {ch='j', attcycle}
    clearwindow{0,maxy-2*height,maxx,maxy);
    xptextxy('Taste',70*width,maxy-2*height,xarr,7);
    ch:=readkey;
    end;  {'i'}