xref: /phasta/phSolver/common/new_interface.c (revision 54d5a0aae21e69d70baf84f3685cd7c844f5d326)
1 /* This file provides interface functions for 'partial ' random
2    access into the PHASTA input files
3 
4    Anil Karanam March 2001 */
5 
6 #include <stdio.h>
7 #include <string.h>
8 #include <ctype.h>
9 #include <stdlib.h>
10 #include <time.h>
11 #include <math.h>
12 #include "mpi.h"
13 #include "phastaIO.h"
14 #include "rdtsc.h"
15 #include <FCMangle.h>
16 #include "new_interface.h"
17 #include "phIO.h"
18 
19 //MR CHANGE
20 #include "common_c.h"
21 //MR CHANGE END
22 
23 #ifdef intel
24 #include <winsock2.h>
25 #else
26 #include <unistd.h>
27 #include <strings.h>
28 #endif
29 
30 //extern double cpu_speed = 2600000000.0; //for Jaguar XT5
31 //extern double cpu_speed = 2100000000.0;  //for Jaguar xt4
32 //extern double cpu_speed =  850000000.0;  //for Intrepid
33 
34 void igetMinMaxAvg(int *ivalue, double *stats, int *statRanks) {
35   int isThisRank;
36 
37   double *value = (double*)malloc(sizeof(double));
38   *value = 1.0*(*ivalue);
39 
40   rgetMinMaxAvg(value,stats,statRanks);
41 
42   /* MPI_Allreduce(value,&stats[0],1,MPI_DOUBLE,MPI_MIN,MPI_COMM_WORLD);
43   isThisRank=workfc.numpe+1;
44   if(*value==stats[0])
45     isThisRank=workfc.myrank;
46   MPI_Allreduce(&isThisRank,&statRanks[0],1,MPI_INT,MPI_MIN,MPI_COMM_WORLD);
47 
48   MPI_Allreduce(value,&stats[1],1,MPI_DOUBLE,MPI_MAX,MPI_COMM_WORLD);
49   isThisRank=workfc.numpe+1;
50   if(*value==stats[1])
51     isThisRank=workfc.myrank;
52   MPI_Allreduce(&isThisRank,&statRanks[1],1,MPI_INT,MPI_MIN,MPI_COMM_WORLD);
53 
54   MPI_Allreduce(value,&stats[2],1,MPI_DOUBLE,MPI_SUM,MPI_COMM_WORLD);
55   stats[2] /= workfc.numpe; */
56 
57   free(value);
58 }
59 
60 void rgetMinMaxAvg(double *value, double *stats, int *statRanks) {
61   int isThisRank;
62 
63   MPI_Allreduce(value,&stats[0],1,MPI_DOUBLE,MPI_MIN,MPI_COMM_WORLD);
64   isThisRank=workfc.numpe+1;
65   if(*value==stats[0])
66     isThisRank=workfc.myrank;
67   MPI_Allreduce(&isThisRank,&statRanks[0],1,MPI_INT,MPI_MIN,MPI_COMM_WORLD);
68 
69   MPI_Allreduce(value,&stats[1],1,MPI_DOUBLE,MPI_MAX,MPI_COMM_WORLD);
70   isThisRank=workfc.numpe+1;
71   if(*value==stats[1])
72     isThisRank=workfc.myrank;
73   MPI_Allreduce(&isThisRank,&statRanks[1],1,MPI_INT,MPI_MIN,MPI_COMM_WORLD);
74 
75   MPI_Allreduce(value,&stats[2],1,MPI_DOUBLE,MPI_SUM,MPI_COMM_WORLD);
76   stats[2] /= workfc.numpe;
77 
78   double sqValue = (*value)*(*value), sqValueAvg = 0.;
79   MPI_Allreduce(&sqValue,&sqValueAvg,1,MPI_DOUBLE,MPI_SUM,MPI_COMM_WORLD);
80   sqValueAvg /= workfc.numpe;
81   // stats[3] = sqValueAvg;
82 
83   stats[3] = sqrt(sqValueAvg-stats[2]*stats[2]);
84 }
85 
86 void print_mesh_stats(void) {
87   int statRanks[2];
88   double iStats[4], rStats[4];
89 
90   igetMinMaxAvg(&conpar.nshg,iStats,statRanks);
91   if(workfc.myrank==workfc.master)
92     printf("nshg    : min [%d,%d], max[%d,%d] and avg[.,%d] (rms=%d)\n",statRanks[0],(int)iStats[0],statRanks[1],(int)iStats[1],(int)iStats[2],(int)iStats[3]);
93   igetMinMaxAvg(&conpar.numel,iStats,statRanks);
94   if(workfc.myrank==workfc.master)
95     printf("numel   : min [%d,%d], max[%d,%d] and avg[.,%d] (rms=%d)\n",statRanks[0],(int)iStats[0],statRanks[1],(int)iStats[1],(int)iStats[2],(int)iStats[3]);
96   igetMinMaxAvg(&conpar.numelb,iStats,statRanks);
97   if(workfc.myrank==workfc.master)
98     printf("numelb  : min [%d,%d], max[%d,%d] and avg[.,%d] (rms=%d)\n",statRanks[0],(int)iStats[0],statRanks[1],(int)iStats[1],(int)iStats[2],(int)iStats[3]);
99   igetMinMaxAvg(&conpar.nnz_tot,iStats,statRanks);
100   if(workfc.myrank==workfc.master) {
101     printf("nnz_tot : min [%d,%d], max[%d,%d] and avg[.,%d] (rms=%d)\n",statRanks[0],(int)iStats[0],statRanks[1],(int)iStats[1],(int)iStats[2],(int)iStats[3]);
102     printf("\n");
103   }
104 }
105 
106 void print_mpi_stats(void) {
107   int statRanks[2];
108   double iStats[4], rStats[4];
109 
110 // NS equations
111   igetMinMaxAvg(&mpistats.iISend,iStats,statRanks);
112   if(workfc.myrank==workfc.master)
113     printf("iISend : min [%d,%d], max[%d,%d] and avg[.,%d] (rms=%d)\n",statRanks[0],(int)iStats[0],statRanks[1],(int)iStats[1],(int)iStats[2],(int)iStats[3]);
114   igetMinMaxAvg(&mpistats.iIRecv,iStats,statRanks);
115   if(workfc.myrank==workfc.master)
116     printf("iIRecv : min [%d,%d], max[%d,%d] and avg[.,%d] (rms=%d)\n",statRanks[0],(int)iStats[0],statRanks[1],(int)iStats[1],(int)iStats[2],(int)iStats[3]);
117   igetMinMaxAvg(&mpistats.iWaitAll,iStats,statRanks);
118   if(workfc.myrank==workfc.master)
119     printf("iWtAll : min [%d,%d], max[%d,%d] and avg[.,%d] (rms=%d)\n",statRanks[0],(int)iStats[0],statRanks[1],(int)iStats[1],(int)iStats[2],(int)iStats[3]);
120   igetMinMaxAvg(&mpistats.iAllR,iStats,statRanks);
121   if(workfc.myrank==workfc.master)
122     printf("iAllR  : min [%d,%d], max[%d,%d] and avg[.,%d] (rms=%d)\n",statRanks[0],(int)iStats[0],statRanks[1],(int)iStats[1],(int)iStats[2],(int)iStats[3]);
123 
124   rgetMinMaxAvg(&mpistats.rISend,rStats,statRanks);
125   if(workfc.myrank==workfc.master)
126     printf("rISend : min [%d,%2.5f], max[%d,%2.5f] and avg[.,%2.5f] (rms=%2.5f)\n",statRanks[0],rStats[0],statRanks[1],rStats[1],rStats[2],rStats[3]);
127   rgetMinMaxAvg(&mpistats.rIRecv,rStats,statRanks);
128   if(workfc.myrank==workfc.master)
129     printf("rIRecv : min [%d,%2.5f], max[%d,%2.5f] and avg[.,%2.5f] (rms=%2.5f)\n",statRanks[0],rStats[0],statRanks[1],rStats[1],rStats[2],rStats[3]);
130   rgetMinMaxAvg(&mpistats.rWaitAll,rStats,statRanks);
131   if(workfc.myrank==workfc.master)
132     printf("rWtAll : min [%d,%2.5f], max[%d,%2.5f] and avg[.,%2.5f] (rms=%2.5f)\n",statRanks[0],rStats[0],statRanks[1],rStats[1],rStats[2],rStats[3]);
133   rgetMinMaxAvg(&mpistats.rCommu,rStats,statRanks);
134   if(workfc.myrank==workfc.master)
135     printf("rCommu : min [%d,%2.5f], max[%d,%2.5f] and avg[.,%2.5f] (rms=%2.5f)\n",statRanks[0],rStats[0],statRanks[1],rStats[1],rStats[2],rStats[3]);
136   rgetMinMaxAvg(&mpistats.rAllR,rStats,statRanks);
137   if(workfc.myrank==workfc.master) {
138     printf("rAllR  : min [%d,%2.5f], max[%d,%2.5f] and avg[.,%2.5f] (rms=%2.5f)\n",statRanks[0],rStats[0],statRanks[1],rStats[1],rStats[2],rStats[3]);
139     printf("\n");
140   }
141 // Scalars
142   igetMinMaxAvg(&mpistats.iISendScal,iStats,statRanks);
143   if(workfc.myrank==workfc.master)
144     printf("iISendScal : min [%d,%d], max[%d,%d] and avg[.,%d] (rms=%d)\n",statRanks[0],(int)iStats[0],statRanks[1],(int)iStats[1],(int)iStats[2],(int)iStats[3]);
145   igetMinMaxAvg(&mpistats.iIRecvScal,iStats,statRanks);
146   if(workfc.myrank==workfc.master)
147     printf("iIRecvScal : min [%d,%d], max[%d,%d] and avg[.,%d] (rms=%d)\n",statRanks[0],(int)iStats[0],statRanks[1],(int)iStats[1],(int)iStats[2],(int)iStats[3]);
148   igetMinMaxAvg(&mpistats.iWaitAllScal,iStats,statRanks);
149   if(workfc.myrank==workfc.master)
150     printf("iWtAllScal : min [%d,%d], max[%d,%d] and avg[.,%d] (rms=%d)\n",statRanks[0],(int)iStats[0],statRanks[1],(int)iStats[1],(int)iStats[2],(int)iStats[3]);
151   igetMinMaxAvg(&mpistats.iAllRScal,iStats,statRanks);
152   if(workfc.myrank==workfc.master)
153     printf("iAllRScal : min [%d,%d], max[%d,%d] and avg[.,%d] (rms=%d)\n",statRanks[0],(int)iStats[0],statRanks[1],(int)iStats[1],(int)iStats[2],(int)iStats[3]);
154 
155   rgetMinMaxAvg(&mpistats.rISendScal,rStats,statRanks);
156   if(workfc.myrank==workfc.master)
157     printf("rISendScal : min [%d,%2.5f], max[%d,%2.5f] and avg[.,%2.5f] (rms=%2.5f)\n",statRanks[0],rStats[0],statRanks[1],rStats[1],rStats[2],rStats[3]);
158   rgetMinMaxAvg(&mpistats.rIRecvScal,rStats,statRanks);
159   if(workfc.myrank==workfc.master)
160     printf("rIRecvScal : min [%d,%2.5f], max[%d,%2.5f] and avg[.,%2.5f] (rms=%2.5f)\n",statRanks[0],rStats[0],statRanks[1],rStats[1],rStats[2],rStats[3]);
161   rgetMinMaxAvg(&mpistats.rWaitAllScal,rStats,statRanks);
162   if(workfc.myrank==workfc.master)
163     printf("rWtAllScal : min [%d,%2.5f], max[%d,%2.5f] and avg[.,%2.5f] (rms=%2.5f)\n",statRanks[0],rStats[0],statRanks[1],rStats[1],rStats[2],rStats[3]);
164   rgetMinMaxAvg(&mpistats.rCommuScal,rStats,statRanks);
165   if(workfc.myrank==workfc.master)
166     printf("rCommuScal : min [%d,%2.5f], max[%d,%2.5f] and avg[.,%2.5f] (rms=%2.5f)\n",statRanks[0],rStats[0],statRanks[1],rStats[1],rStats[2],rStats[3]);
167   rgetMinMaxAvg(&mpistats.rAllRScal,rStats,statRanks);
168   if(workfc.myrank==workfc.master)
169     printf("rAllRScal  : min [%d,%2.5f], max[%d,%2.5f] and avg[.,%2.5f] (rms=%2.5f)\n",statRanks[0],rStats[0],statRanks[1],rStats[1],rStats[2],rStats[3]);
170 
171 
172 }
173 
174 //void print_system_stats(double tcorecp[2]) {
175 void print_system_stats(double *tcorecp, double *tcorecpscal) {
176   int statRanks[2];
177   double iStats[4], rStats[4];
178   double syst_assembly, syst_solve;
179 
180 // NS equations
181   syst_assembly = tcorecp[0];
182   syst_solve = tcorecp[1];
183 
184   rgetMinMaxAvg(&syst_assembly,rStats,statRanks);
185   if(workfc.myrank==workfc.master)
186     printf("Elm. form. : min [%d,%2.5f], max[%d,%2.5f] and avg[.,%2.5f] (rms=%2.5f)\n",statRanks[0],rStats[0],statRanks[1],rStats[1],rStats[2],rStats[3]);
187 
188   rgetMinMaxAvg(&syst_solve,rStats,statRanks);
189   if(workfc.myrank==workfc.master)
190     printf("Lin. alg. sol : min [%d,%2.5f], max[%d,%2.5f] and avg[.,%2.5f] (rms=%2.5f)\n",statRanks[0],rStats[0],statRanks[1],rStats[1],rStats[2],rStats[3]);
191 
192 // Scalars
193   syst_assembly = tcorecpscal[0];
194   syst_solve = tcorecpscal[1];
195 
196   rgetMinMaxAvg(&syst_assembly,rStats,statRanks);
197   if(workfc.myrank==workfc.master)
198     printf("Elm. form. Scal. : min [%d,%2.5f], max[%d,%2.5f] and avg[.,%2.5f] (rms=%2.5f)\n",statRanks[0],rStats[0],statRanks[1],rStats[1],rStats[2],rStats[3]);
199 
200   rgetMinMaxAvg(&syst_solve,rStats,statRanks);
201   if(workfc.myrank==workfc.master) {
202     printf("Lin. alg. sol Scal. : min [%d,%2.5f], max[%d,%2.5f] and avg[.,%2.5f] (rms=%2.5f)\n",statRanks[0],rStats[0],statRanks[1],rStats[1],rStats[2],rStats[3]);
203     printf("\n");
204   }
205   //printf("rank %d - syst_assembly %f - syst_solve %f\n",workfc.myrank,syst_assembly,syst_solve);
206 }
207 
208 
209 
210 void countfieldstowriterestart()
211 {
212   int nfields;
213 
214 //     printf("TEST: %d %d %d %d %d\n",timdat.istep,timdat.itseq,inpdat.nstep[0],inpdat.nstep[1],timdat.lstep);
215 
216   nfields = 2; //solution, time derivatives
217 
218   if(outpar.ivort == 1){
219     nfields++; //vorticity
220   }
221 
222   if(abs(turbvar.itwmod) != 1 && outpar.iowflux == 1) {
223     nfields++; //instantaneous wss in bflux.f
224   }
225 
226 //   if(ideformwall.eq.1) not handled yet
227 
228   if(timdat.istep == inpdat.nstep[timdat.itseq-1]){ //Last time step of the computation
229 
230     //projection vectors and pressure projection vectors (call saveLesRestart in itrdrv)
231     nfields = nfields +2;
232 
233     //if Print Error Indicators = true (call write_error in itrdrv)
234     if(turbvar.ierrcalc == 1){
235       nfields++;
236     }
237 
238     //if Print ybar = True (call write_field(myrank,'a','ybar',4,... in itrdrv)
239     if(outpar.ioybar == 1){
240       nfields++;  //ybar
241 
242       //phase average fields
243       if(outpar.nphasesincycle >0) {
244         nfields = nfields + outpar.nphasesincycle;
245       }
246 
247       if(abs(turbvar.itwmod) != 1 && outpar.iowflux == 1) {
248         nfields++; //wssbar
249       }
250 
251     }
252 
253     if(turbvari.irans < 0) {
254       nfields++; //dwal
255     }
256 
257   }
258 
259   outpar.nsynciofieldswriterestart = nfields;
260 
261   if(workfc.myrank == 0) {
262     printf("Number of fields to write in restart files: %d\n", nfields);
263   }
264 }
265 
266 
267 void
268 Write_Restart(  int* pid,
269                 int* stepno,
270                 int* nshg,
271                 int* numVars,
272                 double* array1,
273                 double* array2 ) {
274 
275     char fname[255];
276     char rfile[60];
277     char existingfile[30], linkfile[30];
278     int irstou;
279     int magic_number = 362436;
280     int* mptr = &magic_number;
281 //    time_t timenow = time ( &timenow);
282     double version=0.0;
283     int isize, nitems;
284     int iarray[10];
285 
286     /*sprintf(rfile,"restart.%d.%d",*stepno,*pid+1);
287     openfile_(rfile,"write", &irstou);
288 
289     // writing the top ascii header for the restart file
290 
291     writestring_( &irstou,"# PHASTA Input File Version 2.0\n");
292     writestring_( &irstou,
293                   "# format \"keyphrase : sizeofnextblock usual headers\"\n");
294 
295     bzero( (void*)fname, 255 );
296     sprintf(fname,"# Output generated by phasta version (NOT YET CURRENT): %lf \n", version);
297     writestring_( &irstou, fname );
298 
299     bzero( (void*)fname, 255 );
300     gethostname(fname,255);
301     writestring_( &irstou,"# This result was produced on: ");
302     writestring_( &irstou, fname );
303     writestring_( &irstou,"\n");
304 
305     bzero( (void*)fname, 255 );
306     sprintf(fname,"# %s\n", ctime( &timenow ));
307     writestring_( &irstou, fname );
308 
309     isize = 1;
310     nitems = 1;
311     iarray[ 0 ] = 1;
312     writeheader_( &irstou, "byteorder magic number ",
313                   (void*)iarray, &nitems, &isize, "integer", phasta_iotype );
314 
315     nitems = 1;
316     writedatablock_( &irstou, "byteorder magic number ",
317                      (void*)mptr, &nitems, "integer", phasta_iotype );
318 
319 
320     bzero( (void*)fname, 255 );
321     sprintf(fname,"number of modes : < 0 > %d\n", *nshg);
322     writestring_( &irstou, fname );
323 
324     bzero( (void*)fname, 255 );
325     sprintf(fname,"number of variables : < 0 > %d\n", *numVars);
326     writestring_( &irstou, fname );
327 
328 
329     isize = (*nshg)*(*numVars);
330     nitems = 3;
331     iarray[ 0 ] = (*nshg);
332     iarray[ 1 ] = (*numVars);
333     iarray[ 2 ] = (*stepno);
334     writeheader_( &irstou, "solution ",
335                   (void*)iarray, &nitems, &isize, "double", phasta_iotype );
336 
337 
338     nitems = (*nshg)*(*numVars);
339     writedatablock_( &irstou, "solution ",
340                      (void*)(array1), &nitems, "double", phasta_iotype );
341 
342 
343 
344     nitems = 3;
345     writeheader_( &irstou, "time derivative of solution ",
346                   (void*)iarray, &nitems, &isize, "double", phasta_iotype );
347 
348 
349     nitems = (*nshg)*(*numVars);
350     writedatablock_( &irstou, "time derivative of solution ",
351                      (void*)(array2), &nitems, "double", phasta_iotype );
352 
353 
354     closefile_( &irstou, "write" );
355     */
356     //MPI_Barrier(MPI_COMM_WORLD);
357 
358     /////////////////////////////// Start of writing using new-lib ////////////////////////////
359 
360 //MR CHANGE
361     int nfiles;
362     int nfields;
363     int numparts;
364     int irank;
365     int nprocs;
366 
367     //  First, count the number of fields to write and store the result in
368     countfieldstowriterestart();
369 
370     //  Retrieve and compute the parameters required for SyncIO
371     nfiles = outpar.nsynciofiles;
372     nfields = outpar.nsynciofieldswriterestart;
373     numparts = workfc.numpe;
374     irank = *pid; //workfc.myrank;
375     nprocs = workfc.numpe;
376 //MR CHANGE END
377     int nppf = numparts/nfiles;
378     int GPID;
379 
380     // Calculate number of parts each proc deal with and where it start and end ...
381     int nppp = numparts/nprocs;// nppp : Number of parts per proc ...
382     int startpart = irank * nppp +1;// Part id from which I (myrank) start ...
383     int endpart = startpart + nppp - 1;// Part id to which I (myrank) end ...
384 
385     int descriptor;
386     char filename[255],path[255],fieldtag_s[255];
387     bzero((void*)filename,255);
388     bzero((void*)fieldtag_s,255);
389 
390     phio_restartname(stepno, filename);
391     phio_openfile_write(filename, &nfiles, &nfields, &nppf, &f_descriptor);
392 
393 //MR CHANGE
394 //  Measure the time - End of timer
395 //    MPI_Barrier(MPI_COMM_WORLD);
396 //    timer_end = rdtsc();
397 //    time_span=(double)((timer_end-timer_start)/cpu_speed);
398 //    if (*pid==0) {
399 //      printf("Time: 'openfile' of %s with %d fields and %d files is:    %f s\n",filename,nfields,nfiles,time_span);
400 //      printf("*****************************\n");
401 //    }
402 //MR CHANGE END
403 
404     field_flag=0;
405 
406      int i;
407      for ( i = 0; i < nppp; i++) { //This loop is useful only if several parts per processor
408      // GPID : global part id, corresponds to rank ...
409         // e.g : (in this example)
410         // proc 0 : 1--4
411         // proc 1 : 5--8 ...
412         GPID = startpart + i;
413 
414         // Write solution field ...
415         sprintf(fieldtag_s,"solution@%d",GPID);
416 
417         isize = (*nshg)*(*numVars);
418         nitems = 3;
419         iarray[ 0 ] = (*nshg);
420         iarray[ 1 ] = (*numVars);
421         iarray[ 2 ] = (*stepno);
422 
423 //MR CHANGE
424 //  Measure the time - Start the timer
425 //        MPI_Barrier(MPI_COMM_WORLD);
426 //        timer_start = rdtsc();
427 //MR CHANGE END
428 
429         phio_writeheader( &f_descriptor, "solution", (void*)iarray, &nitems,
430             &isize, "double", phasta_iotype);
431 
432 //MR CHANGE
433 //  Measure the time - End of timer
434 //        MPI_Barrier(MPI_COMM_WORLD);
435 //        timer_end = rdtsc();
436 //        time_span=(double)((timer_end-timer_start)/cpu_speed);
437 //        if (*pid==0) {
438 //          printf("\n*****************************\n");
439 //          printf("Time: header for 'Solution':    %f s\n",time_span);
440 //        }
441 //MR CHANGE END
442 
443         nitems = (*nshg)*(*numVars);
444 
445 //MR CHANGE
446 //  Measure the time - Start the timer
447 //        MPI_Barrier(MPI_COMM_WORLD);
448 //        timer_start = rdtsc();
449 //MR CHANGE END
450 
451         writedatablock( &f_descriptor, fieldtag_s, (void*)(array1), &isize, "double", phasta_iotype );
452 
453 //MR CHANGE
454 //  Measure the time - End of timer
455 //        MPI_Barrier(MPI_COMM_WORLD);
456 //        timer_end = rdtsc();
457 //        time_span=(double)((timer_end-timer_start)/cpu_speed);
458 
459 //        int isizemin,isizemax,isizetot;
460 //        double sizemin,sizemax,sizeavg,sizetot,rate;
461 
462 //        MPI_Allreduce(&isize,&isizemin,1,MPI_INT,MPI_MIN,MPI_COMM_WORLD);
463 //        MPI_Allreduce(&isize,&isizemax,1,MPI_INT,MPI_MAX,MPI_COMM_WORLD);
464 //        MPI_Allreduce(&isize,&isizetot,1,MPI_INT,MPI_SUM,MPI_COMM_WORLD);
465 
466 //        sizemin=(double)(8.0*isizemin/1024.0/1024.0);
467 //        sizemax=(double)(8.0*isizemax/1024.0/1024.0);
468 //        sizetot=(double)(8.0*isizetot/1024.0/1024.0);
469 //        sizeavg=(double)(1.0*sizetot/workfc.numpe);
470 //        rate=(double)(1.0*sizetot/time_span);
471 
472 //        if (*pid==0) {
473 //          printf("Time: block for 'Solution':    %f s\n",time_span);
474 //          printf("Time: block:   Min= %f MB; Max= %f MB; Avg= %f MB; Tot= %f MB; Rate= %f MB/s; \n",sizemin,sizemax,sizeavg,sizetot,rate);
475 
476 //        }
477 //MR CHANGE END
478 
479     }
480     field_flag++;
481 
482     for ( i = 0; i < nppp; i++) {
483 
484         // GPID : global part id, corresponds to rank ...
485         // e.g : (in this example)
486         // proc 0 : 1--4
487         // proc 1 : 5--8 ...
488         GPID = startpart + i;
489 
490         // Write solution field ...
491         sprintf(fieldtag_s,"time derivative of solution@%d",GPID);
492 
493         isize = (*nshg)*(*numVars);
494         nitems = 3;
495         iarray[ 0 ] = (*nshg);
496         iarray[ 1 ] = (*numVars);
497         iarray[ 2 ] = (*stepno);
498 
499 //MR CHANGE
500 //  Measure the time - Start the timer
501 //        MPI_Barrier(MPI_COMM_WORLD);
502 //        timer_start = rdtsc();
503 //MR CHANGE END
504 
505         phio_writeheader( &f_descriptor, "time derivative of solution",
506             (void*)iarray, &nitems, &isize, "double", phasta_iotype);
507 
508 //MR CHANGE
509 //  Measure the time - End of timer
510 //        MPI_Barrier(MPI_COMM_WORLD);
511 //        timer_end = rdtsc();
512 //        time_span=(double)((timer_end-timer_start)/cpu_speed);
513 //        if (*pid==0) {
514 //          printf("Time: header for 'Time Derivative of solution':    %f s\n",time_span);
515 //        }
516 //MR CHANGE END
517 
518         nitems = (*nshg)*(*numVars);
519 
520 //MR CHANGE
521 //  Measure the time - Start the timer
522 //        MPI_Barrier(MPI_COMM_WORLD);
523 //        timer_start = rdtsc();
524 //MR CHANGE END
525 
526         writedatablock( &f_descriptor, fieldtag_s, (void*)(array2), &isize, "double", phasta_iotype );
527 
528 //MR CHANGE
529 //  Measure the time - End of timer
530 //        MPI_Barrier(MPI_COMM_WORLD);
531 //        timer_end = rdtsc();
532 //        time_span=(double)((timer_end-timer_start)/cpu_speed);
533 
534 //        int isizemin,isizemax,isizetot;
535 //        double sizemin,sizemax,sizeavg,sizetot,rate;
536 
537 //        MPI_Allreduce(&isize,&isizemin,1,MPI_INT,MPI_MIN,MPI_COMM_WORLD);
538 //        MPI_Allreduce(&isize,&isizemax,1,MPI_INT,MPI_MAX,MPI_COMM_WORLD);
539 //        MPI_Allreduce(&isize,&isizetot,1,MPI_INT,MPI_SUM,MPI_COMM_WORLD);
540 
541 //        sizemin=(double)(8.0*isizemin/1024.0/1024.0);
542 //        sizemax=(double)(8.0*isizemax/1024.0/1024.0);
543 //        sizetot=(double)(8.0*isizetot/1024.0/1024.0);
544 //        sizeavg=sizetot/workfc.numpe;
545 //        rate=sizetot/time_span;
546 
547 //        if (*pid==0) {
548 //          printf("Time: block for 'Time Derivative of Solution':    %f s\n",time_span);
549 //          printf("Time: block:   Min= %f MB; Max= %f MB; Avg= %f MB; Tot= %f MB; Rate= %f MB/s; \n",sizemin,sizemax,sizeavg,sizetot,rate);
550 //          printf("*****************************\n");
551 
552 //        }
553 //MR CHANGE END
554 
555     }
556     field_flag++;
557 
558     if (field_flag==nfields){
559 
560 //MR CHANGE
561 //  Measure the time - Start the timer
562 //      MPI_Barrier(MPI_COMM_WORLD);
563 //      timer_start = rdtsc();
564 //MR CHANGE END
565 
566 //MR CHANGE
567 //    Measure the time - End of timer
568 //      MPI_Barrier(MPI_COMM_WORLD);
569 //      timer_end = rdtsc();
570 //      time_span=(double)((timer_end-timer_start)/cpu_speed);
571 //      if (*pid==0) {
572 //        printf("\n*****************************\n");
573 //        printf("Time: 'closefile' is:    %f s\n",time_span);
574 //      }
575 //MR CHANGE END
576 
577 //MR CHANGE
578 //  Measure the time - Start the timer
579 //      MPI_Barrier(MPI_COMM_WORLD);
580 //      timer_start = rdtsc();
581 //MR CHANGE END
582 
583       phio_closefile_write(&f_descriptor);
584 
585 //MR CHANGE
586 //    Measure the time - End of timer
587 //      MPI_Barrier(MPI_COMM_WORLD);
588 //      timer_end = rdtsc();
589 //      time_span=(double)((timer_end-timer_start)/cpu_speed);
590       if (*pid==0) {
591 //        printf("Time: 'finalizephmpiio' is:    %f s\n",time_span);
592 //        printf("Last field %d '%s' finished! \n",nfields, fieldtag_s);
593         printf("\n");
594 //        printf("*****************************\n");
595       }
596     }
597 //MR CHANGE END
598 
599 
600 
601     ///////////////////////////////////////////////////////////////////////////////////////////
602 
603     /* create a soft link of the restart we just wrote to restart.latest
604      this is the file the next run will always try to start from */
605 
606 /*    sprintf( linkfile, "restart.latest.%d", *pid+1 );
607     unlink( linkfile );
608     sprintf( existingfile, "restart.%d.%d", *stepno, *pid+1 );
609     link( existingfile, linkfile );
610 */
611 }
612 
613 void
614 Write_Error(  int* pid,
615               int* stepno,
616               int* nshg,
617               int* numVars,
618               double* array1 ) {
619 
620 
621     char fname[255];
622     char rfile[60];
623     int irstou;
624     int magic_number = 362436;
625     int* mptr = &magic_number;
626     //printf("Time is commented\n");
627     //time_t timenow = time ( &timenow);
628     //printf("Yes\n");
629     double version=0.0;
630     int isize, nitems;
631     int iarray[10];
632 
633     /*sprintf(rfile,"restart.%d.%d",*stepno,*pid+1);
634     openfile_(rfile,"append", &irstou);
635 
636     isize = (*nshg)*(*numVars);
637     nitems = 3;
638     iarray[ 0 ] = (*nshg);
639     iarray[ 1 ] = (*numVars);
640     iarray[ 2 ] = (*stepno);
641     writeheader_( &irstou, "errors", (void*)iarray, &nitems, &isize, "double", phasta_iotype );
642 
643 
644     nitems = (*nshg)*(*numVars);
645     writedatablock_( &irstou, "errors ", (void*)(array1), &nitems, "double", phasta_iotype );
646 
647     closefile_( &irstou, "append" );*/
648 
649     /////////////////////////////// Start of writing using new-lib ////////////////////////////
650 
651     int nfiles;
652     int nfields;
653     int numparts;
654     int irank;
655     int nprocs;
656 
657 //    unsigned long long timer_start;
658 //    unsigned long long timer_end;
659 //    double time_span;
660 
661     nfiles = outpar.nsynciofiles;
662     nfields = outpar.nsynciofieldswriterestart;
663     numparts = workfc.numpe;
664     irank = *pid; //workfc.myrank;
665     nprocs = workfc.numpe;
666 
667     int nppf = numparts/nfiles;
668     int GPID;
669 
670     // Calculate number of parts each  proc deal with and where it start and end ...
671     int nppp = numparts/nprocs;// nppp : Number of parts per proc ...
672     int startpart = irank * nppp +1;// Part id from which I (myrank) start ...
673     int endpart = startpart + nppp - 1;// Part id to which I (myrank) end ...
674 
675     field_flag++;
676 
677     char fieldtag[255];
678 
679     int i;
680     for ( i = 0; i < nppp; i++  ) {
681         GPID = startpart + i;
682         sprintf(fieldtag,"errors@%d",GPID);
683 
684         if(*pid==0) {
685 //          printf("\n*****************************\n");
686           printf("\n");
687           printf("The %d/%d th field to be written is '%s'\n",field_flag,nfields,fieldtag);
688         }
689 
690         isize = (*nshg)*(*numVars);
691         nitems = 3;
692         iarray[ 0 ] = (*nshg);
693         iarray[ 1 ] = (*numVars);
694         iarray[ 2 ] = (*stepno);
695 
696 //MR CHANGE
697 //  Measure the time - Start the timer
698 //        MPI_Barrier(MPI_COMM_WORLD);
699 //        timer_start = rdtsc();
700 //MR CHANGE END
701 
702         phio_writeheader( &f_descriptor, "errors", (void*)iarray, &nitems,
703             &isize, "double", phasta_iotype);
704 
705 //MR CHANGE
706 //  Measure the time - End of timer
707 //        MPI_Barrier(MPI_COMM_WORLD);
708 //        timer_end = rdtsc();
709 //        time_span=(double)((timer_end-timer_start)/cpu_speed);
710 //        if (*pid==0) {
711 //          printf("Time: header for 'error':    %f s\n",time_span);
712 //        }
713 //MR CHANGE END
714 
715 //MR CHANGE
716 //  Measure the time - Start the timer
717 //        MPI_Barrier(MPI_COMM_WORLD);
718 //        timer_start = rdtsc();
719 //MR CHANGE END
720 
721         writedatablock( &f_descriptor, fieldtag, (void*)array1, &isize, "double", phasta_iotype );
722 
723 //MR CHANGE
724 //  Measure the time - End of timer
725 //        MPI_Barrier(MPI_COMM_WORLD);
726 //        timer_end = rdtsc();
727 //        time_span=(double)((timer_end-timer_start)/cpu_speed);
728 
729 //        int isizemin,isizemax,isizetot;
730 //        double sizemin,sizemax,sizeavg,sizetot,rate;
731 
732 //        MPI_Allreduce(&isize,&isizemin,1,MPI_INT,MPI_MIN,MPI_COMM_WORLD);
733 //        MPI_Allreduce(&isize,&isizemax,1,MPI_INT,MPI_MAX,MPI_COMM_WORLD);
734 //        MPI_Allreduce(&isize,&isizetot,1,MPI_INT,MPI_SUM,MPI_COMM_WORLD);
735 
736 //        sizemin=(double)(8.0*isizemin/1024.0/1024.0);
737 //        sizemax=(double)(8.0*isizemax/1024.0/1024.0);
738 //        sizetot=(double)(8.0*isizetot/1024.0/1024.0);
739 //        sizeavg=sizetot/workfc.numpe;
740 //        rate=sizetot/time_span;
741 
742 //        if (*pid==0) {
743 //          printf("Time: block for 'error':    %f s\n",time_span);
744 //          printf("Time: block:   Min= %f MB; Max= %f MB; Avg= %f MB; Tot= %f MB; Rate= %f MB/s; \n",sizemin,sizemax,sizeavg,sizetot,rate);
745 //          printf("*****************************\n");
746 //        }
747 //MR CHANGE END
748 
749     }
750 
751 //     MPI_Barrier(MPI_COMM_WORLD);
752 //     timer_end = rdtsc();
753 //     time_span=(double)(timer_end-timer_start)/cpu_speed;
754 
755 //     if (*pid==0) {
756 //         printf("Field 'error' written in:     %f s\n",time_span);
757 //         printf("Write field '%s' finished! \n",fieldtag);
758 //     }
759 
760     if (field_flag==nfields){
761 
762 //MR CHANGE
763 //  Measure the time - Start the timer
764 //      MPI_Barrier(MPI_COMM_WORLD);
765 //      timer_start = rdtsc();
766 //MR CHANGE END
767 
768 //MR CHANGE
769 //    Measure the time - End of timer
770 //      MPI_Barrier(MPI_COMM_WORLD);
771 //      timer_end = rdtsc();
772 //      time_span=(double)((timer_end-timer_start)/cpu_speed);
773 //      if (*pid==0) {
774 //        printf("\n*****************************\n");
775 //        printf("Time: 'closefile' is:    %f s\n",time_span);
776 //      }
777 //MR CHANGE END
778 
779 //MR CHANGE
780 //  Measure the time - Start the timer
781 //      MPI_Barrier(MPI_COMM_WORLD);
782 //      timer_start = rdtsc();
783 //MR CHANGE END
784 
785       phio_closefile_write(&f_descriptor);
786 
787 //MR CHANGE
788 //    Measure the time - End of timer
789 //      MPI_Barrier(MPI_COMM_WORLD);
790 //      timer_end = rdtsc();
791 //      time_span=(double)((timer_end-timer_start)/cpu_speed);
792       if (*pid==0) {
793 //        printf("Time: 'finalizephmpiio' is:    %f s\n",time_span);
794         printf("Last field %d '%s' finished! \n",nfields, fieldtag);
795         printf("\n");
796 //        printf("*****************************\n");
797       }
798     }
799 //MR CHANGE END
800 
801     ///////////////////////////////////////////////////////////////////////////////////////////
802 
803 
804 }
805 
806 
807 void
808 Write_Displ(  int* pid,
809               int* stepno,
810               int* nshg,
811               int* numVars,
812               double* array1 ) {
813   fprintf(stderr, "This function is dead...exiting\n");
814   exit(1);
815 }
816 
817 void
818 Write_Field(  int *pid,
819               char* filemode,
820               char* fieldtag,
821               int* tagsize,
822               void* array,
823               char* arraytype,
824               int* nshg,
825               int* numvars,
826               int* stepno) {
827 
828     //printf("Rank is %d, field is %s, tagsize is %d, nshg is %d, numvars is %d\n",*pid,fieldtag,*tagsize,*nshg,*numvars);
829 
830 //     char rfile[32];
831     // assuming restart.sn.(pid+1)
832 //     sprintf(rfile,"restart.%d.%d",*stepno,*pid+1);
833 
834     char *fieldlabel = (char *)malloc((*tagsize+1)*sizeof(char));
835     strncpy(fieldlabel, fieldtag, *tagsize);
836     fieldlabel[*tagsize] = '\0';
837 
838     int irstou;
839     int magic_number = 362436;
840     int* mptr = &magic_number;
841     double version=0.0;
842     int isize, nitems;
843     int iarray[10];
844 
845     char fmode[10];
846     if(!strncmp(filemode,"w",1))
847       strcpy(fmode,"write");
848     else // default is append
849       strcpy(fmode,"append");
850 
851     char datatype[10];
852     if(!strncmp(arraytype,"i",1))
853       strcpy(datatype,"int");
854     else // default is double
855       strcpy(datatype,"double");
856 
857 /*     openfile_(rfile, fmode, &irstou);
858 
859      nitems = 3; // assuming field will write 3 items in iarray
860      iarray[ 0 ] = (*nshg);
861      iarray[ 1 ] = (*numvars);
862      iarray[ 2 ] = (*stepno);
863 
864      isize = (*nshg)*(*numvars);
865      writeheader_( &irstou, fieldlabel, (void*)iarray, &nitems, &isize, datatype, phasta_iotype );
866 
867      nitems = (*nshg)*(*numvars);
868      writedatablock_( &irstou, fieldlabel, array, &nitems, datatype, phasta_iotype );
869      closefile_( &irstou, fmode);
870 */
871     /////////////////////////////// Start of writing using new-lib ////////////////////////////
872 
873     int nfiles;
874     int nfields;
875     int numparts;
876     int irank;
877     int nprocs;
878 
879 //    unsigned long long timer_start;
880 //    unsigned long long timer_end;
881 //    double time_span;
882 
883     nfiles = outpar.nsynciofiles;
884     nfields = outpar.nsynciofieldswriterestart;
885     numparts = workfc.numpe;
886     irank = *pid; //workfc.myrank;
887     nprocs = workfc.numpe;
888 
889     int nppf = numparts/nfiles;
890     int GPID;
891 
892     // Calculate number of parts each  proc deal with and where it start and end ...
893     int nppp = numparts/nprocs;// nppp : Number of parts per proc ...
894     int startpart = irank * nppp +1;// Part id from which I (myrank) start ...
895     int endpart = startpart + nppp - 1;// Part id to which I (myrank) end ...
896 
897     char filename[255],path[255],fieldtag_s[255];
898     bzero((void*)filename,255);
899     bzero((void*)fieldtag_s,255);
900 
901     strncpy(fieldlabel, fieldtag, *tagsize);
902 
903     field_flag++;
904     if(*pid==0) {
905 //      printf("\n*****************************\n");
906       printf("\n");
907       printf("The %d/%d th field to be written is '%s'\n",field_flag,nfields,fieldlabel);
908     }
909 
910     sprintf(filename,"restart-dat.%d.%d",*stepno,((int)(irank/(nprocs/nfiles))+1));
911 
912 //     MPI_Barrier(MPI_COMM_WORLD);
913 //     timer_start = rdtsc();
914 
915     int i;
916     for ( i = 0; i < nppp; i++  ) {
917         GPID = startpart + i;
918 
919         // Write solution field ...
920         sprintf(fieldtag_s,"%s@%d",fieldlabel,GPID);
921 
922         isize = (*nshg)*(*numvars);
923         nitems = 3;
924         iarray[ 0 ] = (*nshg);
925         iarray[ 1 ] = (*numvars);
926         iarray[ 2 ] = (*stepno);
927 
928 //MR CHANGE
929 //  Measure the time - Start the timer
930 //        MPI_Barrier(MPI_COMM_WORLD);
931 //        timer_start = rdtsc();
932 //MR CHANGE END
933 
934         phio_writeheader( &f_descriptor, fieldlabel, (void*)iarray, &nitems,
935             &isize, datatype, phasta_iotype);
936 
937 //MR CHANGE
938 //  Measure the time - End of timer
939 //        MPI_Barrier(MPI_COMM_WORLD);
940 //        timer_end = rdtsc();
941 //        time_span=(double)((timer_end-timer_start)/cpu_speed);
942 //        if (*pid==0) {
943 //          printf("Time: header for '%s':    %f s\n",fieldtag_s,time_span);
944 //        }
945 //MR CHANGE END
946 
947         nitems = (*nshg)*(*numvars);
948 
949 //MR CHANGE
950 //  Measure the time - Start the timer
951 //        MPI_Barrier(MPI_COMM_WORLD);
952 //        timer_start = rdtsc();
953 //MR CHANGE END
954 
955         writedatablock( &f_descriptor, fieldtag_s, array, &isize, datatype, phasta_iotype );
956 
957 //MR CHANGE
958 //  Measure the time - End of timer
959 //        MPI_Barrier(MPI_COMM_WORLD);
960 //        timer_end = rdtsc();
961 //        time_span=(double)((timer_end-timer_start)/cpu_speed);
962 
963 //        int isizemin,isizemax,isizetot;
964 //        double sizemin,sizemax,sizeavg,sizetot,rate;
965 
966 //        MPI_Allreduce(&isize,&isizemin,1,MPI_INT,MPI_MIN,MPI_COMM_WORLD);
967 //        MPI_Allreduce(&isize,&isizemax,1,MPI_INT,MPI_MAX,MPI_COMM_WORLD);
968 //        MPI_Allreduce(&isize,&isizetot,1,MPI_INT,MPI_SUM,MPI_COMM_WORLD);
969 
970 //        sizemin=(double)(8.0*isizemin/1024.0/1024.0);
971 //        sizemax=(double)(8.0*isizemax/1024.0/1024.0);
972 //        sizetot=(double)(8.0*isizetot/1024.0/1024.0);
973 //        sizeavg=sizetot/workfc.numpe;
974 //        rate=sizetot/time_span;
975 
976 //        if (*pid==0) {
977 //          printf("Time: block for '%s':    %f s\n",fieldtag_s,time_span);
978 //          printf("Time: block:   Min= %f MB; Max= %f MB; Avg= %f MB; Tot= %f MB; Rate= %f MB/s; \n",sizemin,sizemax,sizeavg,sizetot,rate);
979 //          printf("*****************************\n");
980 //        }
981 //MR CHANGE END
982 
983     }
984 
985 //     MPI_Barrier(MPI_COMM_WORLD);
986 //     timer_end = rdtsc();
987 //     time_span=(double)(timer_end-timer_start)/cpu_speed;
988 
989 //     if (*pid==0) {
990 //         printf("Field '%s' written in:     %f s\n",fieldtag,time_span);
991 //         printf("Write field '%s' finished! \n",fieldtag_s);
992 //     }
993 
994 //     if (field_flag==nfields){
995 //       closefile_(&f_descriptor, "write");
996 //       finalizephmpiio_(&f_descriptor);
997 //       if(*pid==0) {
998 //         printf("Last field %d '%s' finished! \n",nfields, fieldtag_s);
999 //         printf("\n*****************************\n");
1000 //       }
1001 //     }
1002 
1003     if (field_flag==nfields){
1004 
1005 //MR CHANGE
1006 //  Measure the time - Start the timer
1007 //      MPI_Barrier(MPI_COMM_WORLD);
1008 //      timer_start = rdtsc();
1009 //MR CHANGE END
1010 
1011 //MR CHANGE
1012 //    Measure the time - End of timer
1013 //      MPI_Barrier(MPI_COMM_WORLD);
1014 //      timer_end = rdtsc();
1015 //      time_span=(double)((timer_end-timer_start)/cpu_speed);
1016 //      if (*pid==0) {
1017 //        printf("\n*****************************\n");
1018 //        printf("Time: 'closefile' is:    %f s\n",time_span);
1019 //      }
1020 //MR CHANGE END
1021 
1022 //MR CHANGE
1023 //  Measure the time - Start the timer
1024 //      MPI_Barrier(MPI_COMM_WORLD);
1025 //      timer_start = rdtsc();
1026 //MR CHANGE END
1027 
1028       phio_closefile_write(&f_descriptor);
1029 
1030 //MR CHANGE
1031 //    Measure the time - End of timer
1032 //      MPI_Barrier(MPI_COMM_WORLD);
1033 //      timer_end = rdtsc();
1034 //      time_span=(double)((timer_end-timer_start)/cpu_speed);
1035       if (*pid==0) {
1036 //        printf("Time: 'finalizephmpiio' is:    %f s\n",time_span);
1037         printf("Last field %d '%s' finished! \n",nfields, fieldtag);
1038         printf("\n");
1039 //        printf("*****************************\n");
1040       }
1041     }
1042 //MR CHANGE END
1043 
1044     ///////////////////////////////////////////////////////////////////////////////////////////
1045 
1046     free(fieldlabel);
1047 }
1048 
1049 //MR CHANGE
1050 
1051 //MR CHANGE
1052 void
1053 Write_PhAvg2( int* pid,
1054               char* filemode,
1055               char* fieldtag,
1056               int* tagsize,
1057               int* iphase,
1058               int* nphasesincycle,
1059               void* array,
1060               char* arraytype,
1061               int* nshg,
1062               int* numvars,
1063               int* stepno) {
1064 
1065 //     char rfile[32];
1066     // assuming restart.sn.(pid+1)
1067 //     sprintf(rfile,"restart.%d.%d",*stepno,*pid+1);
1068 
1069     int addtagsize=0; // phase number is added to the name of the field
1070     if(*iphase<10)
1071       addtagsize=1;
1072     else if(*iphase<100)
1073       addtagsize=2;
1074     else if(*iphase<1000)
1075       addtagsize=3;
1076 
1077     int tagsize2;
1078     tagsize2=*tagsize+addtagsize;
1079 
1080 //     char *fieldlabel = (char *)malloc((*tagsize+1)*sizeof(char));
1081 //     strncpy(fieldlabel, fieldtag, *tagsize);
1082 //     fieldlabel[*tagsize] = '\0';
1083 
1084     char *fieldlabel = (char *)malloc((tagsize2+1)*sizeof(char));
1085     strncpy(fieldlabel, fieldtag, *tagsize);
1086     fieldlabel[tagsize2] = '\0';
1087 
1088     char straddtagsize[10];
1089     sprintf(straddtagsize,"%d",*iphase);
1090 
1091     if(*iphase<10) {
1092       fieldlabel[tagsize2-1]=straddtagsize[0];
1093     }
1094     else if(*iphase<100) {
1095       fieldlabel[tagsize2-2]=straddtagsize[0];
1096       fieldlabel[tagsize2-1]=straddtagsize[1];
1097     }
1098     else if(*iphase<1000) {
1099       fieldlabel[tagsize2-3]=straddtagsize[0];
1100       fieldlabel[tagsize2-2]=straddtagsize[1];
1101       fieldlabel[tagsize2-1]=straddtagsize[2];
1102     }
1103 
1104     int irstou;
1105     int magic_number = 362436;
1106     int* mptr = &magic_number;
1107     double version=0.0;
1108     int isize, nitems;
1109     int iarray[10];
1110 
1111     char fmode[10];
1112     if(!strncmp(filemode,"w",1))
1113       strcpy(fmode,"write");
1114     else // default is append
1115       strcpy(fmode,"append");
1116 
1117     char datatype[10];
1118     if(!strncmp(arraytype,"i",1))
1119       strcpy(datatype,"int");
1120     else // default is double
1121       strcpy(datatype,"double");
1122 
1123 //
1124 // //     if(*iphase==1) //open the file but then keep it open for the remaining cycles
1125 //     openfile_(rfile, fmode, &irstou);
1126 //
1127 // //     printf("iphase: %d - pid: %d - irstou %d\n",*iphase,*pid,irstou);
1128 //
1129 //
1130 //     nitems = 3; // assuming field will write 3 items in iarray
1131 //     iarray[ 0 ] = (*nshg);
1132 //     iarray[ 1 ] = (*numvars);
1133 //     iarray[ 2 ] = (*stepno);
1134 //
1135 //     isize = (*nshg)*(*numvars);
1136 //     writeheader_( &irstou, fieldlabel, (void*)iarray, &nitems, &isize, datatype, phasta_iotype );
1137 //
1138 //     nitems = (*nshg)*(*numvars);
1139 //     writedatablock_( &irstou, fieldlabel, array, &nitems, datatype, phasta_iotype );
1140 //
1141 // //     if(*iphase==*nphasesincycle) //close the file after nphasesincycle
1142 //       closefile_( &irstou, fmode);
1143 //
1144 
1145     /////////////////////////////// Start of writing using new-lib ////////////////////////////
1146 
1147     int nfiles;
1148     int nfields;
1149     int numparts;
1150     int irank;
1151     int nprocs;
1152 //    unsigned long long timer_start;
1153 //    unsigned long long timer_end;
1154 //    double time_span;
1155 
1156     nfiles = outpar.nsynciofiles;
1157     nfields = outpar.nsynciofieldswriterestart;
1158     numparts = workfc.numpe;
1159     irank = *pid; //workfc.myrank;
1160     nprocs = workfc.numpe;
1161 
1162     int nppf = numparts/nfiles;
1163     int GPID;
1164 
1165     // Calculate number of parts each  proc deal with and where it start and end ...
1166     int nppp = numparts/nprocs;// nppp : Number of parts per proc ...
1167     int startpart = irank * nppp +1;// Part id from which I (myrank) start ...
1168     int endpart = startpart + nppp - 1;// Part id to which I (myrank) end ...
1169 
1170     //int descriptor;
1171     char filename[255],path[255],fieldtag_s[255];
1172     bzero((void*)filename,255);
1173     bzero((void*)fieldtag_s,255);
1174 
1175 //     char * namer;
1176 //     namer = strtok(fieldlabel," ");
1177 //     strncpy(fieldlabel, fieldtag, *tagsize);
1178 
1179     field_flag++;
1180     if(*pid==0) {
1181 //      printf("\n*****************************\n");
1182       printf("\n");
1183       printf("The %d/%d th field to be written is '%s'\n",field_flag,nfields,fieldlabel);
1184     }
1185 
1186     sprintf(filename,"restart-dat.%d.%d",*stepno,((int)(irank/(nprocs/nfiles))+1));
1187 
1188     int i;
1189     for ( i = 0; i < nppp; i++  ) {
1190         GPID = startpart + i;
1191 
1192         // Write solution field ...
1193         sprintf(fieldtag_s,"%s@%d",fieldlabel,GPID);
1194 
1195         //printf("This is %d and fieldtag_s is %s \n",myrank,fieldtag_s);
1196 
1197         isize = (*nshg)*(*numvars);
1198         nitems = 3;
1199         iarray[ 0 ] = (*nshg);
1200         iarray[ 1 ] = (*numvars);
1201         iarray[ 2 ] = (*stepno);
1202 
1203 //MR CHANGE
1204 //  Measure the time - Start the timer
1205 //        MPI_Barrier(MPI_COMM_WORLD);
1206 //        timer_start = rdtsc();
1207 //MR CHANGE END
1208 
1209         phio_writeheader( &f_descriptor, fieldlabel, (void*)iarray, &nitems,
1210             &isize, "double", phasta_iotype);
1211 
1212 //MR CHANGE
1213 //  Measure the time - End of timer
1214 //        MPI_Barrier(MPI_COMM_WORLD);
1215 //        timer_end = rdtsc();
1216 //        time_span=(double)((timer_end-timer_start)/cpu_speed);
1217 //        if (*pid==0) {
1218 //          printf("Time: header for '%s':    %f s\n",fieldtag_s,time_span);
1219 //        }
1220 //MR CHANGE END
1221 
1222         nitems = (*nshg)*(*numvars);
1223 
1224 //MR CHANGE
1225 //  Measure the time - Start the timer
1226 //        MPI_Barrier(MPI_COMM_WORLD);
1227 //        timer_start = rdtsc();
1228 //MR CHANGE END
1229 
1230         writedatablock( &f_descriptor, fieldtag_s, array, &isize, "double", phasta_iotype );
1231 
1232 //MR CHANGE
1233 //  Measure the time - End of timer
1234 //        MPI_Barrier(MPI_COMM_WORLD);
1235 //        timer_end = rdtsc();
1236 //        time_span=(double)((timer_end-timer_start)/cpu_speed);
1237 
1238 //        int isizemin,isizemax,isizetot;
1239 //        double sizemin,sizemax,sizeavg,sizetot,rate;
1240 
1241 //        MPI_Allreduce(&isize,&isizemin,1,MPI_INT,MPI_MIN,MPI_COMM_WORLD);
1242 //        MPI_Allreduce(&isize,&isizemax,1,MPI_INT,MPI_MAX,MPI_COMM_WORLD);
1243 //        MPI_Allreduce(&isize,&isizetot,1,MPI_INT,MPI_SUM,MPI_COMM_WORLD);
1244 
1245 //        sizemin=(double)(8.0*isizemin/1024.0/1024.0);
1246 //        sizemax=(double)(8.0*isizemax/1024.0/1024.0);
1247 //        sizetot=(double)(8.0*isizetot/1024.0/1024.0);
1248 //        sizeavg=sizetot/workfc.numpe;
1249 //        rate=sizetot/time_span;
1250 
1251 //        if (*pid==0) {
1252 //          printf("Time: block for '%s':    %f s\n",fieldtag_s,time_span);
1253 //          printf("Time: block:   Min= %f MB; Max= %f MB; Avg= %f MB; Tot= %f MB; Rate= %f MB/s; \n",sizemin,sizemax,sizeavg,sizetot,rate);
1254 //          printf("*****************************\n");
1255 //        }
1256 //MR CHANGE END
1257 
1258     }
1259 
1260 //     if (*pid==0) {
1261 //         printf("Field '%s' written in:     %f s\n",fieldtag,time_span);
1262 //         printf("Write field '%s' finished! \n",fieldtag_s);
1263 //     }
1264 
1265 //
1266 //     if (field_flag==nfields){
1267 //       closefile_(&f_descriptor, "write");
1268 //       finalizephmpiio_(&f_descriptor);
1269 //       if(*pid==0) {
1270 //         printf("Last field %d '%s' finished! \n",nfields, fieldtag_s);
1271 //         printf("\n*****************************\n");
1272 //       }
1273 //     }
1274 
1275     if (field_flag==nfields){
1276 
1277 //MR CHANGE
1278 //  Measure the time - Start the timer
1279 //      MPI_Barrier(MPI_COMM_WORLD);
1280 //      timer_start = rdtsc();
1281 //MR CHANGE END
1282 
1283 //MR CHANGE
1284 //    Measure the time - End of timer
1285 //      MPI_Barrier(MPI_COMM_WORLD);
1286 //      timer_end = rdtsc();
1287 //      time_span=(double)((timer_end-timer_start)/cpu_speed);
1288 //     if (*pid==0) {
1289 //        printf("\n*****************************\n");
1290 //        printf("Time: 'closefile' is:    %f s\n",time_span);
1291 //      }
1292 //MR CHANGE END
1293 
1294 //MR CHANGE
1295 //  Measure the time - Start the timer
1296 //      MPI_Barrier(MPI_COMM_WORLD);
1297 //      timer_start = rdtsc();
1298 //MR CHANGE END
1299 
1300       phio_closefile_write(&f_descriptor);
1301 
1302 //MR CHANGE
1303 //    Measure the time - End of timer
1304 //      MPI_Barrier(MPI_COMM_WORLD);
1305 //      timer_end = rdtsc();
1306 //      time_span=(double)((timer_end-timer_start)/cpu_speed);
1307       if (*pid==0) {
1308 //        printf("Time: 'finalizephmpiio' is:    %f s\n",time_span);
1309 //        printf("Last field %d '%s' finished! \n",nfields, fieldtag);
1310         printf("\n");
1311 //        printf("*****************************\n");
1312       }
1313     }
1314 //MR CHANGE END
1315 
1316     ///////////////////////////////////////////////////////////////////////////////////////////
1317 
1318     free(fieldlabel);
1319 }
1320 
1321 
1322 void
1323 Write_d2wall(   int* pid,
1324                 int* numnp,
1325                 double* array1 ) {
1326 
1327 //    time_t timenow = time ( &timenow);
1328     int isize, nitems;
1329     int iarray[10];
1330 
1331 //    MPI_Barrier(MPI_COMM_WORLD);
1332 
1333     /////////////////////////////// Start of writing using new-lib ////////////////////////////
1334 
1335     int nfiles;
1336     int nfields;
1337     int numparts;
1338     int irank;
1339     int nprocs;
1340 
1341     //  First, count the number of fields to write and store the result in
1342     //countfieldstowriterestart();
1343 
1344     //  Retrieve and compute the parameters required for SyncIO
1345     nfiles = outpar.nsynciofiles;
1346     nfields = 1; //outpar.nsynciofieldswriterestart;  // Only the distance to the walls in d2wall
1347     numparts = workfc.numpe;
1348     irank = *pid; //workfc.myrank;
1349     nprocs = workfc.numpe;
1350     int nppf = numparts/nfiles;
1351     int GPID;
1352 
1353     // Calculate number of parts each proc deal with and where it start and end ...
1354     int nppp = numparts/nprocs;// nppp : Number of parts per proc ...
1355     int startpart = irank * nppp +1;// Part id from which I (myrank) start ...
1356     int endpart = startpart + nppp - 1;// Part id to which I (myrank) end ...
1357 
1358     int descriptor;
1359     char filename[255],path[255],fieldtag_s[255];
1360     bzero((void*)filename,255);
1361     bzero((void*)fieldtag_s,255);
1362 
1363     phio_openfile_write("d2wall.", &nfiles, &nfields, &nppf, &f_descriptor);
1364 
1365     field_flag=0;
1366 
1367      int i;
1368      for ( i = 0; i < nppp; i++) { //This loop is useful only if several parts per processor
1369      // GPID : global part id, corresponds to rank ...
1370         // e.g : (in this example)
1371         // proc 0 : 1--4
1372         // proc 1 : 5--8 ...
1373         GPID = startpart + i;
1374 
1375         // Write solution field ...
1376         sprintf(fieldtag_s,"d2wall@%d",GPID);
1377 
1378         isize = (*numnp);
1379         nitems = 2;
1380         iarray[ 0 ] = (*numnp);
1381         iarray[ 1 ] = 1; //numVars = 1
1382 
1383         phio_writeheader( &f_descriptor, "d2wall", (void*)iarray, &nitems,
1384             &isize, "double", phasta_iotype);
1385 
1386         //nitems = (*nshg)*(*numVars);
1387         //nitems = (*numnp);
1388 
1389         writedatablock( &f_descriptor, fieldtag_s, (void*)(array1), &isize, "double", phasta_iotype );
1390 
1391 
1392     }
1393     field_flag++;
1394 
1395     if (field_flag==nfields){
1396       phio_closefile_write(&f_descriptor);
1397       if (irank==0) {
1398         printf("\n");
1399       }
1400     }
1401 }
1402 
1403