# A Random Walk with MPI

We would like to start with a simple example that has a simple work distribution. To make it even easier, we will use no communications.

## Two-Dimensional Lattice Random Walk

A particle is moving through a two-dimensional finite grid. At each step, the “walker” can move left, right, up, or down, with equal probability. We seek the distance from the origin after N steps.

Theoretically, for N steps the distance is $\sqrt{N}$. We want to test this empirically, but one trial does not give us very good results. We want to run a lot of trials and average the results.

### Agglomeration and Mapping

This algorithm has the properties that:

• Each trial is independent of the others
• There is a fixed number of tasks
• No communications is needed between tasks to obtain independent results.

### Serial Code

Download the serial code for your language of choice. Compile it (if appropriate), then run a test case.

C++

/* Solves a random\-walk problem by a brute\-force method.
*
*/

#include <iostream>
#include <string>
#include <sstream>
#include <random>
#include <cmath>

using namespace std;

int main(int argc, char **argv) {

random_device rd;
mt19937 rng(rd());
uniform_int_distribution<int> choice(1,4);

long long N;
if (argc != 2) {
cout<<"0,0,0\n";
return 1;
}
else {
string steps=argv;
stringstream ssteps;
ssteps<<steps;
ssteps>>N;
}

double x=0;
double y=0;
int direction;

for (long long i=1; i<=N; ++i) {
direction=choice(rng);
switch (direction) {
case 1:
x+=1.;
break;
case 2:
x-=1.;
break;
case 3:
y+=1.;
break;
case 4:
y-=1.;
break;
}
}

double eucDist=sqrt(x*x+y*y);

cout<<N<<","<<sqrt((double)N)<<","<<eucDist<<"\n";
return 0;

}



Fortran

module random
implicit none

! Author:  Katherine Holcomb
! Shortened version of a longer module
! Module incorporated into this file for convenience; typically it would be in its
! own file.
! Comment out one or the other when the module is incorporated into a code.

! Single precision
!integer, parameter   :: rk = kind(1.0)
! Double precision
integer, parameter   :: rk = kind(1.0d0)

contains
subroutine set_random_seed(seed)
! Sets all elements of the seed array
integer, optional, intent(in)       :: seed

integer                             :: isize
integer, dimension(:), allocatable  :: iseed
integer, dimension(8)               :: idate
integer                             :: icount, i

call random_seed(size=isize)

if ( .not. present(seed) ) then
call system_clock(icount)
allocate(iseed(isize),source=icount+370*[(i,i=0,isize-1)])
else
allocate(iseed(isize))
iseed=seed
endif

call random_seed(put=iseed)

end subroutine set_random_seed

function urand(lb,ub,seed)
! Returns a uniformly-distributed random number in the range lb to ub.
real(rk)                        :: urand
real(rk), optional, intent(in)  :: lb,ub
real(rk), optional, intent(in)  :: seed

integer                         :: iseed
real(rk)                        :: rnd
real(rk)                        :: lower,upper

if ( present(seed) ) then
iseed=int(seed)
call set_random_seed(iseed)
endif

if ( present(lb) ) then
lower=lb
else
lower=0.0_rk
endif

if ( present(ub) ) then
upper = ub
else
upper = 1.0_rk
endif

call random_number(rnd)
urand = lower+(upper-lower)*rnd

return
end function urand

function randint(n,m)
! Returns a random integer between n and m.
integer                 :: randint
integer, intent(in)     :: n,m

randint=ceiling(urand(real(n-1,rk),real(m,rk)))

end function randint

end module random
program random_walk
use random
implicit none

integer, parameter   :: ik = selected_int_kind(15)
integer              :: nargs
character(len=10)    :: ns,nexp
integer              :: m, step
integer(ik)          :: n, nsteps
real, dimension(2)   :: pos, move
real                 :: distance, avg

nargs=command_argument_count()
if ( nargs .ne. 1 ) then
write(*,*) 0, 0., 0.
stop
else
call get_command_argument(1,ns)
endif

call set_random_seed()
pos=[0.,0.]

avg=0.

do n=1,nsteps

step=randint(1,4)
if (step==1) then
move=[0.,1.]
else if (step==2) then
move=[0.,-1.]
else if (step==3) then
move=[1.,0.]
else
move=[-1.,0.]
endif
pos=pos+move

enddo

distance=sqrt(pos(1)**2+pos(2)**2)

write(*,*) nsteps, sqrt(real(nsteps)), distance

end program



Python

import numpy as np
import random
import sys

if len(sys.argv)>1:
Nsteps=int(sys.argv)
else:
print(0, 0., 0.)
sys.exit()

pos=np.array([0,0])

moves=[np.array([0,1]),np.array([0,-1]),np.array([1,0]),np.array([-1,0])]

for n in range(Nsteps):
move=random.choice(moves)
pos+=move

distance=np.sqrt(pos**2+pos**2)
print("{:.2f},{:.2f},{:.2f}".format(Nsteps, np.sqrt(Nsteps), distance))



Using the mpi1 program in your language as an example, add the lines to run this code with MPI. In your print statement make the first output value be the rank.

## Run it

mpicxx -o randc mpirandom_walk.cxx
mpiexec -np 4 ./randc 1000000

0:1000000,1000,898.186
1:1000000,1000,189.589
2:1000000,1000,1235.87
3:1000000,1000,391.479


We have to compute the average manually

678.781


Try with -np 8

0:1000000,1000,1011.24
5:1000000,1000,1569.2
3:1000000,1000,1753.23
2:1000000,1000,967.042
1:1000000,1000,1212.17
6:1000000,1000,418.708
7:1000000,1000,881.862
4:1000000,1000,744.641


Why is the rank order jumbled?

MPI output is nondeterministic unless the programmer forces it to be ordered, using a barrier.

Previous