Fractals

The Sierpinski carpet is a straightforward example of a fractal that you can code conveniently using recursion. We converted the c code in Wikipedia to Pascal. Program SierpinskiCarpet recursively divides each blank square into a grid of nine squares and colours the central square. It produced the following output.

Output from program SierpinskiCarpet

Output from program SierpinskiCarpet

Program SierpinskiCarpet uses the graph unit which is available in Lazarus but not in Delphi. See also the Smart Pascal version.

program SierpinskiCarpet;
  {$mode objfpc}{$H+}
uses
 SysUtils, Graph;
var
  Gd, Gm : smallint;
  CountX, CountY, MaxX, MaxY : integer;

function PixelToBeFilled (x, y, Width, Height : integer) : Boolean;
var
  x2, y2 : integer;
begin
  // First stopping condition
  if x < 1  then
    result := False
  else
    begin
      //If the grid was split in 9 parts, what part(x2,y2) would x,y fit into?
      x2 := (x * 3) DIV Width; //an integer from 0..2 inclusive
      y2 := (y * 3) DIV Height; //an integer from 0..2 inclusive
      //Second stopping condition
      //if in the centre square, it should be filled.
      if (x2 = 1) and (y2 = 1) then 
        result := True
      else
        begin
          // General case
          {Offset x and y so it becomes bounded by 0 .. width DIV 3 and
          0.. height DIV 3 and prepares for recursive call.}
          x := x - (x2 * Width) DIV 3;
          y := y - (y2 * Height) DIV 3;
          //Recursive call
          Result := PixelToBeFilled(x, y, Width DIV 3, Height DIV 3);
        end;
    end;
end;

begin
  write('When the graphics window opens, you may need to restore'+
        ' down or minimise'#13#10'the window then maximise it' +
        ' to see all of the graphic.'#13#10 +
          'Please press Return to continue.');
  readln;
  Gd := Detect;
  Gm := 0;
  InitGraph(Gd, Gm, '');
  //For best results length of side should be a power of 3.
  MaxY :=  3;
  repeat
    MaxY := MaxY * 3;
  until MaxY * 3 > GetMaxY;
  MaxX := MaxY;
  //Draw carpet
  for CountX := 1 to MaxX do
    for CountY := 1 to MaxY do
      if PixelToBeFilled (CountX, CountY, MaxX, MaxY) then
        putPixel(CountX, CountY, Red);
  readln;
end.         

Sierpinski Curves

Program SierpinskiCurve provides good examples of mutual recursion. The code is based on the outline given in Recursion via Pascal by J. S. Rohl. The code was designed to enable continuous drawing by a plotter rather than for the console, which had poor resolution back in 1984! We show the curves of order 1 and 3. Try to predict the output of other orders before running them. We now have an online Smart Pascal demonstration.

Sierpinski Curve (order 1)

Sierpinski Curve (order 1)

Sierpinski Curve (order 3)

Sierpinski Curve (order 3)

Lines and procedures are named by the letters of the compass directions that the output is following. For example, LineSE draws a diagonal line from top left to bottom right and procedure N draws a sequence of lines starting one diagonal line above bottom of the graphic and finishing one diagonal line below the top. The following graphic shows the increasing complexity of the output of procedure N as the order increases from 1 to 3.

Output from procedure N with orders 1, 2 and 3

Output from procedure N with orders 1, 2 and 3

We include forward declarations of procedures S, E and W because they call each other and otherwise the calls to procedures not yet coded would not be accepted by the compiler.

program SierpinskiCurve;
  {$mode objfpc}{$H+}
uses
 SysUtils, Graph, Math;
type
  order = 1..7;
var
  SleepTime, LeftMargin, BottomMargin, GraphicWidth, x, y, h: integer;
  Gd, Gm : smallint;
  UserOrder : order;   

procedure S(i: order); forward;
procedure E(i: order); forward;
procedure W(i: order); forward;

procedure LineN;
begin
  y := y - 2 * h;
  lineTo(x, y);
