Double for loop in a random walk. How do I loop over all particles simultaneously instead of looping the amount of steps for each separate particle one at a time?

I currently have a script that utilizes a double for loop in order to create a random walk in 2D. It however is looping over all the steps of the random walk one particle at a time. That is, it loops over N steps for one particle, then it moves on to the M:th particle and loops it over all N. I want to loop all particles M simultaneously instead of one at a time. I don't know how to do this, tried implementing a zero vector and inserting all values of x_t and y_t in it and plotting it but it didn't work, the figure looked weird. Here is the code:

 clear
 clc

 N = 300; % Length of the x-axis, also known as the length of the random walks.
 M = 40; % The amount of random walks.
 x_t(1) = 0;
 y_t(1) = 0;
 radie=50;
 x=0;
 y=0;
 theta=linspace(0,2*pi);

 for m=1:M
    for n = 1:N % Looping all values of N into x_t(n).
        A = sign(randn); % Generates either +1/-1 depending on the SIGN of RAND.
        x_t(n+1) = x_t(n) + A;

        A = sign(randn);
        y_t(n+1) = y_t(n) + A;

        distance = sqrt(x_t(n+1)^2 + y_t(n+1)^2);

        while distance > radie % så att de ej passerar cirkeln.
            A = sign(randn); % Generates either +1/-1 depending on the SIGN of RAND.
            x_t(n+1) = x_t(n) + A;

            A = sign(randn);
            y_t(n+1) = y_t(n) + A;
            % Compute distance again with new values.
            distance = sqrt(x_t(n+1)^2 + y_t(n+1)^2);
        end

    end
    plot(x_t, y_t,'Markersize',10)
    hold on
    plot(x+radie*cos(theta),y+radie*sin(theta));
    axis equal
 end
 grid on

I want all of the M particles to start moving N steps at the same time instead of looping all of N over each of the inputs M contains.

[Guillaume]

I don't fully understand your question, is it related to the plotting or to the calculation of x_t or y_t? What should be different from what you have already?

Unrelated to your question:

A = sign(randn); % Generates either +1/-1 depending on the SIGN of RAND.

Technically, this can also generate 0 although very unlikely. A guaranteed way to just generate +-1:

A = (-1).^randi(2);

Also,

somevalue = somelongcalculation
while somevalue > somethreshold
    somevalue = somelongcalculationlc;
end

is more elegantly written as:

somevalue = inf;
while somevalue > somethreshold
    somevalue = somelongcalculationlc;
end

[Delshad Ayoubi]

Basically, when running the script, it will run the for loop that contains N amount of "steps", that is the length of the random walk, how many steps each particle takes.

Then there's the double for loop which represents each particle, it will run M amount of times, that is M amount of particles.

What I want to do is run EVERY particle at once with N steps, instead of having the double for loop put all of the values of N in one particle, then to particle nr 2, then nr 3 to M. I want all particles M to start at one point, and then move all at once instead of one at a time.

Answer 1

[Guillaume]

If you want to change the order of the loops, then your x_t and y_t need to be 2d matrices

x_t = zeros(M, N+1);
y_t = zeros(M, N+1);
for n = 1:N
   for m = 1:M
      distance = inf;
      while distance > radie
          x_t(m, n+1) = x_t(m, n) + (-1).^randi(2);
          y_t(m, n+1) = y_t(m, n) + (-1).^randi(2);
          distance = hypot(x_t(m, n+1), y_t(m, n+1));
      end
   end
end
plot(x_t, y_t, 'Markersize', 10);

[Delshad Ayoubi]

Does it work? It runs very slowly for some reason and when I vary m the amount of particles changes whereas n changes the steps, I thought it would be the opposite.

[Guillaume]

With

N = 300;
M = 40;
radie = 50;

the above takes 0.208821 seconds on my machine. I wouldn't call that very slow.

I'm using the exact same code you were using, only simplified. n is still the same as it was in your original code.

Note that with my code, the order of the n and m loop does not matter.

Of course, if your radie is very small, my code (and your original code) will spend a lot of time retrying new movements.

[Delshad Ayoubi]

