Appendix 1

/*                                            */
/*   Program ufmac to generate random         */
/*    iterative formulas and to compute       */
/*    their Julia set                         */
/*                                            */
#include <stdio.h>
#include <math.h>

main(argc,argv) unsigned int argc; char *argv[]; {          
 int Xscreen = 200,Yscreen=200;
 float Left= -4.,Right= 4.,Bottom= -4.,Top= 4.;
 int numcol = 256;
 int maxiter = 255; 
 float x,y,x0,y0,xd,yd,Deltax,Deltay;
 int i,j,ib,iter,ic,l,j1;
 float stack[50];
 int istack,k,numpull;
 int iexe[2][50];
 float p[2],p0[2];
/* nsel,isel and c contain the randomly generated formula*/
 int nsel[2],isel[2][50];
 float c[2][50];

/* Generate here random formula */

/* rndm() should return a random*/
/* float in the range 0-1 always*/
/* different every time the */
/* program is started */
 for(i=0; i < 2; i++) {
    nsel[i] = (int) (rndm()*49.+2.);
    for(j=0; j < nsel[i]; j++) {
       isel [i] [j] = (int) (rndm()*22. + 1.);
       c [i] [j] = rndm()*2. -1.;
       }
       }

/*Force first command to be a set accumulator */
/* command; not necessary but will speed up  */
/* formula computation                       */
 for(i=0; i < 2; i++) {
       isel [i] [0] = (int)(rndm()*3. + 13.);
       }
  


/* This is to speed up formula computation */
/* The accumulator value is pushed in the stack only  */
/* if there is a subsequent command that uses the stack */
/* iexe[i][j]=0 means that the command isel [i][j] must */
/* be ignored.                                          */
 for(i=0; i < 2; i++)
    for(j=0; j < nsel[i]; j++) iexe[i][j] = 1;
          for (j=0; j < 2; j++) {
            	for (i=0; i < nsel[j]; i++) {
                  if(isel[j][i] > 15 && isel[j][i] < 19) {
                 numpull = 0;
                 if((i+1) < nsel[j]) {
                 	for(k = i+1; k < nsel[j]; k++) {
                   if(isel[j][k] > 15 && isel[j][k] < 19) 
                    numpull = numpull - 1;
                   if(isel[j][k] > 18 && isel[j][k] < 23)
               numpull = numpull + 1;
                   if(numpull == 1) break;
                   }
                 }
                 if(numpull != 1) iexe[j][i] = 0;
                 }
                  }
                 }

/* Initialize here graphics */
/* Left,Right,Bottom,Top define the world coordinates,
   Xscreen,Yscreen the number of X and Y pixel */

/* Main loop to compute Julia set */
      Deltax = (Right - Left) / Xscreen;
      Deltay = (Top - Bottom) / Yscreen;
      for (j=0; j < Yscreen; j++) {
         for (l=0;l < Xscreen; l++) {
            x0 = (l-1)*Deltax + Left;
            y0 = (j-1)*Deltay + Bottom;
            x = x0 ;
            y = y0;
            iter = 1;
            while(iter <= maxiter){

/* Iteration starts here */
            p0[0] = x0;
            p0[1] = y0;
            for(i=0; i< 2; i++) {
            istack = 0;
            p[i] = p0[i];
            for (k=0; k < nsel[i]; k++) {
            	switch (isel[i][k]) {
   	    case  1: p[i] = p[i] + x0;break;
            case  2: p[i] = p[i] + y0;break;
            case  3: p[i] = p[i] + c[i][k];break;
            case  4: p[i] = p[i] - x0;break;
            case  5: p[i] = p[i] - y0;break;
            case  6: p[i] = p[i] - c[i][k];break;
            case  7: p[i] = p[i] * x0;break;
            case  8: p[i] = p[i] * y0;break;
            case  9: p[i] = p[i] * c[i][k];break;
            case 10: p[i] = p[i] / x0;break;
            case 11: p[i] = p[i] / y0;break;
            case 12: p[i] = p[i] /c[i][k];break;
            case 13: if(k == 0)p[i] = x0;break;
            case 14: if(k == 0)p[i] = y0;break;
            case 15: if(k == 0)p[i] = c[i][k];break;
            case 16: if(iexe[i][k] == 1){
              istack = istack + 1; stack[istack] = p[i];
              p[i] = x0; }break;
            case 17: if(iexe[i][k] == 1){
              istack = istack + 1; stack[istack] = p[i];
              p[i] = y0; }break;
            case 18: if(iexe[i][k] == 1){
              istack = istack + 1; stack[istack] = p[i];
                 p[i] = c[i][k]; }break;
            case 19: if(istack > 0) {p[i] = p[i] + stack[istack];
                     istack = istack - 1;}break;
            case 20: if(istack > 0) {p[i] = p[i] - stack[istack];
                     istack = istack - 1;}break;
            case 21: if(istack > 0) {p[i] = p[i] * stack[istack];
                     istack = istack - 1;}break;
            case 22: if(istack > 0) {
            	p[i] = p[i] / stack[istack];
                     istack = istack - 1;}break;
                }
                }
                }
            x = p[0];
            y = p[1];
/* end of iteration */

       if(fabs((double)x) < 1.e-10)x = 0.;
       if(fabs((double)y) < 1.e-10)y = 0.;
       if(fabs((double)x) > 1.e5 || fabs((double)y) > 1.e5)break;
       iter = iter + 1;
       x0 = x;
       y0 = y;
        }
       ic = (iter % numcol) + 1;
/*  Here set pixel l,j with color ic   */
       setc(&l,&j,&ic);
       }
       }
}