end;

procedure LineS;
begin
  y := y + 2 * h;
  lineTo(x, y);
end;

procedure LineE;
begin
  x := x + 2 * h;
  lineTo(x, y);
end;

procedure LineW;
begin
  x := x - 2 * h;
  lineTo(x, y);
end;

procedure LineNW;
begin
  x := x - h;
  y := y - h;
  lineTo(x, y);
end;
procedure LineNE;
begin
  x := x + h;
  y := y - h;
  lineTo(x, y);
end;

procedure LineSE;
begin
  x := x + h;
  y := y + h;
  lineTo(x, y);
end;

procedure LineSW;
begin
  x := x - h;
  y := y + h;
  lineTo(x, y);
end;

 procedure N(i: order);
 begin
   if i = 1 then
     begin
       LineNE;
       LineN;
       LineNW;
       sleep(SleepTime);
     end
   else
     begin
       N(i - 1);
       LineNE;
       E(i - 1);
       LineN;
       W(i - 1);
       LineNW;
       N(i - 1);
     end;
 end;

procedure E(i: order);
begin
  if i = 1 then
    begin
      LineSE;
      LineE;
      LineNE;
      sleep(SleepTime);
    end
  else
    begin
      E(i-1);
      LineSE;
      S(i-1);
      LineE;
      N(i-1);
      LineNE;
      E(i-1);
    end;
end;

procedure S(i: order);
begin
  if i = 1 then
    begin
      LineSW;
      LineS;
      LineSE;
      sleep(SleepTime);
    end
  else
    begin
      S(i - 1);
      LineSW;
      W(i - 1);
      LineS;
      E(i - 1);
      LineSE;
      S(i - 1);
    end;
end;

procedure W(i: order);
begin
  if i = 1 then
    begin
      LineNW;
      LineW;
      LineSW;
      sleep(SleepTime);
    end
  else
    begin
      W(i - 1);
      LineNW;
      N(i - 1);
      LineW;
      S(i - 1);
      LineSW;
      W(i - 1);
    end;
end;

procedure Sierpinski(i: order);
begin
  x := LeftMargin;
  y := GetMaxY - (h + BottomMargin); 
  moveTo(x, y );
  N(i);
  LineNE;
  E(i);
  LineSE;
  S(i);
  LineSW;
  W(i);
  LineNW;
end;

begin
  writeln('When the graphics window opens, you may need to restore'+
        ' down or minimise'#13#10'the window then maximise it' +
        ' to see all of the graphic.'#13#10);
  write( 'Enter complexity integer (1 for basic shape, 7 for most complex). ');
  readln(UserOrder);
  if UserOrder < 4 then
    SleepTime := 600 DIV UserOrder
  else
    SleepTime := 0;
  Gd := Detect;
  Gm := 0;
  InitGraph(Gd, Gm, '');
  //GraphicWidth is number of units
  GraphicWidth := trunc(power(2, UserOrder + 2) - 2);
  //h is size of 1 unit
  h:= GetMaxY DIV GraphicWidth;
  LeftMargin := (GetMaxX - h * GraphicWidth) DIV 2;
  BottomMargin := (GetMaxY - h * GraphicWidth) DIV 2;
  Sierpinski(UserOrder);
  writeln('Please press Return to exit.');  
  readln;
end.
 

To run the program in Delphi 7, download and add wingraph.pas to the project as described in a section of our graphics tutorial. The first four lines should be as follows, with the rest of the program as above for Lazarus.

program SierpinskiCurve;
  {$APPTYPE CONSOLE}
uses
  SysUtils, Math, wingraph in 'wingraph.pas';  

Experimenting with Fractals

You should be tempted to ring the changes on these programs, for example by changing line colours, fill colours, and line drawing procedures. Try to superimpose fractals with different orders. You can find plenty of other fractals to code and join the growing number of fractal artists.

Programming - a skill for life!

First textbook examples of recursion, binary search trees and fractals