1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
|
! $ mpifort -o mandelbrot.exe mandelbrot.f90
program main
use MPI
implicit none
integer master, nprocs, myrank, ierr
real :: rmsg(4), x1, x2, y1, y2
integer :: imsg(3), w, h, maxit
integer :: N
integer, allocatable :: I(:)
integer, allocatable :: C(:,:)
integer :: empty(0)
integer :: j, k
real :: x, dx, y, dy
call MPI_Init(ierr)
call MPI_Comm_get_parent(master, ierr)
if (master == MPI_COMM_NULL) then
print *, "parent communicator is MPI_COMM_NULL"
call MPI_Abort(MPI_COMM_WORLD, 1, ierr)
end if
call MPI_Comm_size(master, nprocs, ierr)
call MPI_Comm_rank(master, myrank, ierr)
! receive parameters and unpack
call MPI_Bcast(rmsg, 4, MPI_REAL, 0, master, ierr)
call MPI_Bcast(imsg, 3, MPI_INTEGER, 0, master, ierr)
x1 = rmsg(1); x2 = rmsg(2)
y1 = rmsg(3); y2 = rmsg(4)
w = imsg(1); h = imsg(2); maxit = imsg(3)
dx = (x2-x1)/real(w)
dy = (y2-y1)/real(h)
! number of lines to compute here
N = h / nprocs
if (modulo(h, nprocs) > myrank) then
N = N + 1
end if
! indices of lines to compute here
allocate( I(0:N-1) )
I = (/ (k, k=myrank, h-1, nprocs) /)
! compute local lines
allocate( C(0:w-1, 0:N-1) )
do k = 0, N-1
y = y1 + real(I(k)) * dy
do j = 0, w-1
x = x1 + real(j) * dx
C(j, k) = mandelbrot(x, y, maxit)
end do
end do
! send number of lines computed here
call MPI_Gather(N, 1, MPI_INTEGER, &
MPI_BOTTOM, 0, MPI_BYTE, &
0, master, ierr)
! send indices of lines computed here
call MPI_Gatherv(I, N, MPI_INTEGER, &
MPI_BOTTOM, empty, empty, MPI_BYTE, &
0, master, ierr)
! send data of lines computed here
call MPI_Gatherv(C, N*w, MPI_INTEGER, &
MPI_BOTTOM, empty, empty, MPI_BYTE, &
0, master, ierr)
deallocate(C)
deallocate(I)
! we are done
call MPI_Comm_disconnect(master, ierr)
call MPI_Finalize(ierr)
contains
function mandelbrot(x, y, maxit) result (it)
implicit none
real, intent(in) :: x, y
integer, intent(in) :: maxit
integer :: it
complex :: z, c
z = cmplx(0, 0)
c = cmplx(x, y)
it = 0
do while (abs(z) < 2.0 .and. it < maxit)
z = z*z + c
it = it + 1
end do
end function mandelbrot
end program main
|
