9 #include "instr/instr.h" //TRACE_
21 //---NOTE : all the timers function have been modified to
22 // avoid global timers (privatize these).
23 // ----------------------- timers ---------------------
24 void timer_clear(double *onetimer) {
29 void timer_start(double *onetimer) {
30 *onetimer = MPI_Wtime();
33 void timer_stop(int n,double *elapsed,double *start) {
41 double timer_read(int n, double *elapsed) { /* ok, useless, but jsut to keep function call */
44 /********************************************************************
45 ***************** V R A N L C ******************
46 ***************** *****************/
47 double vranlc(int n, double x, double a, double *y)
50 long i246m1=0x00003FFFFFFFFFFF;
54 // This doesn't work, because the compiler does the calculation in 32
55 // bits and overflows. No standard way (without f90 stuff) to specify
56 // that the rhs should be done in 64 bit arithmetic.
57 // parameter(i246m1=2**46-1)
61 // c Note that the v6 compiler on an R8000 does something stupid with
62 // c the above. Using the following instead (or various other things)
63 // c makes the calculation run almost 10 times as fast.
66 // c data d2m46/0.0d0/
67 // c if (d2m46 .eq. 0.0d0) then
68 // c d2m46 = 0.5d0**46
73 //fprintf(stdout,("================== Vranlc ================");
74 //fprintf(stdout,("Before Loop: Lx = " + Lx + ", La = " + La);
76 for (i=0; i< n; i++) {
79 y[i] = d2m46 * (double)LLx;
82 fprintf(stdout,("After loop 0:");
83 fprintf(stdout,("Lx = " + Lx + ", La = " + La);
84 fprintf(stdout,("d2m46 = " + d2m46);
85 fprintf(stdout,("LLX(Lx) = " + LLX.doubleValue());
86 fprintf(stdout,("Y[0]" + y[0]);
93 fprintf(stdout,("Change: Lx = " + Lx);
94 fprintf(stdout,("=============End Vranlc ================");
101 //-------------- the core (unique function) -----------
102 void doTest(int argc, char **argv) {
103 double dum[3] = {1.,1.,1.};
104 double x1, x2, sx, sy, tm, an, tt, gc;
106 double epsilon=1.0E-8, a = 1220703125., s=271828183.;
107 double t1, t2, t3, t4;
108 double sx_verify_value, sy_verify_value, sx_err, sy_err;
110 #include "npbparams.h"
112 // --> set by make : in npbparams.h
113 //m=28, // for CLASS=A
114 //m=30, // for CLASS=B
117 nn = (int)(pow(2,mm)),
118 nk = (int)(pow(2,mk)),
127 k, nit, no_large_nodes,
129 int me, nprocs, root=0, dp_type;
132 char size[500]; // mind the size of the string to represent a big number
136 double R23, R46, T23, T46;
138 double *qq = (double *) malloc (10000*sizeof(double));
139 double *start = (double *) malloc (64*sizeof(double));
140 double *elapsed = (double *) malloc (64*sizeof(double));
142 double *x = (double *) malloc (2*nk*sizeof(double));
143 double *q = (double *) malloc (nq*sizeof(double));
145 TRACE_smpi_set_category ("start");
147 MPI_Init( &argc, &argv );
148 MPI_Comm_size( MPI_COMM_WORLD, &no_nodes);
149 MPI_Comm_rank( MPI_COMM_WORLD, &node);
157 /* Because the size of the problem is too large to store in a 32-bit
158 * integer for some classes, we put it into a string (for printing).
159 * Have to strip off the decimal point put in there by the floating
160 * point print statement (internal file)
162 fprintf(stdout," NAS Parallel Benchmarks 3.2 -- EP Benchmark");
163 sprintf(size,"%d",pow(2,m+1));
164 //size = size.replace('.', ' ');
165 fprintf(stdout," Number of random numbers generated: %s\n",size);
166 fprintf(stdout," Number of active processes: %d\n",no_nodes);
171 /* c Compute the number of "batches" of random number pairs generated
172 c per processor. Adjust if the number of processors does not evenly
173 c divide the total number
177 no_large_nodes = nn % no_nodes;
178 if (node < no_large_nodes) np_add = 1;
183 fprintf(stdout,"Too many nodes: %d %d",no_nodes,nn);
184 MPI_Abort(MPI_COMM_WORLD,1);
188 /* c Call the random number generator functions and initialize
189 c the x-array to reduce the effects of paging on the timings.
190 c Also, call all mathematical functions that are used. Make
191 c sure these initializations cannot be eliminated as dead code.
194 //call vranlc(0, dum[1], dum[2], dum[3]);
195 // Array indexes start at 1 in Fortran, 0 in Java
196 vranlc(0, dum[0], dum[1], &(dum[2]));
198 dum[0] = randlc(&(dum[1]),&(dum[2]));
199 /////////////////////////////////
200 for (i=0;i<2*nk;i++) {
203 Mops = log(sqrt(abs(1)));
206 c---------------------------------------------------------------------
207 c Synchronize before placing time stamp
208 c---------------------------------------------------------------------
210 MPI_Barrier( MPI_COMM_WORLD );
213 TRACE_smpi_set_category ("ep");
215 timer_clear(&(elapsed[1]));
216 timer_clear(&(elapsed[2]));
217 timer_clear(&(elapsed[3]));
218 timer_start(&(start[1]));
221 //fprintf(stdout,("(ep.f:160) t1 = " + t1);
222 t1 = vranlc(0, t1, a, x);
223 //fprintf(stdout,("(ep.f:161) t1 = " + t1);
226 /* c Compute AN = A ^ (2 * NK) (mod 2^46). */
229 //fprintf(stdout,("(ep.f:165) t1 = " + t1);
230 for (i=1; i <= mk+1; i++) {
231 t2 = randlc(&t1, &t1);
232 //fprintf(stdout,("(ep.f:168)[loop i=" + i +"] t1 = " + t1);
235 //fprintf(stdout,("(ep.f:172) s = " + s);
240 for (i=0; i < nq ; i++) {
245 Each instance of this loop may be performed independently. We compute
246 the k offsets separately to take into account the fact that some nodes
247 have more numbers to generate than others
251 k_offset = node * np -1;
253 k_offset = no_large_nodes*(np+1) + (node-no_large_nodes)*np -1;
256 for(k = 1; k <= np; k++) {
260 //fprintf(stdout,("(ep.f:193) t1 = " + t1);
263 // Find starting seed t1 for this kk.
265 for (i=1;i<=100 && !stop;i++) {
267 //fprintf(stdout,("(ep.f:199) ik = " +ik+", kk = " + kk);
269 t3 = randlc(&t1, &t2);
270 //fprintf(stdout,("(ep.f:200) t1= " +t1 );
275 t3 = randlc(&t2, &t2);
279 // Compute uniform pseudorandom numbers.
281 //if (timers_enabled) timer_start(3);
282 timer_start(&(start[3]));
283 //call vranlc(2 * nk, t1, a, x) --> t1 and y are modified
285 //fprintf(stdout,">>>>>>>>>>>Before vranlc(l.210)<<<<<<<<<<<<<");
286 //fprintf(stdout,"2*nk = " + (2*nk));
287 //fprintf(stdout,"t1 = " + t1);
288 //fprintf(stdout,"a = " + a);
289 //fprintf(stdout,"x[0] = " + x[0]);
290 //fprintf(stdout,">>>>>>>>>>>>>>>>>>>>>><<<<<<<<<<<<<<<<<<<<<<");
292 t1 = vranlc(2 * nk, t1, a, x);
294 //fprintf(stdout,(">>>>>>>>>>>After Enter vranlc (l.210)<<<<<<");
295 //fprintf(stdout,("2*nk = " + (2*nk));
296 //fprintf(stdout,("t1 = " + t1);
297 //fprintf(stdout,("a = " + a);
298 //fprintf(stdout,("x[0] = " + x[0]);
299 //fprintf(stdout,(">>>>>>>>>>>>>>>>>>>>>><<<<<<<<<<<<<<<<<<<<<<");
301 //if (timers_enabled) timer_stop(3);
302 timer_stop(3,elapsed,start);
304 /* Compute Gaussian deviates by acceptance-rejection method and
305 * tally counts in concentric square annuli. This loop is not
308 //if (timers_enabled) timer_start(2);
309 timer_start(&(start[2]));
310 for(i=1; i<=nk;i++) {
311 x1 = 2. * x[2*i-2] -1.0;
312 x2 = 2. * x[2*i-1] - 1.0;
315 t2 = sqrt(-2. * log(t1) / t1);
318 l = (int)(abs(t3) > abs(t4) ? abs(t3) : abs(t4));
325 fprintf(stdout,"x1 = " + x1);
326 fprintf(stdout,"x2 = " + x2);
327 fprintf(stdout,"t1 = " + t1);
328 fprintf(stdout,"t2 = " + t2);
329 fprintf(stdout,"t3 = " + t3);
330 fprintf(stdout,"t4 = " + t4);
331 fprintf(stdout,"l = " + l);
332 fprintf(stdout,"q[l] = " + q[l]);
333 fprintf(stdout,"sx = " + sx);
334 fprintf(stdout,"sy = " + sy);
338 //if (timers_enabled) timer_stop(2);
339 timer_stop(2,elapsed,start);
342 TRACE_smpi_set_category ("finalize");
344 //int MPI_Allreduce(void *sbuf, void *rbuf, int count, MPI_Datatype dtype, MPI_Op op, MPI_Comm comm)
345 MPI_Allreduce(&sx, x, 1, MPI_DOUBLE, MPI_SUM, MPI_COMM_WORLD);
346 sx = x[0]; //FIXME : x[0] or x[1] => x[0] because fortran starts with 1
347 MPI_Allreduce(&sy, x, 1, MPI_DOUBLE, MPI_SUM, MPI_COMM_WORLD);
349 MPI_Allreduce(q, x, nq, MPI_DOUBLE, MPI_SUM, MPI_COMM_WORLD);
351 for(i = 0; i < nq; i++) {
354 for(i = 0; i < nq; i++) {
358 timer_stop(1,elapsed,start);
359 tm = timer_read(1,elapsed);
360 MPI_Allreduce(&tm, x, 1, MPI_DOUBLE, MPI_MAX, MPI_COMM_WORLD);
368 sx_verify_value = -3.247834652034740E3;
369 sy_verify_value = -6.958407078382297E3;
371 sx_verify_value = -2.863319731645753E3;
372 sy_verify_value = -6.320053679109499E3;
374 sx_verify_value = -4.295875165629892E3;
375 sy_verify_value = -1.580732573678431E4;
377 sx_verify_value = 4.033815542441498E4;
378 sy_verify_value = -2.660669192809235E4;
380 sx_verify_value = 4.764367927995374E4;
381 sy_verify_value = -8.084072988043731E4;
383 sx_verify_value = 1.982481200946593E5;
384 sy_verify_value = -1.020596636361769E5;
390 fprintf(stdout,("sx = " + sx);
391 fprintf(stdout,("sx_verify = " + sx_verify_value);
392 fprintf(stdout,("sy = " + sy);
393 fprintf(stdout,("sy_verify = " + sy_verify_value);
396 sx_err = abs((sx - sx_verify_value)/sx_verify_value);
397 sy_err = abs((sy - sy_verify_value)/sy_verify_value);
399 fprintf(stdout,("sx_err = " + sx_err);
400 fprintf(stdout,("sy_err = " + sx_err);
401 fprintf(stdout,("epsilon= " + epsilon);
403 verified = ((sx_err < epsilon) && (sy_err < epsilon));
406 Mops = (pow(2.0, m+1))/tm/1000;
408 fprintf(stdout,"EP Benchmark Results:\n");
409 fprintf(stdout,"CPU Time=%d\n",tm);
410 fprintf(stdout,"N = 2^%d\n",m);
411 fprintf(stdout,"No. Gaussain Pairs =%d\n",gc);
412 fprintf(stdout,"Sum = %lf %ld\n",sx,sy);
413 fprintf(stdout,"Count:");
414 for(i = 0; i < nq; i++) {
415 fprintf(stdout,"%d\t %ld\n",i,q[i]);
419 print_results("EP", _class, m+1, 0, 0, nit, npm, no_nodes, tm, Mops,
420 "Random numbers generated", verified, npbversion,
421 compiletime, cs1, cs2, cs3, cs4, cs5, cs6, cs7) */
422 fprintf(stdout,"\nEP Benchmark Completed\n");
423 fprintf(stdout,"Class = %s\n", _class);
424 fprintf(stdout,"Size = %s\n", size);
425 fprintf(stdout,"Iteration = %d\n", nit);
426 fprintf(stdout,"Time in seconds = %lf\n",(tm/1000));
427 fprintf(stdout,"Total processes = %d\n",no_nodes);
428 fprintf(stdout,"Mops/s total = %lf\n",Mops);
429 fprintf(stdout,"Mops/s/process = %lf\n", Mops/no_nodes);
430 fprintf(stdout,"Operation type = Random number generated\n");
432 fprintf(stdout,"Verification = SUCCESSFUL\n");
434 fprintf(stdout,"Verification = UNSUCCESSFUL\n");
436 fprintf(stdout,"Total time: %lf\n",(timer_read(1,elapsed)/1000));
437 fprintf(stdout,"Gaussian pairs: %lf\n",(timer_read(2,elapsed)/1000));
438 fprintf(stdout,"Random numbers: %lf\n",(timer_read(3,elapsed)/1000));
441 MPE_Finish_log(argv[0]);
447 int main(int argc, char **argv) {