0

I have created a function that generates N particles to do a random walk, 1 at a time.

I have tracked the end location of each particle ('result'). However, I would like to keep the particle in the end position (which is when y=99). When the new particle is generated the end position should be when y=99 OR when a particle lands on top of the end location of the previous particle.

Can anyone please point me in the direction to achieve this?

At the moment my code stops the particle when y=99 but does not consider if it is landing on top of another particle (ie.if there is already a particle at (50,99) the new particle should stop at (50,98) not (50,99))

function [result]=random_walk_gravity(N,s,w,e)

%N (the number of particles),  and s, w, e (the probabilities of moving south, west and east).

result = zeros(1,N);

for k=1:N

    %start particle i at x0, y0
    %In this case, all particles  will start in column 50 in the top row
    %step with stepsize 1 either South, West or East

    %sample u(0,1) with the equal probability of moving in the east, west or
    %south direction

    i=1;


    % x(i)=50;
    % y(i)=1;


    y(i)=1;
    x(i)=50;

    %for k=1:n
    %generate random number between 0 and 1, N times
    u=rand();


    while y(i)<99 %we want it to stop at y=99 and look at the x-value

        if (0<=u)&&(u<=e) && (x(i)<99) %Check that the x position is less than 99
            %Move east, therefore xposition is increased by 1 unit
            x(i+1)= x(i)+1;
            y(i+1)= y(i);
            i=i+1;
            u=rand();

        elseif (0<=u)&&(u<=e) && (x(i)>=99) %Check if the x position is greater than or equal to 99
            %particle hits the side boundary it moves back in the direction
            %it came from (west)
            x(i+1)= x(i)-1;
            y(i+1)= y(i);
            i=i+1;
            u=rand();

        elseif (e<u)&&(u<=(e+w))&& (x(i)>0) %Check that the x position is greater than 0
            %Move west, therefore xposition is increased by 1 unit
            x(i+1)= x(i)-1;
            y(i+1)= y(i);
            i=i+1;
            u=rand();

        elseif (e<u)&&(u<=(e+w))&& (x(i)<0)
            %Check if the x position is less than or equal to 0
            %particle hits the side boundary it moves back in the direction
            %it came from (east)
            x(i+1)= x(i)+1;
            y(i+1)= y(i);
            i=i+1;
            u=rand();

        else %Move south, therefore yposition is increased by 1 unit
            x(i+1)= x(i);
            y(i+1)= y(i)+1;
            i=i+1;
            u=rand();

        end

    end

    result(k) = x(i);
    histogram(result)
end
Lauren
  • 9
  • 2

1 Answers1

0

I'm slightly unclear on a couple of things but the code below should stop a particle if it moves into a position directly 1 step above the resting position of any previous particle.

  • Note I fixed what seemed like a typo on the line elseif (e<u)&&(u<=(e+w))&& (x(i)<0)

  • I split your results array into a resultsX and resultsY. This isn't strictly necessary if it risks causing problems elsewhere, but seemed to make the calculation easier.

  • The main addition is that at the end of each while iteration (after the particle has moved), the code checks whether there are any previous particles with final resting coordinates that are directly 1 step below this particle's current position (i.e. same x coordinate, and 1 more than the current y coordinate). If there is, then the while loop breaks and the particle stops there.

function [resultX, resultY]=random_walk_gravity(N,s,w,e)

%N (the number of particles),  and s, w, e (the probabilities of moving south, west and east).

resultX = nan(1,N); % I'd recommend using NaN rather than 0 here.
resultY = nan(1,N);

for k=1:N

    %start particle i at x0, y0
    %In this case, all particles  will start in column 50 in the top row
    %step with stepsize 1 either South, West or East

    %sample u(0,1) with the equal probability of moving in the east, west or
    %south direction

    i=1;


    % x(i)=50;
    % y(i)=1;


    y(i)=1;
    x(i)=50;

    %for k=1:n
    %generate random number between 0 and 1, N times
    u=rand();


    while y(i)<99 %we want it to stop at y=99 and look at the x-value

        if (0<=u)&&(u<=e) && (x(i)<99) %Check that the x position is less than 99
            %Move east, therefore xposition is increased by 1 unit
            x(i+1)= x(i)+1;
            y(i+1)= y(i);
            i=i+1;
            u=rand();

        elseif (0<=u)&&(u<=e) && (x(i)>=99) %Check if the x position is greater than or equal to 99
            %particle hits the side boundary it moves back in the direction
            %it came from (west)
            x(i+1)= x(i)-1;
            y(i+1)= y(i);
            i=i+1;
            u=rand();

        elseif (e<u)&&(u<=(e+w))&& (x(i)>0) %Check that the x position is greater than 0
            %Move west, therefore xposition is increased by 1 unit
            x(i+1)= x(i)-1;
            y(i+1)= y(i);
            i=i+1;
            u=rand();

        elseif (e<u)&&(u<=(e+w))&& (x(i)<=0) % You had just <0 rather than <=0
            %Check if the x position is less than or equal to 0
            %particle hits the side boundary it moves back in the direction
            %it came from (east)
            x(i+1)= x(i)+1;
            y(i+1)= y(i);
            i=i+1;
            u=rand();

        else %Move south, therefore yposition is increased by 1 unit
            x(i+1)= x(i);
            y(i+1)= y(i)+1;
            i=i+1;
            u=rand();

        end

        % If there are any particles with the same X coordinate that
        % stopped 1 Y coordinate below this particle's current position,
        % then end the walk now.
        if any(x(i)==resultX & (y(i)+1)==resultY)
            break
        end

    end

    resultX(k) = x(i);
    resultY(k) = y(i);
    histogram(resultX)
end

Obviously say if this isn't quite what you meant. However, in particular there are a couple things I was unsure of:

  • Do you want the particle to stop immediately if it moves above the resting point of a previous particle, or only if it tries to move south into the resting position of a previous particle. For example, if there is a resting particle at (X=51,Y=99) and a new particle moves west from (X=52,Y=98) to (X=51,Y=98) so that it is directly above the previous particle, do you want to the particle to stop immediately (as in the above code) or only if it then tries to move south on the next step? The difference here is whether you want a particle to be able to move horizontally just above a previous particle without getting 'stuck' or not.

  • The code above will let multiple particles stack on top of each other. So if three particles have already come to rest at X=51, their respective Y positions would be 99, 98 and 97. Another particle would then 'stop' if it moves into coordinate (X=51,Y=96). However, what do you want to happen if a particle tried to move west from coordinate (X=52,Y=97) into (X=51,Y=97). Here, a new particle is trying to move into the resting position of a previous particle, but from the side (rather than above). What would you want to happen in a sitation like this? For example, you might want particles to be able to (a) move through previous particles in this way, (b) to reverse direction (as when it hits the boundaries at 0 and 99), (c) to fall straight down and stop, (d) to stop where it is (i.e. gets 'stuck' to the side of a stack of previous particles, even if there aren't any particles directly below it)?

Craig Arnold
  • 58
  • 5