PROGRAM dotProductMPI ! ! This program computes the dot product of two vectors X,Y ! (each of size N) with component i having value i ! in parallel using P processes. ! Vectors are initialized in the code by the root process, ! then statically distributed in blocks to all processes. ! It is not assumed N is divisible by P. ! INCLUDE 'mpif.h' ! variable declarations INTEGER, PARAMETER :: N = 100 REAL, PARAMETER :: ROOT = 0 INTEGER :: P, NBAR INTEGER :: RANK, I, EXTRA, INDEX, OFFSET = 0 INTEGER :: IERR REAL :: X(N), Y(N) REAL :: DOT, DOT_LOC = 0.0 ! initialize MPI CALL MPI_INIT(IERR) IF (IERR.NE.MPI_SUCCESS) THEN PRINT*, "ERROR: MPI not initialized." STOP ENDIF ! Get the number of processors: CALL MPI_COMM_SIZE(MPI_COMM_WORLD, P, IERR) IF (IERR.NE.MPI_SUCCESS) THEN PRINT*, "ERROR: MPI processors not established." STOP ENDIF ! Get ranks of processors: CALL MPI_COMM_RANK(MPI_COMM_WORLD, RANK, IERR) IF (IERR.NE.MPI_SUCCESS) THEN PRINT*, "ERROR: MPI ranks not established." STOP ENDIF ! Root process initializes vectors x,y and distributes them IF (RANK.EQ.ROOT) THEN DO 10 I=1,N X(I) = I Y(I) = I 10 END DO ENDIF ! this could probably be done more efficiently by packing x and y ! into one entity and broadcasting it CALL MPI_BCAST(X, N, MPI_REAL, ROOT, MPI_COMM_WORLD, IERR) IF (IERR.NE.MPI_SUCCESS) THEN PRINT*, "ERROR: MPI_BCAST not successful." STOP ENDIF CALL MPI_BCAST(Y, N, MPI_REAL, ROOT, MPI_COMM_WORLD, IERR) IF (IERR.NE.MPI_SUCCESS) THEN PRINT*, "ERROR: MPI_BCAST not successful." STOP ENDIF ! determine which block of data to work on and compute dot product NBAR = N/P EXTRA = MOD(N,P) IF (RANK < EXTRA) OFFSET = 1 DO 20 I=1,NBAR+OFFSET INDEX = RANK*NBAR + I + MIN(EXTRA,RANK) DOT_LOC = DOT_LOC + X(INDEX)*Y(INDEX) 20 END DO ! gather and reduce the data and print the result CALL MPI_REDUCE(DOT_LOC, DOT, 1, MPI_REAL, MPI_SUM, ROOT, & MPI_COMM_WORLD, IERR) IF (RANK.EQ.ROOT) THEN IF (IERR.NE.MPI_SUCCESS) THEN PRINT*, "ERROR: MPI_REDUCE not successful." STOP ENDIF PRINT*, 'The dot product is: ', DOT PRINT*, 'The answer should be: ', N*(N+1)*(2*N+1)/6 ENDIF ! Finalize MPI: CALL MPI_FINALIZE(IERR) IF (IERR.NE.MPI_SUCCESS) THEN PRINT*, "ERROR: MPI not finalized." STOP ENDIF END PROGRAM dotProductMPI