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
98
99
100
101
102
103
104
105
106
|
function MPIPi(N,mod)
// Pi_seq: Compute PI (sequential) by num.integr. of arctan(x) in [0..1]
//
// Pi_seq [ (N) ]
//
// N [1E7] #subdivisions of the [0, 1] interval
// mod ['s'] communication modality: (s)end (r)educe
// printed results struct contains
// pi estimated pi value
// err error
// time from argument xmit to pi computed
//
// Examples:
//
// octave:1> Pi_seq
// results =
// {
// pi = 3.1416
// err = -6.4837e-14
// time = 2.7761
// }
//
////////////////////
// ArgChk //
////////////////////
if argn(2)<1, N=1E7; end
if argn(2)<2, mod="s"; end
if argn(2)>2, error("usage MPIPi(N,mod)"); end
flag=0;
[%H,%ierr] = evstr(string(N));
flag=flag | size(N,"*")<>1 | %ierr<>0;
flag=flag | fix(N)~=N | N<1;
if flag, error("usage MPIPi( <int> N>0, <char> mod==''s|r'' )"), end
////////////////////////////////////
// Results struct //
////////////////////////////////////
results.pi =0;
results.err =0;
results.time=0;
////////////
// PARALLEL / initialization, include MPI_Init time in measurement
////////////
// T=clock; /
////////////
MPI_Init(); // should have lambooted first
rnk = MPI_Comm_rank(); // let it abort if it fails
siz = MPI_Comm_size();
SLV = rnk; // handy shortcuts, master is rank 0
MST = ~ SLV; // slaves are all other
if MST
disp("Master here")
else
disp("Slave here")
end
////////////
// PARALLEL / computation (depends on rank/size)
//////////// // vectorized code, equivalent to
width=1/N; lsum=0; // for i=rnk:siz:N-1
i=rnk:siz:N-1; // x=(i+0.5)*width;
x=(i+0.5)*width; // lsum=lsum+4/(1+x^2);
lsum=sum(4 ./(1+x.^2)); // end
// mem-bound, N=1E7 => size(i)=8E7/siz (80MB!!!)
////////////
// PARALLEL / reduction and finish
////////////
select mod
case "s", TAG=7; // Any tag would do
if SLV // All slaves send result back
MPI_Send(lsum, 0,TAG);
else // Here at master
Sum =lsum; // save local result
for slv=1:siz-1 // collect in any order
// MPI_Recv(lsum,MPI_ANY_SOURCE,TAG,MPI_COMM_WORLD);
MPI_Recv(lsum,TAG);
Sum = Sum + lsum; // and accumulate
end // order: slv or MPI_ANY_SOURCE
end
case "r", Sum=0; // reduction master = rank 0 @ WORLD
MPI_Reduce(lsum,Sum, MPI_SUM, 0,MPI_COMM_WORLD);
end
MPI_Finalize;
//////////////////////////////////////////////////////////////////
//results.time = etime(clock,T); //
//////////////////////////////////////////////////////////////////
results.err = Sum-%pi;
results.pi = Sum // ;
//////////////////////////////
// FINAL CHECK //
//////////////////////////////
if abs(results.err)>5E-10
warning("pretty nice pi value! go fix it")
end
disp(results)
endfunction
|