% DIFFUSION LIMITED AGGREGATION along a sticky wall using Brownian motion
% 
% Read "The Science of Fractal Images", Ed. Peitgen and Saupe, p. 37 (1988)

% Author : B. Militzer, University of California, Berkeley 
% Date   : Sep 21, 2009

nParticles = 20000;
n = 0;
maxX = 500;
maxY = 500;

nNewParticlesPerFrame = 150;


A = zeros(maxX, maxY);

centerx = maxX/2;
centery = maxY/2;
count = 0;
for i = 1:size(A,1)
    for j = 1:size(A,1)
        if (i-centerx)^2 + (j-centery)^2 < 200^2;
            A(i,j) = 100;
        end
        if (i-centerx)^2 + (j-centery)^2 < 195^2;
            count = count+1;
            A(i,j) = 0;
            i_in_circle(count) = i;
            j_in_circle(count) = j;
        end
    end
end

yBuffer = 5;
yStart = 250+yBuffer;
count = 0;
for i = 1:nParticles
	
    x = centerx;
    y = centery;
    while 1
        xOld = x; 
        yOld = y;
        
        r = rand();
     
        if 0<=r && r<0.25; 
            x = x + 1;
        elseif 0.25<=r && r<0.5; 
            x = x - 1;
        elseif 0.5<=r && r<0.75; 
            y = y + 1;
        elseif 0.75<=r && r<=1; 
            y = y - 1;
        end
        
        if x>maxX; x=1; end
        if x<1; x=maxX; end
                
        if (y>yStart)   
            y = yStart; 
        end
        if y < 1; y=1; end 
        if (A(x,y)>0) 
            x=xOld;
            y=yOld;
            continue;
        end
        
        xR = x+1; 
        xL = x-1; 
        yU = y+1; 
        yD = y-1; 
        if xR>maxX; xR=1; end
        if xL<1; xL=maxX; end
        if yD < 1; yD = 1; end 
        
        if (A(xR,y) > 0 || A(xL,y) > 0 || A(x,yU) > 0 || A(x,yD) > 0)
            count = count+1;
            A(x,y) = count; 
            if (y+yBuffer>yStart && y+yBuffer<maxY) 
                yStart = y+yBuffer;
            end
            fprintf('i= %5d x=%4d y=%4d yStart=%4d\n',i,x,y,yStart);
            break;
        end

    end

    if (i==nParticles || mod(i,nNewParticlesPerFrame)==0 || yStart+1==maxY) 
        imagesc(A');
        set(gca,'YDir','normal');
        colormap('hot')
        axis equal;
        axis tight;
        drawnow
        n = n+1;
        filename = ['2ndinwanimation_frame_number_',sprintf('%03d',n),'.png']
        print(filename,'-r600','-dpng')
    end

    if (yStart+1==maxY) 
        fprintf('Growth reached Y limit after only %d particles\n',i);
        break;
    end
    if A(centerx,centery) > 0
        fprintf('Growth reached center limit after only %d particles\n',i);
        break;
    end  
end
p=1;
title(['Partical only sticks with p = ' num2str(p)])
