MPI Project Set 2
Project 4
A “token ring” is a circular messaging system. Think of a relay, with a message or “baton” being passed from one runner to the next, but around a circle so that the last “runner” passes it back to the first. Write a program that implements this with MPI. Hint: separate the ranks into root and everybody else. You may wish to use MPI_ANY_TAG and MPI_ANY_SOURCE.
Example Solutions
C++
#include <iostream>
#include <mpi.h>
using namespace std;
int main(int argc, char **argv)
{
int nProcs;
int rank;
int baton = 42;
MPI_Status status;
double startTime;
double stopTime;
MPI_Init(&argc, &argv);
MPI_Comm_size(MPI_COMM_WORLD, &nProcs);
MPI_Comm_rank(MPI_COMM_WORLD, &rank);
if (nProcs > 1)
{
if (rank == 0) {
startTime = MPI_Wtime();
MPI_Send(&baton, 1, MPI_INT, 1, 0, MPI_COMM_WORLD);
MPI_Recv(&baton, 1, MPI_INT, MPI_ANY_SOURCE,
MPI_ANY_TAG, MPI_COMM_WORLD, &status);
cout<<"Process "<<rank<<" has the baton.\n";
stopTime = MPI_Wtime();
cout<<"Elapsed time to pass a token around a ring of size "<<nProcs<<" was "<<stopTime-startTime<<"seconds\n";
fflush(stdout);
} else {
MPI_Recv(&baton, 1, MPI_INT, MPI_ANY_SOURCE,
MPI_ANY_TAG, MPI_COMM_WORLD, &status);
cout<<"Process "<<rank<<" has the baton\n";
fflush(stdout);
MPI_Send(&baton, 1, MPI_INT, (rank + 1) % nProcs,
0, MPI_COMM_WORLD);
}
}
MPI_Finalize();
return 0;
}
Fortran
program ring
use mpi
implicit none
integer :: baton
integer :: my_rank, npes
integer, parameter:: tag = 0
integer :: err, errcode
integer :: i
double precision :: start, runtime
integer, dimension(MPI_STATUS_SIZE) :: status
call MPI_INIT(err)
call MPI_COMM_RANK(MPI_COMM_WORLD, my_rank, err)
call MPI_COMM_SIZE(MPI_COMM_WORLD, npes, err)
if ( my_rank .eq. 0 ) then
baton = 42
endif
if ( npes .eq. 1 ) then
print *, 'Rank ',my_rank, 'has the baton ', baton
call MPI_Finalize(err)
stop
endif
if ( my_rank .eq. 0 ) then
start = MPI_Wtime()
call MPI_Send(baton,1,MPI_INTEGER,1,tag,MPI_COMM_WORLD,err)
call MPI_Recv(baton,1,MPI_INTEGER,MPI_ANY_SOURCE,MPI_ANY_TAG, &
MPI_COMM_WORLD,status,err)
print *, 'Rank ',my_rank, 'has the baton ', baton
call flush(6)
runtime = MPI_Wtime()-start
write(*,'(a,f0.4,a)') 'Total elapsed time was ', runtime, ' seconds'
call flush(6)
else
call MPI_Recv(baton,1,MPI_INTEGER,MPI_ANY_SOURCE,MPI_ANY_TAG, &
MPI_COMM_WORLD,status,err)
call MPI_Send(baton,1,MPI_INTEGER,mod(my_rank+1,npes),tag, &
MPI_COMM_WORLD,err)
write(*,*) 'Rank ',my_rank, 'has the baton ', baton
call flush(6)
endif
call MPI_Finalize(err)
end program ring
Python
from mpi4py import MPI
import numpy as np
import sys
comm=MPI.COMM_WORLD
my_rank=comm.Get_rank()
npes =comm.Get_size()
if my_rank==0:
baton=42*np.array([1],dtype='int')
else:
baton=np.zeros((1,),dtype='int')
if npes==1:
print('Rank ',my_rank,'has the baton ',baton)
MPI.Finalize()
sys.exit()
if my_rank==0:
start=MPI.Wtime()
comm.Send([baton,MPI.INT],1,0)
comm.Recv([baton,MPI.INT],MPI.ANY_SOURCE,MPI.ANY_TAG)
print('Rank ',my_rank,'has the baton ',baton)
runtime=MPI.Wtime()-start
print('Total elapsed time was %.4f seconds'%runtime)
else:
comm.Recv([baton,MPI.INT],MPI.ANY_SOURCE,MPI.ANY_TAG)
comm.Send([baton,MPI.INT],(my_rank+1)%npes,0)
print('Rank ',my_rank,'has the baton ',baton)
Project 5
Write a program in which each process determines a unique partner to exchange messages. One way to do this is to use
if rank < npes//2:
partner=npes//2 + rank
else
partner=rank-npes//2
where the // indicates integer division (no fractional part).
Each tasks sends its rank to its partner. Each task receives the partner’s rank. Print the message received when done. You may assume that the number of processes is even, but check that this is the case before proceeding.
Hints: do not overwrite the receiver’s rank. As always, Python programmers should take care that NumPy arrays are declared with the correct type.
Example Solutions
C++
#include <iostream>
#include "mpi.h"
using namespace std;
int main(int argc, char *argv[]) {
int rank, nprocs;
MPI_Status status;
MPI_Init(&argc, &argv);
MPI_Comm_size(MPI_COMM_WORLD, &nprocs);
MPI_Comm_rank(MPI_COMM_WORLD, &rank);
if ( nprocs%2 != 0 ) {
cout<<"This program runs only on an even number of processes\n";
MPI_Finalize();
exit(1);
}
int message;
int partner;
int half=nprocs/2;
if ( rank < half ) {
partner=half+rank;
}
else {
partner=rank-half;
}
if ( rank < half ) {
MPI_Recv(&message,1,MPI_INT,partner,0,MPI_COMM_WORLD,&status);
MPI_Send(&rank,1,MPI_INT,partner,0,MPI_COMM_WORLD);
}
else {
MPI_Send(&rank,1,MPI_INT,partner,0,MPI_COMM_WORLD);
MPI_Recv(&message,1,MPI_INT,partner,0,MPI_COMM_WORLD,&status);
}
cout<<rank<<" "<<message<<"\n";
MPI_Finalize();
exit(0);
}
Fortran
program send_receive
use mpi
implicit none
integer :: rank, partner, nprocs, ierr
integer, dimension(MPI_STATUS_SIZE) :: status
integer :: half, message
call MPI_Init(ierr)
call MPI_Comm_size(MPI_COMM_WORLD, nprocs, ierr);
call MPI_Comm_rank(MPI_COMM_WORLD, rank, ierr);
if ( mod(nprocs,2) /= 0 ) then
write(*,*) "This program runs only on an even number of processes"
call MPI_Finalize(ierr)
stop
endif
half = nprocs/2
if ( rank < half ) then
partner= half+rank
else
partner=rank-half
endif
if ( rank < nprocs/2 ) then
call MPI_Recv(message,1,MPI_INTEGER,partner,0,MPI_COMM_WORLD,status,ierr);
call MPI_Send(rank,1,MPI_INTEGER,partner,0,MPI_COMM_WORLD,ierr)
else
call MPI_Send(rank,1,MPI_INTEGER,partner,0,MPI_COMM_WORLD,ierr)
call MPI_Recv(message,1,MPI_INTEGER,partner,0,MPI_COMM_WORLD,status,ierr)
endif
write(*,*) rank, message
call MPI_Finalize(ierr)
end program
Python
import sys
import numpy as np
from mpi4py import MPI
comm=MPI.COMM_WORLD
nprocs=comm.Get_size()
rank=comm.Get_rank()
if nprocs%2 != 0:
print("This program works only for an even number of processes.")
sys.exit()
message=np.zeros(1,dtype='int')
rank_val=rank*np.ones(1,dtype='int')
half = nprocs//2
if rank < half:
partner=half+rank
else:
partner=rank-half
if rank<half:
comm.Recv([message,MPI.INT],source=partner)
comm.Send([rank_val,MPI.INT],dest=partner)
else:
comm.Send([rank_val,MPI.INT],dest=partner)
comm.Recv([message,MPI.INT],source=partner)
print(rank,message)
Project 6
A common pattern in parallel programming is the manager-worker structure. A “manager” process distributes work to “worker” processes. Sometimes manager code is separated by if rank==0 (the manager should nearly always be rank 0) statements, while the other ranks execute the “worker” code. Sometimes the manager spawns distinct worker processes, but that requires using the more advanced MPI spawn capability. In this project we will use a single code. Usually the manager distributes work to the workers, which return results; the manager then hands more work to those processes that are ready, until all is completed. For this example the workers will do only one round of work.
Starting from the stub for your language, complete the send and receive calls.
Starting Codes
C++
#include <iostream>
#include <iomanip>
#include <string>
#include <sstream>
#include <random>
#include <cmath>
#include <vector>
#include "mpi.h"
using namespace std;
int random_int(int n,int m) {
//quick and stupid way, using full C++ machinery is better but complex
return n+rand()%m;
}
double do_work() {
//hardcoded bounds for convenience
int nsteps=random_int(10000,30000);
float result=0.;
for (int i=0; i<nsteps; ++i) {
result+=i;
}
return result;
}
int main(int argc, char **argv) {
int nprocs, rank;
MPI_Status status;
MPI_Init(&argc, &argv);
MPI_Comm_size(MPI_COMM_WORLD, &nprocs);
MPI_Comm_rank(MPI_COMM_WORLD, &rank);
//Don't use rand and srand if good random numbers are needed
unsigned int seed = (unsigned) time(NULL);
//taken from some University course site
unsigned int my_seed = (seed & 0xFFFFFFF0) | (rank + 1);
srand(my_seed);
vector<float> results(nprocs);
int done=0;
float result;
int sender;
if (rank==0) {
for (int i=1; i<nprocs; ++i) {
MPI_Recv(
sender=
results[sender]=result;
done=1;
MPI_Send(&done,1,MPI_INT,sender,0,MPI_COMM_WORLD);
}
} else {
for(int n=1; n<nprocs; ++n) {
if (rank==n) {
result=do_work();
MPI_Send(
MPI_Recv(
}
}
}
float total=0;
if (rank==0) {
for (int i=1; i<nprocs; ++i ) {
total+=results[i];
}
cout << "The final result is "<<total<<"\n";
}
MPI_Finalize();
}
Fortran
module Random
implicit none
! 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 get_random_seed(seed)
! To use this subroutine, seed must be declared allocatable in the
! calling unit.
integer, dimension(:), allocatable, intent(out) :: seed
integer :: isize
call random_seed(size=isize)
if (.not. allocated(seed)) allocate(seed(isize))
call random_seed(get=seed)
end subroutine get_random_seed
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, parameter :: default_seed=2345678
call get_random_seed(iseed)
if ( .not. present(seed) ) then
call date_and_time(values=idate)
! idate(8) contains millisecond
if ( all(iseed .ne. 0) ) then
iseed = iseed * (idate(8))
else
iseed = default_seed * (idate(8))
endif
else
iseed=int(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 manager_worker
use random
use mpi
integer :: i, n
integer :: nprocs, rank, sender, ierr
integer, dimension(MPI_STATUS_SIZE) :: status;
integer :: done=0
real :: my_result, total=0.
real, dimension(:), allocatable :: results
integer, dimension(:), allocatable :: seed
interface
function do_work() result(my_result)
use random
real :: my_result
end function do_work
end interface
call MPI_Init(ierr)
call MPI_Comm_size(MPI_COMM_WORLD, nprocs, ierr)
call MPI_Comm_rank(MPI_COMM_WORLD, rank, ierr)
allocate(results(nprocs))
results=0.
! The topic of handling random-number generation in parallel programs is
! an issue in applied mathematics, and beyond our scope here. We just
! use the rank to spread out the seeds some.
call get_random_seed(seed)
seed=seed*(rank+1)
call set_random_seed(seed(1))
if (rank==0) then
do i=1,nprocs-1
call MPI_Recv(
sender=
results(sender)=my_result
done=1;
call MPI_Send(done,1,MPI_INTEGER,sender,0,MPI_COMM_WORLD,ierr)
enddo
else
do n=1,nprocs-1
if (rank==n) then
my_result=do_work();
call MPI_Send(my_result,1,MPI_REAL,0,rank,MPI_COMM_WORLD,ierr)
call MPI_Recv(done,1,MPI_INTEGER,0,MPI_ANY_TAG,MPI_COMM_WORLD,status,ierr)
endif
enddo
endif
total=sum(results)
if (rank==0) then
write(*,*) "The final result is",total
endif
call MPI_Finalize(ierr)
end program
function do_work() result(my_result)
use random
real :: my_result
integer nsteps
!hardcoded bounds for convenience
nsteps=randint(10000,30000);
my_result=0.
do i=1,nsteps
my_result=my_result+i;
enddo
end function
Python
import sys
import numpy as np
from mpi4py import MPI
def do_work():
nsteps=rng.integers(low=10000, high=30000)
result=0
for i in range(nsteps):
result+=i
return np.array([result],dtype='float')
comm=MPI.COMM_WORLD
nprocs=comm.Get_size()
rank=comm.Get_rank()
status=MPI.Status()
# The topic of handling random-number generation in parallel programs is
# an issue in applied mathematics, and beyond our scope here. This should
# generally space out the seeds a bit.
# rng is global
rng = np.random.default_rng()
done=np.array([0])
results=np.zeros(nprocs)
if rank==0:
result=np.empty(1)
for i in range(1,nprocs):
comm.Recv([result,MPI.DOUBLE],source=MPI.ANY_SOURCE,tag=MPI.ANY_TAG,status=status)
sender=status.Get_source()
results[sender]=result
done[:]=1
comm.Send([done,MPI.INT],dest=sender)
else:
for n in range(1,nprocs):
if (rank==n):
result=do_work()
comm.Send([result,MPI.DOUBLE],dest=0,tag=rank)
comm.Recv([done,MPI.INT],source=0,status=status)
total=np.sum(results)
if rank==0:
print(f"The final result is {total:e}")
#Disconnect all processes from communicator
comm.Disconnect()
Example Solutions
C++
#include <iostream>
#include <iomanip>
#include <string>
#include <sstream>
#include <random>
#include <cmath>
#include <vector>
#include "mpi.h"
using namespace std;
int random_int(int n,int m) {
//quick and stupid way, using full C++ machinery is better but complex
return n+rand()%m;
}
double do_work() {
//hardcoded bounds for convenience
int nsteps=random_int(10000,30000);
float result=0.;
for (int i=0; i<nsteps; ++i) {
result+=i;
}
return result;
}
int main(int argc, char **argv) {
int nprocs, rank;
MPI_Status status;
MPI_Init(&argc, &argv);
MPI_Comm_size(MPI_COMM_WORLD, &nprocs);
MPI_Comm_rank(MPI_COMM_WORLD, &rank);
//Don't use rand and srand if good random numbers are needed
unsigned int seed = (unsigned) time(NULL);
//taken from some University course site
unsigned int my_seed = (seed & 0xFFFFFFF0) | (rank + 1);
srand(my_seed);
vector<float> results(nprocs);
int done=0;
float result;
int sender;
if (rank==0) {
for (int i=1; i<nprocs; ++i) {
MPI_Recv(&result,1,MPI_FLOAT,MPI_ANY_SOURCE,MPI_ANY_TAG,MPI_COMM_WORLD,&status);
sender=status.MPI_SOURCE;
results[sender]=result;
done=1;
MPI_Send(&done,1,MPI_INT,sender,0,MPI_COMM_WORLD);
}
} else {
for(int n=1; n<nprocs; ++n) {
if (rank==n) {
result=do_work();
MPI_Send(&result,1,MPI_FLOAT,0,rank,MPI_COMM_WORLD);
MPI_Recv(&done,1,MPI_INT,0,MPI_ANY_TAG,MPI_COMM_WORLD,&status);
}
}
}
float total=0;
if (rank==0) {
for (int i=1; i<nprocs; ++i ) {
total+=results[i];
}
cout << "The final result is "<<total<<"\n";
}
MPI_Finalize();
}
Fortran
module Random
implicit none
! 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 get_random_seed(seed)
! To use this subroutine, seed must be declared allocatable in the
! calling unit.
integer, dimension(:), allocatable, intent(out) :: seed
integer :: isize
call random_seed(size=isize)
if (.not. allocated(seed)) allocate(seed(isize))
call random_seed(get=seed)
end subroutine get_random_seed
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, parameter :: default_seed=2345678
call get_random_seed(iseed)
if ( .not. present(seed) ) then
call date_and_time(values=idate)
! idate(8) contains millisecond
if ( all(iseed .ne. 0) ) then
iseed = iseed * (idate(8))
else
iseed = default_seed * (idate(8))
endif
else
iseed=int(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 manager_worker
use random
use mpi
integer :: i, n
integer :: nprocs, rank, sender, ierr
integer, dimension(MPI_STATUS_SIZE) :: status;
integer :: done=0
real :: my_result, total=0.
real, dimension(:), allocatable :: results
integer, dimension(:), allocatable :: seed
interface
function do_work() result(my_result)
use random
real :: my_result
end function do_work
end interface
call MPI_Init(ierr)
call MPI_Comm_size(MPI_COMM_WORLD, nprocs, ierr)
call MPI_Comm_rank(MPI_COMM_WORLD, rank, ierr)
allocate(results(nprocs))
results=0.
! The topic of handling random-number generation in parallel programs is
! an issue in applied mathematics, and beyond our scope here. We just
! use the rank to spread out the seeds some.
call get_random_seed(seed)
seed=seed*(rank+1)
call set_random_seed(seed(1))
if (rank==0) then
do i=1,nprocs-1
call MPI_Recv(my_result,1,MPI_REAL,MPI_ANY_SOURCE,MPI_ANY_TAG,MPI_COMM_WORLD,status,ierr);
sender=status(MPI_SOURCE)
results(sender)=my_result
done=1;
call MPI_Send(done,1,MPI_INTEGER,sender,0,MPI_COMM_WORLD,ierr)
enddo
else
do n=1,nprocs-1
if (rank==n) then
my_result=do_work();
call MPI_Send(my_result,1,MPI_REAL,0,rank,MPI_COMM_WORLD,ierr)
call MPI_Recv(done,1,MPI_INTEGER,0,MPI_ANY_TAG,MPI_COMM_WORLD,status,ierr)
endif
enddo
endif
total=sum(results)
if (rank==0) then
write(*,*) "The final result is",total
endif
call MPI_Finalize(ierr)
end program
function do_work() result(my_result)
use random
real :: my_result
integer nsteps
!hardcoded bounds for convenience
nsteps=randint(10000,30000);
my_result=0.
do i=1,nsteps
my_result=my_result+i;
enddo
end function
Python
import sys
import numpy as np
from mpi4py import MPI
def do_work():
nsteps=rng.integers(low=10000, high=30000)
result=0
for i in range(nsteps):
result+=i
return np.array([result],dtype='float')
comm=MPI.COMM_WORLD
nprocs=comm.Get_size()
rank=comm.Get_rank()
status=MPI.Status()
# The topic of handling random-number generation in parallel programs is
# an issue in applied mathematics, and beyond our scope here. This should
# generally space out the seeds a bit.
# rng is global
rng = np.random.default_rng()
done=np.array([0])
results=np.zeros(nprocs)
if rank==0:
result=np.empty(1)
for i in range(1,nprocs):
comm.Recv([result,MPI.DOUBLE],source=MPI.ANY_SOURCE,tag=MPI.ANY_TAG,status=status)
sender=status.Get_source()
results[sender]=result
done[:]=1
comm.Send([done,MPI.INT],dest=sender)
else:
for n in range(1,nprocs):
if (rank==n):
result=do_work()
comm.Send([result,MPI.DOUBLE],dest=0,tag=rank)
comm.Recv([done,MPI.INT],source=0,status=status)
total=np.sum(results)
if rank==0:
print(f"The final result is {total:e}")
#Disconnect all processes from communicator
comm.Disconnect()