100 REM ITERATED FUNCTION SYSTEMS Ver 1.0 (c) 1993 by J. C. Sprott 110 DEFDBL A-Z: DIM A(12) 120 RANDOMIZE TIMER 'Reseed random number generator 130 SCREEN 12 'Assume VGA graphics 140 GOSUB 300 'Set parameters 150 GOSUB 400 'Iterate equations 160 GOSUB 500 'Display results 170 GOSUB 600 'Test results 180 ON T% GOTO 130, 140, 150 190 END 300 REM Set parameters 310 X = .05: Y = .05: XE = X + .00001: YE = Y 320 FOR I% = 1 TO 12: A(I%) = .1 * (INT(25 * RND) - 12): NEXT I% 330 J0 = ABS(A(1) * A(4) - A(2) * A(3)) 340 J1 = ABS(A(7) * A(10) - A(8) * A(9)) 350 IF J0 + J1 = 0 OR J0 > 1 OR J1 > 1 THEN GOTO 320 'Not contracting 360 P = J0 / (J0 + J1) 370 T% = 3: LSUM = 0: N = 1: N1 = 0: N2 = 0 380 XMIN = 1000000#: XMAX = -XMIN: YMIN = XMIN: YMAX = XMAX 390 RETURN 400 REM Iterate equations 410 IF X <> XE THEN IF RND > P THEN R% = 6 ELSE R% = 0 420 XNEW = A(1 + R%) * X + A(2 + R%) * Y + A(5 + R%) 430 YNEW = A(3 + R%) * X + A(4 + R%) * Y + A(6 + R%) 490 RETURN 500 REM Display results 510 IF N < 100 OR N > 1000 THEN GOTO 560 520 IF X < XMIN THEN XMIN = X 530 IF X > XMAX THEN XMAX = X 540 IF Y < YMIN THEN YMIN = Y 550 IF Y > YMAX THEN YMAX = Y 560 IF N = 1000 THEN GOSUB 800 570 IF X > XL AND X < XH AND Y > YL AND Y < YH AND N > 1000 THEN PSET (X, Y) 590 RETURN 600 REM Test results 610 GOSUB 700 'Calculate Lyapunov exponent 620 GOSUB 900 'Calculate correlation dimension 630 IF N > 21000 THEN T% = 2 'Candidate IFS found 640 IF ABS(XNEW) + ABS(YNEW) > 1000000# THEN T% = 2 'Unbounded 650 IF N > 998 AND (F < 1 OR L > -.2) THEN T% = 2 'Uninteresting 660 IF LEN(INKEY\$) THEN T% = 0 'User key press 670 X = XNEW: Y = YNEW: N = N + 1 690 RETURN 700 REM Calculate Lyapunov exponent 710 XSAVE = XNEW: YSAVE = YNEW: X = XE: Y = YE 720 GOSUB 400 'Reiterate equations 730 DLX = XNEW - XSAVE: DLY = YNEW - YSAVE: DL2 = DLX * DLX + DLY * DLY 740 DF = 10000000000# * DL2: RS = 1# / SQR(DF) 750 XE = XSAVE + RS * (XNEW - XSAVE): YE = YSAVE + RS * (YNEW - YSAVE) 760 XNEW = XSAVE: YNEW = YSAVE 770 LSUM = LSUM + LOG(DF): L = .721347 * LSUM / N 790 RETURN 800 REM Resize the screen (and discard the first thousand iterates) 810 DX = .1 * (XMAX - XMIN): DY = .1 * (YMAX - YMIN) 820 XL = XMIN - DX: XH = XMAX + DX: YL = YMIN - DY: YH = YMAX + DY 830 IF XH - XL < .000001 OR YH - YL < .000001 THEN GOTO 890 840 WINDOW (XL, YL)-(XH, YH): CLS 850 LINE (XL, YL)-(XH, YH), , B 890 RETURN 900 REM Calculate fractal dimension 910 IF N < 200 THEN GOTO 990 'Wait for transient to settle 920 IF N = 200 THEN D2MAX = (XMAX - XMIN) ^ 2 + (YMAX - YMIN) ^ 2 930 IF N = 200 OR RND < .02 THEN XS = X: YS = Y 'New reference point 940 DX = XNEW - XS: DY = YNEW - YS 950 D2 = DX * DX + DY * DY 960 IF D2 < .001 * D2MAX THEN N2 = N2 + 1 970 IF D2 < .00001 * D2MAX THEN N1 = N1 + 1: F = .434294 * LOG(N2 / N1) 990 RETURN