In this case, how do I know how many particles I start out with? And variable changes the length? Seem to be confused by changing the for loop, m and n.

[Guillaume]

Again, N and M are still the same variables as in your original code. So M is the number of particles.

I do not understand "And variable changes the length". Unlike your original code, all variables have their size defined before the loop starts.

However, I did make one mistake with the plot call, which should have been:

plot(x_t.', y_t.', 'MarkerSize', 10);

[Delshad Ayoubi]

For some reason the code is running extremely slowly for me. These are the changes I did

 clear
 clc
 N = 300; 
 M = 40; 
 radie = 50;
 x_t(1) = 0;
 y_t(1) = 0;
 x=0;
 y=0;
 theta=linspace(0,2*pi);
 x_t = zeros(M, N+1);
 y_t = zeros(M, N+1);

 tic
 for n = 1:N
    for m = 1:M
        distance = inf;
        while distance > radie % så att de ej passerar cirkeln.
            x_t(m, n+1) = x_t(m, n) + (-1).^randi(2);
            y_t(m, n+1) = y_t(m, n) + (-1).^randi(2);
            distance = hypot(x_t(m, n+1), y_t(m, n+1));
        end
    end
    plot(x_t.', y_t.','.', 'MarkerSize', 5);
    hold on
    plot(x+radie*cos(theta),y+radie*sin(theta));
    axis equal
 end
 toc
 grid on

I'm not sure why but it sometimes takes more than 30 sec.

[Guillaume]

It's not for some reason, it's because you're replotting all the points at each step of the loop, including the ones you've plotted on the previous steps. Instead of plotting MxN = 1200 points, you're plotting MxNx(N+1)/2 = 2257500 points, so yes it's going to be slow.

You either need to move all the plotting out of the loops, as I'd written or if you want it in the loop so you can see the progress:

   for n = 1:N
      for m = 1:M
          distance = inf;
          while distance > radie % så att de ej passerar cirkeln.
              x_t(m, n+1) = x_t(m, n) + (-1).^randi(2);
              y_t(m, n+1) = y_t(m, n) + (-1).^randi(2);
              distance = hypot(x_t(m, n+1), y_t(m, n+1));
          end
      end
      plot(x_t(:, 1:n+1)', y_t(:, 1:n+1).','.', 'MarkerSize', 5);
      hold on
      plot(x+radie*cos(theta),y+radie*sin(theta));
      axis equal
      hold off 
      drawnow
   end

The hold off and drawnow I've added are both critical. It's still going to be slower than doing the plotting after the loops because at each step you're erasing the previous plot and replotting everything.

Answer 2

[James Tursa]

If you want to walk all of the particles simultaneously, then simply reverse the order of your for-loops:

 for n = 1:N % Looping all values of N into x_t(n).
    for m=1:M

You will also need to modify all of the interior loop code to account for M particles. Currently the x_t and y_t variables are only dimensioned for 1 particle. So maybe turn them into arrays so you are using x_t(n,m) and y_t(n,m) etc within the loop.

Also this test doesn't look correct:

        while distance > radie % så att de ej passerar cirkeln.

It will always fail at the start of your run. Maybe you meant to use < instead?

[Delshad Ayoubi]

How would I change the for loop?

 for n=1:N
    for m = 1:M % Looping all values of N into x_t(n).
        A = sign(randn); % Generates either +1/-1 depending on the SIGN 
 of RAND.
        x_t(n+1) = x_t(n) + A;

        A = sign(randn);
        y_t(n+1) = y_t(n) + A;

        distance = sqrt(x_t(n+1)^2 + y_t(n+1)^2);

        while distance > radie % så att de ej passerar cirkeln.
            A = sign(randn); % Generates either +1/-1 depending on the 
 SIGN of RAND.
            x_t(n+1) = x_t(n) + A;

            A = sign(randn);
            y_t(n+1) = y_t(n) + A;
            % Compute distance again with new values.
            distance = sqrt(x_t(n+1)^2 + y_t(n+1)^2);
        end

    end
    plot(x_t, y_t,'Markersize',10)

I'm specifically talking about changing the order of the two for loops but I'm not sure how I'm supposed to change all of the x_t and y_t inside.