%% Plot the Sierpinski Gasket.

Numits=3001; % Specify the number of iterations (or number of points to graph)
NumFunc=3;   % Specify how many functions.

F=[0.5 0 0 0.5 0 0 
   0.5 0 0 0.5 50 0 
   0.5 0 0 0.5 25 50];

%% If everything above is good, you don't need to change anything below.

X=zeros(2,Numits);
% Initialize the first point randomly:
X(:,1)=randn(2,1);

for k=1:Numits-1
  idx=ceil(NumFunc*rand); % Select a function at random
% Build the matrix and vector for the function:
  A=[F(idx,1) F(idx,2); F(idx,3) F(idx,4)];  b=[F(idx,5); F(idx,6)];
  X(:,k+1)=A*X(:,k)+b;
end
plot(X(1,:),X(2,:),'k.');