Table of Contents Editor's Forum ...................................... 1 Rand_c .............................................. 2 Cellular Automata ................................... 5 Iterated Functions Systems .......................... 23 Editor's Forum By Timothy Swenson This issue is another one that took some time to come together. I have spent part of the last two months beating my head against a single program. I have been working on a C version of my Cellular Automata program. After I finally found a decent random number generator, I made some changes to the generator that caused it to die. Looking in the documentation about the Small C floating point library, I could not find an answer. After a few dents in the old forehead, I mailed a copy of the program to Herb Schaaf. He found the problem in an instant. It seems that the documentation of the Small C library is WRONG. The details about the fmove() function are reversed. The two arguements should be switched to get the desired effects. After the change, the program worked (at least until I took out a few more bugs). Last issue I wrote about Software Patenting and an article in Dr. Doob's Journal. After I read the article, I wrote a letter to the editor (along with numberous others). In the June issue, I noticed that a letter to the editor sounded sort of familiar. Then it dawned on me, "Hey, I wrote that." It's sort of neat to have your name in a national publication. I'm not the only QL'er that has had the same experience. Herb Schaaf had a letter to the editor printed up in the C User's Journal. If you know of any other's (yourself included), please let me know. I'm sure most of you are members of some Sinclair group, but just in case you are not, you might be interested to know that there is a new publication/group in SinclairLand. Don Lambert and Bob Swoger have broken off from the lethargic SNUG group and formed a new TSNUG group. Don was bothered that the newsletter he put together for SNUG was not being distributed. He decided to grab the bull by the horns and page 1 start off on his own, with help from Bob. I suggest we (being QLer's) support the new group and it's energetic "motivator." The group should help to bring the various Sinclair groups together into a tighter community. As time goes on, our community is shrinking and we need to support each other into the coming years. The address for the new group/newsletter is: Don Lambert ZXir QLive Alive Newsletter 1301 Kiblinger Pl. Auburn, IN 46706 That's about all I have to say for now. Hope you enjoy the issue. Rand_c By Tim Swenson and Herb Schaaf In translating one of my pet programs from SuperBasic (first it started in Pascal), I needed a number of random numbers. Not finding a random number generator in Small-C, I decided to write the program in Lattice C. But, Lattice C does not really have any QL graphic functions. So back I go to Small-C to find a random number generator (RNG). I tried a number of RNGs but had no luck. I was limited by Small-C handling only 16 bit integers with a range from -32K to 32K. Herb Schaaf found that floating point numbers go beyond the 16 bit limitation. Soon, he sent me a RNG and a test program. I took the RNG and turned it into a seperate C function and file that can be included into your C programs. The seed for the RNG is taken from the date. This works well if you need a large number of random numbers. But if you run the program again, you will start with the same seed. If you need to have a program run numberous times, either build it into the code so it only initialises the RNG once, or have the user type in a seed at runtime. /* rand_c For QL Small-C with floating point library. Gets seed from date(). To make more random you could call initrand() before every occurance of rand(). USAGE: int num1[3], temp[3], num2; initrand(); * must call at least once * *temp = float(100); * 100 in float * fmove(rand(),num1); * put rand in num1 * *num1 = fmult(num1,temp); * num1 * 100 * num2 = int(num1); */ /* page 2 Global variables used for functions Small-C does not support Static Variables */ int rt, rnd[3], ra[3], rd[3], rm[3], rr[3]; initrand() { int t[2]; *rm = float(125); *rd = atof("2796203"); *t = date(); *ra = float(t[1]); } rand() { *ra = fmult(ra,rm); *ra = fdiv(ra,rd); rt = int(ra); *rr = float(rt); *rnd = fsub(ra,rr); *ra = fmult(rnd,rd); return rnd; } Herb, not one to say "good enough", has also come up with another RNG. This one is based on the algorithm that the ZX81 and T/S 2068 use in thier RNG's. This program is not set up to be included in your program, but with some minor work can be made to do so. /* ZXrand_c April 26, 1991 11am * Using QL Small-C to * Generate Random Numbers insame way as ZX-81 & TS 2068 */ #include int a[3], b[3], m[3], d[3]; int seed[3], random[3], i; char c, strg[6], *ow; main() { initz(); printf("Choose an integer seed in the range from 0 to 65535 \n"); gets(strg); *seed = atof(strg); printf("The seed will be %s\n",ftoa(seed,strg)); printf("touch [SPACE BAR] to begin, q to quit\n"); while((getc(c)) != 'q') { rand(); printf("%6s",ftoa(seed,strg)); } printf("\nPress any key to exit\n"); beep(10,1); for(i=0; i<20480; i++) { ; page 3 } beep(10,2); getc(); } initz() { *a = atof("65536"); *b = atof("65537"); *m = atof("75"); *d = atof("74"); } rand() { *random = fmult(m,seed); *random = fadd(random,d); *random = fmod(); fmove(seed,random); } fmod() { while(fcmp(random,b) >0 ) { *random = fsub(random,b); } return (*random); } _console() { mode(4); fopen(ow,"con","w"); window(ow,512,256,0,0); paper(ow,0); ink(ow,4); border(ow,7,2); cls(ow); } /* End of ZXrand_c Apr 26, 1991 noon ** */ Now that I know how to beat the 16 bit limit with floating point numbers, I have a few other RNGs that could be ported to Small-C. Plus, Herb still has a few more up his sleeves. If there is interest, I will publish these other RNGs. Cellular Automata By Timothy Swenson Back in college, I was short 1/2 a unit from graduating. To make up the unit, I enrolled in a 1 unit Independent Study course. With hours left before the enrollment deadline, I picked what I thought was a weak but interesting subject, Cellular Automata (CA). Well, the more I researched the topic, the more interesting the class became, and the more involved I became. The end result of the class was a paper on CA and a demonstration program. I originally wrote the program in Pascal, because Pascal was the language the college used in most of the undergrad courses. Later I translated the program to SuperBasic so I could utilize the full power of page 4 the QL. Since it ran so slow in SuperBasic, I decided to write it in C for speed. Before I go into the program, let me first explain what Cellular Automata is. CA is comprised of four items; an array (Cellular Space), a number of values for each cell in the array (States), a Neighborhood relation that defines how cells relate to each other, and a Transition Rule that defines how cells change from one state to another over time. One popular version of CA, is the Game of Life by John Conway. Here you have a 2 dimensional CA with 2 states and an 8 cell neighborhood. The Quanta library has a good version of this game (although it's not considered a game by many). I focused on 1 dimensional CA (1DCA). This is where the array is only 1 dimensional and the array is displayed over time by plotting each new array below the previous one. This makes for some rather interesting displays. I could to on for pages about CA, but I hope you have the general idea. If you would like to know more, let me know and I'll send you a copy of my college paper and copies of my references. The SuperBasic version of the program is the most complete. The C version is short and sweet, due to my lack of knowing how to do many things in C. The Pascal version is a little more advanced over the C version. I'm not too sure if the Pascal version is working 100%. My Pascal compiler is on long-term loan to a friend so I have not had the chance to compile it and try it out. It is written for the Metacomco Pascal compiler. 100 LET xmax = 512 110 alldone = 0 120 DIM world(xmax) 130 DIM rule(64,9,2) 140 WINDOW #1,512,256,0,0 150 PAPER #1,0 : INK #1,4 160 LET rule1 = 1 : LET cycle = 15 170 LET bloksiz = 2 180 IF alldone = 1 THEN STOP 190 CLS #1 200 CSIZE 2,1 210 PRINT " ONE DIMENSIONAL" 220 PRINT " CELLULAR AUTOMATA"\\ 230 CSIZE 0,0 240 PRINT \\ 250 LET tab = 20 260 PRINT TO tab;" M E N U " 270 PRINT TO tab;" --------------"\\ 280 PRINT TO tab;" 1] Load Rule File" 290 PRINT TO tab;" 2] Choose Rule To Use" 300 PRINT TO tab;" 3] Display Rule" 310 PRINT TO tab;" 4] Set Number of Cycles" 320 PRINT TO tab;" 5] Set Block Size" 330 PRINT TO tab;" 6] Run Automata" 340 PRINT TO tab;" 7] Instructions"\\ 350 INPUT TO tab;" Enter Choice : ";menu 360 SELect ON menu 370 ON menu = 1 page 5 380 loadrule 390 ON menu = 2 400 chooserule 410 ON menu = 3 420 displayrule 430 ON menu = 4 440 setcycles 450 ON menu = 5 460 setblock 470 ON menu = 6 480 startautomata 490 ON menu = 7 500 instruct 510 END SELect 520 GO TO 180 530 REMark ************************ 540 DEFine PROCedure loadrule 550 LOCal a$,b$,xx,yy 560 CLS #1 570 PRINT\\\\\\\\\\ 580 PRINT "Enter File Name to Load. Include Device Name (i.e MDV1_)"\\ 590 INPUT a$ 600 OPEN_IN #4,a$ 610 INPUT #4,b$ 620 INPUT #4,neighbor 630 INPUT #4,number 640 FOR xx = 1 TO number 650 INPUT #4,b$ 660 INPUT #4,b$ 670 FOR yy = 1 TO neighbor+1 680 rule(xx,yy,1) = b$(yy) 690 NEXT yy 700 INPUT #4,b$ 710 FOR yy = 1 TO neighbor+1 720 rule(xx,yy,2) = b$(yy) 730 NEXT yy 740 NEXT xx 750 CLOSE #4 760 END DEFine 770 REMark ********************** 780 DEFine PROCedure chooserule 790 CLS #1 800 PRINT\\\\\\\\\\ 810 PRINT "Enter Rule to Use. Enter the Nth Rule in the File." 820 INPUT rule1 830 PAUSE 20 840 END DEFine 850 REMark ********************* 860 DEFine PROCedure displayrule 870 LOCal xx 880 CLS #1 890 PRINT \\\\\\\ : PRINT TO 28;"R U L E ";rule1 900 PRINT TO 27;"-------------" 910 PRINT\ : PRINT TO 28; 920 FOR xx = 1 TO neighbor+1 930 PRINT rule(rule1,xx,1);" "; 940 NEXT xx 950 PRINT : PRINT TO 28; 960 FOR xx = 1 TO neighbor+1 970 PRINT rule(rule1,xx,2);" "; 980 NEXT xx page 6 990 PRINT \\\\\ 1000 PRINT TO 20;"Enter Any Key to Continue" 1010 PAUSE 1020 END DEFine 1030 REMark ************************* 1040 DEFine PROCedure setcycles 1050 CLS #1 1060 PRINT\\\\\\\\\\ 1070 PRINT TO 15;"Enter Number of Cycles to Run Through" 1080 PRINT TO 15;" the Program Pauses"\\ 1090 INPUT cycle 1100 PAUSE 10 1110 END DEFine 1120 REMark ********************* 1130 DEFine PROCedure setblock 1140 CLS #1 1150 PRINT\\\\\\\\\\ 1160 PRINT TO 15;"Enter the Size of the Blocks to Use." 1170 PRINT TO 15;" ( 1 to 8 )"\\ 1180 INPUT bloksiz 1190 PAUSE 20 1200 END DEFine 1210 REMark ************************* 1220 DEFine PROCedure instruct 1230 CLS #1 : PRINT \\\\\ 1240 PRINT TO 10;" This program will run a cellular system through a specified rule." 1250 PRINT TO 10;"The program lets you choose the size of the blocks to display." 1260 PRINT TO 10;"This allows you to see in greater detail what is going on. Once" 1270 PRINT TO 10;"the program pauses, you may hit F1 to start to program over, or" 1280 PRINT TO 10;"you may hit F2 to stop the program. The program does not clear" 1290 PRINT TO 10;"the screen so that you may use a screen dump program to print " 1300 PRINT TO 10;"out the results. This will give you a hard copy for later ref-" 1310 PRINT TO 10;"erence." 1320 PRINT \\\\ : PRINT TO 30;"Hit Any Key To Continue" 1330 PAUSE 1340 END DEFine 1350 REMark ************************* 1360 DEFine PROCedure startautomata 1370 setinit 1380 count = 15 1390 CLS #1 1400 PRINT TO 15;"Neighborhood ";neighbor;" Rule ";rule1 1410 display 1420 up = -bloksiz 1430 REPeat done 1440 FOR loop = 1 TO cycle 1450 count = count + bloksiz 1460 IF count > 250 THEN 1470 count = 250 1480 SCROLL up,0 1490 END IF 1500 nextcycle 1510 display 1520 NEXT loop 1530 BEEP 150,4 page 7 1540 PAUSE 1550 LET key = KEYROW(0) 1560 IF key = 2 THEN EXIT done 1570 IF key = 8 THEN 1580 LET alldone = 1 1590 EXIT done 1600 END IF 1610 END REPeat done 1620 END DEFine 1630 REMark ************************ 1640 DEFine PROCedure getdata 1650 PRINT\\\\\\\\\\ 1660 PRINT TO 15;"Enter File Name of Set-up File to Load." 1670 PRINT TO 15;"Include Drive Name (i.e MDV1_ )"\\ 1680 INPUT a$ 1690 OPEN_IN #4,a$ 1700 INPUT #4,b$ 1710 LET world(b$)=1 1720 IF EOF(#4) THEN 1730 CLOSE #4 1740 RETurn 1750 ELSE 1760 GO TO 1700 1770 END IF 1780 END DEFine 1790 REMark ********************* 1800 DEFine PROCedure initialise1 1810 LOCal z, rand 1820 world(1) = 0 : world(xmax) = 0 1830 FOR z = 2 TO INT(xmax/bloksiz)-1 1840 LET rand = RND(1 TO 10) 1850 IF rand <4 THEN 1860 world(z) = 1 1870 ELSE 1880 world(z) = 0 1890 END IF 1900 NEXT z 1910 END DEFine 1920 REMark ************************ 1930 DEFine PROCedure initialise2 1940 LOCal xx 1950 FOR xx = 1 TO xmax 1960 world(xx) = 0 1970 NEXT xx 1980 END DEFine 1990 REMark *********************** 2000 DEFine PROCedure nextcycle 2010 LOCal start, finish, x, number 2020 start = neighbor/2 : finish=INT(xmax/bloksiz)-start 2030 FOR x = start TO finish 2040 number = ruleeval 2050 IF (number=1) AND (world(x)=0) THEN world(x)=2 2060 IF (number=0) AND (world(x)=1) THEN world(x)=3 2070 NEXT x 2080 END DEFine 2090 REMark ***************************** 2100 DEFine PROCedure display 2110 LOCal xx, xxx 2120 FOR xx = 1 TO INT(xmax/bloksiz)-1 2130 xxx = xx * bloksiz 2140 IF (world(xx)=1) OR (world(xx)=2) THEN BLOCK bloksiz,bloksiz,xxx,count,4 2150 IF world(xx) = 2 THEN world(xx) = 1 page 8 2160 IF world(xx) = 3 THEN world(xx) = 0 2170 NEXT xx 2180 END DEFine 2190 REMark ************************** 2200 DEFine FuNction neighborhood2 2210 LOCal count 2220 LET count = 0 2230 IF (world(x-1)=1) OR (world(x-1)=3) THEN count=count+1 2240 IF (world(x+1)=1) OR (world(x+1)=3) THEN count=count+1 2250 RETurn count 2260 END DEFine 2270 REMark ************************* 2280 DEFine FuNction neighborhood4 2290 LOCal count 2300 LET count = 0 2310 IF (world(x-2)=1) OR (world(x-2)=3) THEN count=count+1 2320 IF (world(x-1)=1) OR (world(x-1)=3) THEN count=count+1 2330 IF (world(x+1)=1) OR (world(x+1)=3) THEN count=count+1 2340 IF (world(x+2)=1) OR (world(x+2)=3) THEN count=count+1 2350 RETurn count 2360 END DEFine 2370 REMark ************************* 2380 DEFine FuNction neighborhood6 2390 LOCal count 2400 LET count = 0 2410 IF (world(x-3)=1) OR (world(x-3)=3) THEN count=count+1 2420 IF (world(x-2)=1) OR (world(x-2)=3) THEN count=count+1 2430 IF (world(x-1)=1) OR (world(x-1)=3) THEN count=count+1 2440 IF (world(x+1)=1) OR (world(x+1)=3) THEN count=count+1 2450 IF (world(x+2)=1) OR (world(x+2)=3) THEN count=count+1 2460 IF (world(x+3)=1) OR (world(x+3)=3) THEN count=count+1 2470 RETurn count 2480 END DEFine 2490 REMark ************************* 2500 DEFine FuNction ruleeval 2510 LOCal number, tempx 2520 SELect ON neighbor 2530 ON neighbor = 2 2540 number = neighborhood2 2550 ON neighbor = 4 2560 number = neighborhood4 2570 ON neighbor = 6 2580 number = neighborhood6 2590 END SELect 2600 tempx = world(x) 2610 RETurn rule(rule1, number+1, tempx+1) 2620 END DEFine 2630 REMark ************************* 2640 DEFine PROCedure setinit 2650 CLS #1 2660 PRINT\\\\\\ 2670 PRINT TO 20;"Method of Initialisation" 2680 PRINT TO 20;"------------------------"\\ 2690 PRINT TO 20;" 1] Load From a File" 2700 PRINT TO 20;" 2] Random "\\ 2710 PRINT TO 20;"Enter Way to Initialise"\\ 2720 INPUT init 2730 IF init = 1 THEN 2740 CLS #1 2750 initialise2 2760 getdata 2770 END IF 2780 IF init = 2 THEN initialise1 page 9 2790 END DEFine Example Rule file for Neighborhood of 2. Totalalistic Rules 2 7 Rule 1 100 100 Rule 2 010 010 Rule 3 110 110 Rule 4 001 001 Rule 5 101 101 Rule 6 011 011 Rule 7 111 111 /* cell_c One Dimensional Cellular Automata (1DCA) This program is used to study 1DCA. It will only evaluate one rule at a time. It can be setup to start processing on a particular array or you may choose a random array. This program requires a rule file and optional an array file. The formats are below: Rule File: two lines of 1's or 0's, length equal to MAX_NEIGHBORHOOD. Sample: 011010 - Alive 001100 - Dead */ #include #include "flp2_rand_c" #define MAX_WORLD 516 #define MAX_NEIGHBORHOOD 5 #define MAX_Y 256 #define EOL LF page 10 int world_array[MAX_WORLD], alive_array[MAX_NEIGHBORHOOD], dead_array[MAX_NEIGHBORHOOD], y, fd; main() { int i; y = 0; printf("One Dimensional Cellular Automata\n"); printf(" Analysis Program\n"); printf("\n"); printf(" by Timothy Swenson\n"); printf("\n"); initrand(); load_rule(); init_array(); fd = fopen("scr_512x256a0x0","w"); if (fd == NULL) { printf("Error Opening Screen\n"); abort(1); } cls(fd); display(); ++y; for (i = 0; i < MAX_Y; ++i) { next_cycle(); display(); ++y; } } /*******************************************/ /* This procedure loads the rule from a */ /* file into two arrays; alive_array and */ /* dead_array. */ /*******************************************/ load_rule() { int rule_fd, count, c; char string[30]; puts("Enter Rule File to Load\n"); gets(string); rule_fd = fopen(string,"r"); if (rule_fd == NULL) { printf("Error Loading File"); abort(1); } count = 0; while ( (c = fgetc(rule_fd)) != EOL) { alive_array[count] = c-'0'; count++; } count = 0; page 11 while ( (c = fgetc(rule_fd)) != EOL) { dead_array[count] = c-'0'; count++; } fclose(rule_fd); } /************************************/ /* This procedures initalises the */ /* array randomly */ /************************************/ init_array() { int i,temp1[3], temp2[3]; char c; *temp2 = atof("0.45"); printf("Creating Random Array\n"); for (i=2; i 0) world_array[i] = 1; else world_array[i] = 0; } } /********************************/ /* This procedures diplays the */ /* current array below the */ /* previous array. */ /********************************/ display() { int i; for (i=2; i < MAX_WORLD-2; ++i) { if (world_array[i] == 1 || world_array[i] == 2) { world_array[i] = 1; block(fd,1,1,i-2,y,4); } else world_array[i] = 0; } } /***********************************/ /* This procedures evalutates the */ /* next cycle. */ /***********************************/ next_cycle() { int i, temp, result; for (i=2; i 250 THEN BEGIN COUNT := 250; ERR := SCREEN2(SCREENSCROLL,UP); END; NEXTCYCLE(RULENUM, NEIGHBOR, BLOKSIZ); DISPLAY(COUNT, BLOKSIZ); END; END; {@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@} {@@@@@@ M A I N P R O G R A M @@@@@@} {@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@} BEGIN MAXX := XMAX+1; MAXY :=256; ZERO := 0; CYCLE := 10; RULENUM := 1; BLOKSIZ := 1; ERR := WINDOW7(WINDOWDEFINE,MAXX,MAXY,ZERO,ZERO,ZERO,ZERO); ALLDONE := FALSE; WHILE NOT ALLDONE DO BEGIN ERR := SCREEN1(SCREENCLEAR); WRITELN(' C E L L U L A R A U T O M A T A'); WRITELN(' - O N E D I M E N S I O N A L - '); WRITELN; WRITELN(' By Timothy Swenson '); WRITELN; WRITELN(' When the program pauses between a number of cycles '); WRITELN(' you may hit any key to go on. Hit F1 to restart '); WRITELN(' the program. Hit F2 to end the program.'); WRITELN; WRITELN; WRITELN(' M E N U'); WRITELN; WRITELN(' 1] Load Rule File'); WRITELN(' 2] Choose Rule to use'); WRITELN(' 3] Set number of Cycles'); WRITELN(' 4] Set Block Size'); WRITELN(' 5] Run Automata '); WRITELN; WRITE(' Enter Choice : '); READLN( MENU ); CASE MENU OF 1 : GETRULES( NEIGHBOR ); 2 : BEGIN WRITELN; WRITELN; WRITELN(' Enter a rule to use: by number only '); page 20 WRITELN(' The Nth rule in the Rule file is the Nth Rule '); WRITELN(' and will be entered as such. '); READLN ( RULENUM ); END; 3 : BEGIN WRITELN; WRITELN('How many cycles to run through before a pause.'); READLN( CYCLE ); END; 4 : BEGIN WRITELN; WRITELN(' What size of blocks do you want for the display? '); WRITELN(' 1 is the smallest and 8 is the largest suggested '); READLN(BLOKSIZ); END; 5 : BEGIN SETINIT; ERR := SCREEN1(SCREENCLEAR); DONE := FALSE; COUNT := 15; WRITELN(' R U L E ',RULENUM); DISPLAY(COUNT, BLOKSIZ); UP := -BLOKSIZ; WHILE NOT DONE DO BEGIN MAINLOOP; KEY := READKEY; { once the program has cyclesd through the main loop for the desired number of loops, the program waits for a user input for what it should do next } IF KEY = 1256 THEN DONE := TRUE; { start the program again } IF KEY = 1260 THEN BEGIN { end the progam } ALLDONE := TRUE; DONE := TRUE; END; { If } END; { While done } END; { Case of 5} END; { Case total } END; { while alldone } { The program does not clear the screen. This leaves the results on the screen and ready to be dumped to the printer. This way I have a hard copy of the results. } END. Iterated Function Systems By Herb Schaaf IFS is what Michael Barnsley calls his Iterated Function Systems, and IFS is another part of fractals (remember those Mandelbrot all-nighters ?). Barnsley did IFS a few years ago and has gone on to bigger and better stuff now. He explains how IFS can be done by anyone with a personal computer using page 21 a "collage" approach. Ready for a collage education? Read on. What's the big idea ? - data compression! with the IFS algorithm and just a few well chosen coefficients, you can represent a very complicated picture that otherwise would require many kilobytes to describe. In other words 'one picture is worth a thousand words' but an IFS only needs a dozen or so coefficients. The concept is that a bigger picture is composed of self-similar smaller pictures that just happen to cover the bigger picture nicely. For example, a checkerboard is just a big square that can be neatly covered with 64 smaller squares. Actually we only need 32, since half are blank. You get the idea. Now it can be more complicated, since we can also use rectangles, which can be considered as self-similar to squares where we have changed the scales along the x and/or y axes. Not only that, but we are also allowed to tilt and rotate the little pieces to help us get a collage with better coverage. Not too useful for squares, but real handy for more organic swirly 'natural' shapes like Barnsley's popular fern patterns. Once you have the big outline, and then get a covering of it with a collage of smaller self-similar shapes, you can work out the coefficients. A convenient way is to use three point affine transformations to develop sets of simultaneous equations. Solving these three at a time produces the coefficients needed by the IFS. You can then assign probability values to each transformation set based on how dark you want it to be in the final picture. With the coefficients and probabilities in hand, the drudge work begins, such as calculating and plotting 10,000 points where the coordinates of each point depends on the coordinates of the previous point, (with a little randomness thrown in). This tedious task is the kind of thing computers do well. It takes about 10 minutes for the QL to do 10000 points. As you watch and wait you will see a nebulous scatter of points take on shape and substance. I've written a SuperBasic program to do 2-dimensional IFS, copying some of the better examples from the literature and working out a few on my own. I've also programmmed IFS in Small-C. Barnsley has gone on to 3-D and color, and real time HDTV encoding/decoding etc. with dedicated microchips. One of Barnsley's pet IFS pictures is a fernlike image, a complicated structure that is encoded with only a few coefficients. If you'd like to try it on your computer the data given by Barnsley for a black spleenwort fern is: Translation Rotation Scale Probability 0 0 0 0 0 .16 .005 0 1.6 -2.5 -2.5 .85 .85 .8 0 1.6 49 49 .3 .34 .0975 0 .44 120 -50 .3 .37 .0975 page 22 We can convert these into a simpler form by using trigonometric functions to put the rotations into their equivalent x,y shifts. After doing that we end up with a table (from BYTE Jan '88) such as: W a b c d e f p 1 0 0 0 .16 0 0 .01 2 .2 -.26 .23 .22 0 1.6 .07 3 -.15 .28 .26 .24 0 .44 .07 4 .85 .04 -.04 .85 0 1.6 .85 These values are used in the SuperBasic program 'IFS4_bas' and can be entered at the prompts in Small-C IFS4_c. The algorithm is said to be 'robust', and small changes in the values only change the final result slightly. Enjoy the show by changing one or more of them and allow about 10 minutes (or slightly less using the Small-C version) for the program to do its thing. If there is interest we have coefficients for a few other figures, such as snowflakes, leaves, curlicues, galaxies, pentagons, dragon curves, etc. References: BYTE January and April 1988 Algorithms May/June 1990 Scientific American August 1990 Electronics & Wireless World August 1990 There are a series of books out by Robert L. Devaney, some at the popular level, and others for hard-core mathematicians. He also made some beautiful educational videos. I recommend the books and videos to those who are curious about things like fractals, chaos, Mandelbrot and Julia sets, strange attractors, etc. 100 REMark IFS4_bas Iterative Function System 110 REMark  H. L. Schaaf June 3, 1991 120 n$="fern" :REMark change to suit 130 REMark set up to do fern 140 RESTORE 150 WMON : MODE 4:WINDOW #1,512,256,0,0 160 PAPER #1,4: INK #1,0:CLS#0 170 PRINT "probabilities", 180 REMark m is the number of transformations 190 m=4 : REMark change to suit, fern only uses 4 200 DIM a(m), b(m), c(m), d(m), e(m), f(m), p(m) 210 DATA 0, 0, 0, .16, 0, 0, 1E-2 220 DATA .2, -.26, .23, .22, 0, 1.6, 7E-2 230 DATA -.15, .28, .26, .24, 0, .44, 7E-2 240 DATA .85, 4E-2, -4E-2, .85, 0, 1.6, .85 250 pt=0 :REMark total of probabilities 260 FOR j=1 TO m 270 READ a(j), b(j), c(j), d(j), e(j), f(j), p(j) 280 pt = pt + p(j) 290 END FOR j 300 FOR j = 1 TO m: p(j)= p(j)/pt :PRINT p(j),:END FOR j 310 pt=0 320 FOR i = 1 TO m 330 pt=pt+p(i) : p(i) = pt 340 PRINT\ a(i),b(i),c(i),d(i),e(i),f(i),p(i) 350 END FOR i 360 x=0 : y=0 : maxy=-1E6 : maxx=-1E6 : minx=1E6 : miny=1E6 page 23 370 FOR n = 1 TO 100 380 AT 5,0: PRINT n, 390 pk = RND 400 PRINT pk, 410 k = 1 420 REPeat loopa 430 IF pk <= p(k) THEN EXIT loopa 440 k = k+1 450 END REPeat loopa 460 xnxt=a(k)*x+b(k)*y+e(k) 470 ynxt=c(k)*x+d(k)*y+f(k) 480 x=xnxt : y=ynxt 490 IF y > maxy THEN maxy = y 500 IF x > maxx THEN maxx = x 510 IF x < minx THEN minx = x 520 IF y < miny THEN miny = y 530 PRINT x,y\minx,maxx,miny,maxy 540 END FOR n 550 CLS 560 yscale = 1.1* ABS(maxy-miny) 570 xscale = 1.1* ABS(maxx-minx) 580 PRINT "X & Y & P scales = ",xscale,yscale, 590 IF xscale>yscale : pscale=xscale :ELSE pscale = yscale 600 PRINT pscale\ 610 xoffset = minx: yoffset = miny 620 xfill = ABS(maxx-minx) : xmarg=(1.3*pscale)-xfill 630 PRINT "X range & margin ",xfill,xmarg\ 640 yfill = ABS(maxy-miny) : ymarg=pscale-yfill 650 PRINT "Y range & margin ",yfill,ymarg\ 660 xoffset = minx -.5*xmarg 670 yoffset = miny -.5*ymarg 680 PRINT "Scale = ";pscale;" lower left (x,y)= "; xoffset, yoffset 690 CLS#0: PRINT#0;"any key to begin 5000 iterations"; " = +/- 5 minutes" 700 PAUSE 710 CLS 720 SCALE pscale,xoffset,yoffset 730 iters = 5000 :REMark 5000 iterations change to suit 740 FOR n = 1 TO iters 750 pk = RND 760 k = 1 770 REPeat loopb 780 IF pk<=p(k) THEN EXIT loopb 790 k = k+1 800 END REPeat loopb 810 xnxt=a(k)*x+b(k)*y+e(k) 820 ynxt=c(k)*x+d(k)*y+f(k) 830 x = xnxt : y = ynxt 840 POINT x,y 850 END FOR n 860 BEEP 30000,2 :REMark lets you know it's done 870 SBYTES "ram1_"&n$&"_pix",131072,32768 880 REMark saves it 'clean' 890 PRINT #0\\"Clean screen in ram1_"&n$&"_pix" 900 REMark you can save it to flp, etc. /* IFS4b_c H.L.Schaaf June 3, 1991 * Iterated Function Systems */ #include page 24 int a[21], b[21], c[21], d[21], e[21], f[21], p[21], ps[21], pk[3], fpw[3], fpx[3], fpy[3], fpz[3], rnd[3], ra[3], rd[3], rm[3], rr[3], det[3], tmp[3], w[3], x[3], y[3], z[3], xs[3], ys[3], maxx[3], maxy[3], minx[3], miny[3], yscale[3], xscale[3], pscale[3], xoffset[3], yoffset[3], i, j, k, n, t, iter; char *ow, ch, str[16]; main() { initrand(); prompt("choose number of transforms"); choose_n(); prompt("enter array values"); selarray(); prompt("determine raw probabilities"); rawprob(); prompt("equalize probabilities, sum = 1"); sumprob(); prompt("get ranges of graphics values"); sizeup(); prompt("see IFS develop, 1/2 minute for K points"); forshow(x,y); at(ow,0,0); printf("any key to exit"); getchar(); exit(); } getsit() { printf("?"); gets(str); *tmp = atof(str); return tmp; } selarray() { csize(ow,0,0); cls(ow); printf("\n%7c%14c%14c%14c%14c%14c\n" ,'a','b','c','d','e','f'); for ( j = 0; j < n; j++ ) { k = 2 * j; at(ow,( 5 + k ), 0 ); printf("%u",j); } for ( j = 0; j < n; ++j ) { k = 2 * j; at(ow,( 5 + k),(7+(0*14))); getsit(); fputa(a,j,tmp); at(ow,(5 + k),6); fpica(a,j,tmp); printf("%s \n",ftoa(tmp,str)); at(ow,( 5 + k),(7+(1*14))); getsit(); fputa(b,j,tmp); at(ow,(5 + k),20); fpica(b,j,tmp); printf("%s \n",ftoa(tmp,str)); at(ow,( 5 + k),(7+(2*14))); getsit(); fputa(c,j,tmp); at(ow,(5 + k),34); fpica(c,j,tmp); printf("%s \n",ftoa(tmp,str)); at(ow,( 5 + k),(7+(3*14))); getsit(); fputa(d,j,tmp); at(ow,(5 + k),48); fpica(d,j,tmp); printf("%s \n",ftoa(tmp,str)); at(ow,( 5 + k),(7+(4*14))); getsit(); fputa(e,j,tmp); at(ow,(5 + k),62); fpica(e,j,tmp); printf("%s \n",ftoa(tmp,str)); at(ow,( 5 + k),(7+(5*14))); getsit(); fputa(f,j,tmp); at(ow,(5 + k),76); fpica(f,j,tmp); printf("%s \n",ftoa(tmp,str)); } printf("\nIs this OK ? [y/n]"); if((ch = getchar()) == 'n') { selarray(); } else { csize(ow,1,0); } } page 25 askmore() { at(ow,24,0); printf("%uK want another K ?",iter); ch = getchar(); if(ch != 'n') { at(ow,24,0); printf(" "); doakay(); } } forshow(x,y) int x[3], y[3]; { scale(ow,pscale,xoffset,yoffset); ink(ow,6); cls(ow); point(ow,maxx,maxy); point(ow,maxx,miny); point(ow,minx,maxy); point(ow,minx,miny); ink(ow,4); iter = 0; doakay(); } doakay() { for ( i = 0 ; i < 1024 ; i++) { iterate(x,y); point(ow,x,y); } iter ++; askmore(); } iterate(x,y) int x[3], y[3]; { k = selectk(); fpica(a,k,fpw); fpica(b,k,fpx); fpica(c,k,fpy); fpica(d,k,fpz); *xs = fmult(fpw,x); *tmp = fmult(fpx,y); *xs = fadd(xs,tmp); fpica(e,k,fpw); *xs = fadd(xs,fpw); *ys = fmult(fpy,x); *tmp = fmult(fpz,y); *ys = fadd(ys,tmp); fpica(f,k,fpw); *ys = fadd(ys,fpw); fmove(x,xs); fmove(y,ys); return (x,y); } selectk() { fmove(pk,rand()); for(k=0; k0) { fmove(maxy,y); } if (fcmp(x,maxx) >0) { fmove(maxx,x); } if (fcmp(y,miny) <0) { fmove(miny,y); } if (fcmp(x,minx) <0) { fmove(minx,x); } at(stdout,10,31); printf("%u \n",i); at(stdout,10,8); printf("X min = %s ", ftoa(minx,str)); at(stdout,10,36); printf("X max = %s \n",ftoa(maxx,str)); at(stdout,16,25); printf("Y min = %s ", ftoa(miny,str)); at(stdout,4,25); printf("Y max = %s \n",ftoa(maxy,str)); page 26 } trace("BEGINNING GRAPHICS SCALING"); *w = atof("1.2"); /* 10% margin */ *tmp = fsub(maxy,miny); *tmp = fabs(tmp); *yscale = fmult(tmp,w); *tmp = fsub(maxx,minx); *tmp = fabs(tmp); *xscale = fmult(tmp,w); if(fcmp(xscale,yscale) >0 ) fmove(pscale,xscale); else fmove(pscale,yscale); fmove(xoffset,minx); fmove(yoffset,miny); *fpx = fsub(maxx,minx); *fpx = fabs(fpx); *fpy = fsub(maxy,miny); *fpy = fabs(fpy); *w = atof("1.3"); *fpw = fmult(pscale,w); *fpw = fsub(fpw,fpx); *fpz = fsub(pscale,fpy); *w = atof(".5"); *fpw = fmult(w,fpw); *fpz = fmult(w,fpz); *xoffset = fsub(minx,fpw); *yoffset = fsub(miny,fpz); trace("END OF SCALING\n\n"); trace("Scaling factors based on 200 iterations:"); printf("pscale = %s \n",ftoa(pscale,str)); printf("xoffset = %s \n",ftoa(xoffset,str)); printf("yoffset = %s \n",ftoa(yoffset,str)); printf("x = %s \n",ftoa(x,str)); printf("y = %s \n",ftoa(y,str)); scale(stdout,pscale,xoffset,yoffset); trace("Ready to do a IFS\n"); return(pscale,xoffset,yoffset); } sumprob() { z[0] = float(0); for( i = 0; i < n; i++) { fpica(p,i,fpw); *tmp = fdiv(fpw,pk); *z = fadd(z,tmp); printf("cumulative probability at %u is %s \n" ,i,ftoa(z,str)); fputa(ps,i,z); fputa(p,i,tmp); } printf("sum of probabilities = %s\n",ftoa(z,str)); } rawprob() { *pk = float(0); for( i=0; i(n-1)) { printf(" Which transform 0 to %u ? ",n-1); fgets(ch,1,stdin); iter = atoi(ch); printf(" probability for %u ",iter); getsit(); fputa(p,iter,tmp); page 27 printf(" = %12s \n",ftoa(tmp,str)); } printf("\n Review of probabilities \n\n"); *pk = float(0); for( i = 0; i < n; i++) { fpica(p,i,fpw); printf(" transform %u has probability of %12s\n" ,i,ftoa(fpw,str)); *pk = fadd(pk,fpw); } printf("\nChange probability value(s) ? [y/n]"); } } fputa(array,element,value) int array[21], element, value[3]; { int i; element *= 3; for(i = 0; i < 3 ; i++) { array[i+element] = value[i]; } } fpica(array,element,value) int array[21], element, value[3]; { int i; element *= 3; for(i = 0; i < 3 ; i++) { value[i] = array[i+element]; } } rand() { *ra = fmult(ra,rm); *ra = fdiv(ra,rd); t = int(ra); *rr = float(t); *rnd = fsub(ra,rr); *ra = fmult(rnd,rd); return rnd ; } choose_n() { cls(stdout); n = 0; while( n<2 || n >7 ) { trace("Choose a number between 2 and 6\n"); fgets(ch,1,stdin); n = atoi(ch); printf("%2u was chosen \n",n); } return(n); } prompt(str) char str[]; { printf("\nTouch [ SPACE BAR ] to %s\n",str); getchar(); } trace(str) char str[]; { printf("\ntrace %s\n",str); } initrand() { *rm = float(125); *rd = atof("2796203"); *ra = atof("100001"); for( i = 0; i<21; i++) { a[i]= b[i]= c[i]= d[i]= e[i]= f[i]= p[i]= ps[i]= 0; } } _console() { mode(4); fopen(ow,"con","w"); window(ow,512,256,0,0); paper(ow,0); ink(ow,4); border(ow,1,2); cls(ow); } page 28 /* JUNE 3 1991 10pm * */ page 29