/**************************************************************************/
/* feigenbaum41.c                                               07.01.88  */
/* Version 1.0                                                            */
/*                                                                        */
/*                                                                        */
/* Auswertung  der Folge  p(n+1)=p(n)*(1+k*(1-p(n)))  durch Auftragen von */
/* p(n)  nach  Rechts  und  p(n+1)  nach Unten fuer fest gewaehltes k und */
/* Vergleichsparabel, um die mit Teilungspunkte verstaendlich zu machen.  */
/* (Eingabemoeglichkeit der Pararmeter aus dem CLI.)                      */
/*                                                                        */
/* fuer Markt+Technik geschrieben von Heiko Knappe                        */
/*                                                                        */
/* Komplexe Dynamische Systeme                                            */
/* Feigenbaeume - ein Wachstumsproblem                                    */
/**************************************************************************/


#include <exec/types.h>
#include <intuition/intuition.h>

#define WIDTH  640         /* Screen Modus                     */
#define HEIGHT 500
#define DEPTH  2
#define FLAGS  HIRES+LACE

#define XVAL 1.5 /* Wertebereich */
#define YVAL 1.5
#define XWIDTH 500
#define YHEIGHT 500


struct IntuitionBase *IntuitionBase;
struct GfxBase *GfxBase;
struct Screen *MyScreen;
struct RastPort *rp;
struct NewScreen MyNewScreen=
 {0,0,WIDTH,HEIGHT,DEPTH,1,0,FLAGS,CUSTOMSCREEN,0,0,0,0};

float p0=.1,pneu,palt;     /* Start-Population / Population    */
float k=2;                 /* rel. Wachstums-Rate              */
int n,nmax=10;             /* Anzahl der Iterationen           */
                          

/***********************/
/* Hauptprogramm       */
/***********************/

main(argc,argv)
int argc;
char *argv[];
{
 if((argc!=2)&&(argc!=1))
 {
    printf("Eingabe: %s <nmax,p0,k>\n",argv[0]);
    printf(" nmax : Anzahl der Iterationen\n");
    printf(" p0   : Anfangswert\n");
    printf(" k    : Wachstumsrate\n");
    close_things();
 }
 if(argc==2)
    if(sscanf(argv[1],"%d,%f,%f",&nmax,&p0,&k)!=3)
    {
       printf("falsche Eingabe\n");
       close_things();
    }
 printf("Parameter:\n nmax=%d\n p0=%f\n k=%f\n",nmax,p0,k);
 open_things();
 parabel_(); rech();
 while(MyScreen->TopEdge<1);
 close_things();
}


/***********************/
/* Initialisierungs-   */
/* routinen            */
/***********************/

open_things()
{
 IntuitionBase=OpenLibrary("intuition.library",0);
 if(IntuitionBase==0)
    close_things();
 GfxBase=OpenLibrary("graphics.library",0);
 if(GfxBase==0)
    close_things();
 MyScreen=OpenScreen(&MyNewScreen);
 if(MyScreen==0)
    close_things();
 rp=&(MyScreen->RastPort);
 init_colors();
} 

close_things()
{
 if(IntuitionBase!=0) CloseLibrary(IntuitionBase);
 if(GfxBase!=0) CloseLibrary(GfxBase);
 if(MyScreen!=0) CloseScreen(MyScreen);
 exit();
} 

init_colors()
{
 SetRGB4(&(MyScreen->ViewPort),0,0,0,0);
 SetRGB4(&(MyScreen->ViewPort),1,15,0,15);
 SetAPen(rp,1);
}


/***********************/
/* Berechnung          */
/***********************/

parabel_()
/* zeichnet die Vergleichsparabel und die Winkelhalbierende */
{
 short x,y;

 Move(rp,0,0); /* Parabel zeichnen */
 for(x=0;x<XWIDTH;++x,palt=x*(float)XVAL/XWIDTH)
 {
    pneu=-k*palt*palt+(k+1)*palt;
    y=pneu/YVAL*YHEIGHT;
    if(y<0) y=0;
    if(y>YHEIGHT) y=YHEIGHT;
    Draw(rp,x,y);
 }
 Move(rp,0,0); /* Winkelhalbierende zeichnen */
 Draw(rp,XWIDTH,YHEIGHT);
 Move(rp,0,YHEIGHT); /* Achsen einzeichnen */
 Draw(rp,0,0); Draw(rp,XWIDTH,0);
}

rech()
{
 short x,y;

 palt=p0;
 x=palt/XVAL*XWIDTH;
 if(x<0) x=0;
 if(x>XWIDTH) x=XWIDTH;
 Move(rp,x,0);
 for(y=x,n=0;n<nmax;++n)
 {
    Draw(rp,y,y);
    pneu=palt*(1+k*(1-palt));
    x=palt/XVAL*XWIDTH;
    y=pneu/YVAL*YHEIGHT;
    if(x<0) x=0;
    if(x>XWIDTH) x=XWIDTH;
    if(y>YHEIGHT) y=YHEIGHT;
    if(y<0) y=0;
    Draw(rp,x,y);
    palt=pneu;
 }
}