xref: /phasta/phSolver/common/new_interface.c (revision b63f4427c8540c053341d248402d86b42f257463)
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         writeheader( &f_descriptor, fieldtag_s, (void*)iarray, &nitems, &isize, "double", phasta_iotype);
430 
431 //MR CHANGE
432 //  Measure the time - End of timer
433 //        MPI_Barrier(MPI_COMM_WORLD);
434 //        timer_end = rdtsc();
435 //        time_span=(double)((timer_end-timer_start)/cpu_speed);
436 //        if (*pid==0) {
437 //          printf("\n*****************************\n");
438 //          printf("Time: header for 'Solution':    %f s\n",time_span);
439 //        }
440 //MR CHANGE END
441 
442         nitems = (*nshg)*(*numVars);
443 
444 //MR CHANGE
445 //  Measure the time - Start the timer
446 //        MPI_Barrier(MPI_COMM_WORLD);
447 //        timer_start = rdtsc();
448 //MR CHANGE END
449 
450         writedatablock( &f_descriptor, fieldtag_s, (void*)(array1), &isize, "double", phasta_iotype );
451 
452 //MR CHANGE
453 //  Measure the time - End of timer
454 //        MPI_Barrier(MPI_COMM_WORLD);
455 //        timer_end = rdtsc();
456 //        time_span=(double)((timer_end-timer_start)/cpu_speed);
457 
458 //        int isizemin,isizemax,isizetot;
459 //        double sizemin,sizemax,sizeavg,sizetot,rate;
460 
461 //        MPI_Allreduce(&isize,&isizemin,1,MPI_INT,MPI_MIN,MPI_COMM_WORLD);
462 //        MPI_Allreduce(&isize,&isizemax,1,MPI_INT,MPI_MAX,MPI_COMM_WORLD);
463 //        MPI_Allreduce(&isize,&isizetot,1,MPI_INT,MPI_SUM,MPI_COMM_WORLD);
464 
465 //        sizemin=(double)(8.0*isizemin/1024.0/1024.0);
466 //        sizemax=(double)(8.0*isizemax/1024.0/1024.0);
467 //        sizetot=(double)(8.0*isizetot/1024.0/1024.0);
468 //        sizeavg=(double)(1.0*sizetot/workfc.numpe);
469 //        rate=(double)(1.0*sizetot/time_span);
470 
471 //        if (*pid==0) {
472 //          printf("Time: block for 'Solution':    %f s\n",time_span);
473 //          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);
474 
475 //        }
476 //MR CHANGE END
477 
478     }
479     field_flag++;
480 
481     for ( i = 0; i < nppp; i++) {
482 
483         // GPID : global part id, corresponds to rank ...
484         // e.g : (in this example)
485         // proc 0 : 1--4
486         // proc 1 : 5--8 ...
487         GPID = startpart + i;
488 
489         // Write solution field ...
490         sprintf(fieldtag_s,"time derivative of solution@%d",GPID);
491 
492         isize = (*nshg)*(*numVars);
493         nitems = 3;
494         iarray[ 0 ] = (*nshg);
495         iarray[ 1 ] = (*numVars);
496         iarray[ 2 ] = (*stepno);
497 
498 //MR CHANGE
499 //  Measure the time - Start the timer
500 //        MPI_Barrier(MPI_COMM_WORLD);
501 //        timer_start = rdtsc();
502 //MR CHANGE END
503 
504         writeheader( &f_descriptor, fieldtag_s, (void*)iarray, &nitems, &isize, "double", phasta_iotype);
505 
506 //MR CHANGE
507 //  Measure the time - End of timer
508 //        MPI_Barrier(MPI_COMM_WORLD);
509 //        timer_end = rdtsc();
510 //        time_span=(double)((timer_end-timer_start)/cpu_speed);
511 //        if (*pid==0) {
512 //          printf("Time: header for 'Time Derivative of solution':    %f s\n",time_span);
513 //        }
514 //MR CHANGE END
515 
516         nitems = (*nshg)*(*numVars);
517 
518 //MR CHANGE
519 //  Measure the time - Start the timer
520 //        MPI_Barrier(MPI_COMM_WORLD);
521 //        timer_start = rdtsc();
522 //MR CHANGE END
523 
524         writedatablock( &f_descriptor, fieldtag_s, (void*)(array2), &isize, "double", phasta_iotype );
525 
526 //MR CHANGE
527 //  Measure the time - End of timer
528 //        MPI_Barrier(MPI_COMM_WORLD);
529 //        timer_end = rdtsc();
530 //        time_span=(double)((timer_end-timer_start)/cpu_speed);
531 
532 //        int isizemin,isizemax,isizetot;
533 //        double sizemin,sizemax,sizeavg,sizetot,rate;
534 
535 //        MPI_Allreduce(&isize,&isizemin,1,MPI_INT,MPI_MIN,MPI_COMM_WORLD);
536 //        MPI_Allreduce(&isize,&isizemax,1,MPI_INT,MPI_MAX,MPI_COMM_WORLD);
537 //        MPI_Allreduce(&isize,&isizetot,1,MPI_INT,MPI_SUM,MPI_COMM_WORLD);
538 
539 //        sizemin=(double)(8.0*isizemin/1024.0/1024.0);
540 //        sizemax=(double)(8.0*isizemax/1024.0/1024.0);
541 //        sizetot=(double)(8.0*isizetot/1024.0/1024.0);
542 //        sizeavg=sizetot/workfc.numpe;
543 //        rate=sizetot/time_span;
544 
545 //        if (*pid==0) {
546 //          printf("Time: block for 'Time Derivative of Solution':    %f s\n",time_span);
547 //          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);
548 //          printf("*****************************\n");
549 
550 //        }
551 //MR CHANGE END
552 
553     }
554     field_flag++;
555 
556     if (field_flag==nfields){
557 
558 //MR CHANGE
559 //  Measure the time - Start the timer
560 //      MPI_Barrier(MPI_COMM_WORLD);
561 //      timer_start = rdtsc();
562 //MR CHANGE END
563 
564 //MR CHANGE
565 //    Measure the time - End of timer
566 //      MPI_Barrier(MPI_COMM_WORLD);
567 //      timer_end = rdtsc();
568 //      time_span=(double)((timer_end-timer_start)/cpu_speed);
569 //      if (*pid==0) {
570 //        printf("\n*****************************\n");
571 //        printf("Time: 'closefile' is:    %f s\n",time_span);
572 //      }
573 //MR CHANGE END
574 
575 //MR CHANGE
576 //  Measure the time - Start the timer
577 //      MPI_Barrier(MPI_COMM_WORLD);
578 //      timer_start = rdtsc();
579 //MR CHANGE END
580 
581       phio_closefile_write(&f_descriptor);
582 
583 //MR CHANGE
584 //    Measure the time - End of timer
585 //      MPI_Barrier(MPI_COMM_WORLD);
586 //      timer_end = rdtsc();
587 //      time_span=(double)((timer_end-timer_start)/cpu_speed);
588       if (*pid==0) {
589 //        printf("Time: 'finalizephmpiio' is:    %f s\n",time_span);
590 //        printf("Last field %d '%s' finished! \n",nfields, fieldtag_s);
591         printf("\n");
592 //        printf("*****************************\n");
593       }
594     }
595 //MR CHANGE END
596 
597 
598 
599     ///////////////////////////////////////////////////////////////////////////////////////////
600 
601     /* create a soft link of the restart we just wrote to restart.latest
602      this is the file the next run will always try to start from */
603 
604 /*    sprintf( linkfile, "restart.latest.%d", *pid+1 );
605     unlink( linkfile );
606     sprintf( existingfile, "restart.%d.%d", *stepno, *pid+1 );
607     link( existingfile, linkfile );
608 */
609 }
610 
611 void
612 Write_Error(  int* pid,
613               int* stepno,
614               int* nshg,
615               int* numVars,
616               double* array1 ) {
617 
618 
619     char fname[255];
620     char rfile[60];
621     int irstou;
622     int magic_number = 362436;
623     int* mptr = &magic_number;
624     //printf("Time is commented\n");
625     //time_t timenow = time ( &timenow);
626     //printf("Yes\n");
627     double version=0.0;
628     int isize, nitems;
629     int iarray[10];
630 
631     /*sprintf(rfile,"restart.%d.%d",*stepno,*pid+1);
632     openfile_(rfile,"append", &irstou);
633 
634     isize = (*nshg)*(*numVars);
635     nitems = 3;
636     iarray[ 0 ] = (*nshg);
637     iarray[ 1 ] = (*numVars);
638     iarray[ 2 ] = (*stepno);
639     writeheader_( &irstou, "errors", (void*)iarray, &nitems, &isize, "double", phasta_iotype );
640 
641 
642     nitems = (*nshg)*(*numVars);
643     writedatablock_( &irstou, "errors ", (void*)(array1), &nitems, "double", phasta_iotype );
644 
645     closefile_( &irstou, "append" );*/
646 
647     /////////////////////////////// Start of writing using new-lib ////////////////////////////
648 
649     int nfiles;
650     int nfields;
651     int numparts;
652     int irank;
653     int nprocs;
654 
655 //    unsigned long long timer_start;
656 //    unsigned long long timer_end;
657 //    double time_span;
658 
659     nfiles = outpar.nsynciofiles;
660     nfields = outpar.nsynciofieldswriterestart;
661     numparts = workfc.numpe;
662     irank = *pid; //workfc.myrank;
663     nprocs = workfc.numpe;
664 
665     int nppf = numparts/nfiles;
666     int GPID;
667 
668     // Calculate number of parts each  proc deal with and where it start and end ...
669     int nppp = numparts/nprocs;// nppp : Number of parts per proc ...
670     int startpart = irank * nppp +1;// Part id from which I (myrank) start ...
671     int endpart = startpart + nppp - 1;// Part id to which I (myrank) end ...
672 
673     field_flag++;
674 
675     char fieldtag[255];
676 
677     int i;
678     for ( i = 0; i < nppp; i++  ) {
679         GPID = startpart + i;
680         sprintf(fieldtag,"errors@%d",GPID);
681 
682         if(*pid==0) {
683 //          printf("\n*****************************\n");
684           printf("\n");
685           printf("The %d/%d th field to be written is '%s'\n",field_flag,nfields,fieldtag);
686         }
687 
688         isize = (*nshg)*(*numVars);
689         nitems = 3;
690         iarray[ 0 ] = (*nshg);
691         iarray[ 1 ] = (*numVars);
692         iarray[ 2 ] = (*stepno);
693 
694 //MR CHANGE
695 //  Measure the time - Start the timer
696 //        MPI_Barrier(MPI_COMM_WORLD);
697 //        timer_start = rdtsc();
698 //MR CHANGE END
699 
700         writeheader( &f_descriptor, fieldtag, (void*)iarray, &nitems, &isize, "double", phasta_iotype);
701 
702 //MR CHANGE
703 //  Measure the time - End of timer
704 //        MPI_Barrier(MPI_COMM_WORLD);
705 //        timer_end = rdtsc();
706 //        time_span=(double)((timer_end-timer_start)/cpu_speed);
707 //        if (*pid==0) {
708 //          printf("Time: header for 'error':    %f s\n",time_span);
709 //        }
710 //MR CHANGE END
711 
712 //MR CHANGE
713 //  Measure the time - Start the timer
714 //        MPI_Barrier(MPI_COMM_WORLD);
715 //        timer_start = rdtsc();
716 //MR CHANGE END
717 
718         writedatablock( &f_descriptor, fieldtag, (void*)array1, &isize, "double", phasta_iotype );
719 
720 //MR CHANGE
721 //  Measure the time - End of timer
722 //        MPI_Barrier(MPI_COMM_WORLD);
723 //        timer_end = rdtsc();
724 //        time_span=(double)((timer_end-timer_start)/cpu_speed);
725 
726 //        int isizemin,isizemax,isizetot;
727 //        double sizemin,sizemax,sizeavg,sizetot,rate;
728 
729 //        MPI_Allreduce(&isize,&isizemin,1,MPI_INT,MPI_MIN,MPI_COMM_WORLD);
730 //        MPI_Allreduce(&isize,&isizemax,1,MPI_INT,MPI_MAX,MPI_COMM_WORLD);
731 //        MPI_Allreduce(&isize,&isizetot,1,MPI_INT,MPI_SUM,MPI_COMM_WORLD);
732 
733 //        sizemin=(double)(8.0*isizemin/1024.0/1024.0);
734 //        sizemax=(double)(8.0*isizemax/1024.0/1024.0);
735 //        sizetot=(double)(8.0*isizetot/1024.0/1024.0);
736 //        sizeavg=sizetot/workfc.numpe;
737 //        rate=sizetot/time_span;
738 
739 //        if (*pid==0) {
740 //          printf("Time: block for 'error':    %f s\n",time_span);
741 //          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);
742 //          printf("*****************************\n");
743 //        }
744 //MR CHANGE END
745 
746     }
747 
748 //     MPI_Barrier(MPI_COMM_WORLD);
749 //     timer_end = rdtsc();
750 //     time_span=(double)(timer_end-timer_start)/cpu_speed;
751 
752 //     if (*pid==0) {
753 //         printf("Field 'error' written in:     %f s\n",time_span);
754 //         printf("Write field '%s' finished! \n",fieldtag);
755 //     }
756 
757     if (field_flag==nfields){
758 
759 //MR CHANGE
760 //  Measure the time - Start the timer
761 //      MPI_Barrier(MPI_COMM_WORLD);
762 //      timer_start = rdtsc();
763 //MR CHANGE END
764 
765 //MR CHANGE
766 //    Measure the time - End of timer
767 //      MPI_Barrier(MPI_COMM_WORLD);
768 //      timer_end = rdtsc();
769 //      time_span=(double)((timer_end-timer_start)/cpu_speed);
770 //      if (*pid==0) {
771 //        printf("\n*****************************\n");
772 //        printf("Time: 'closefile' is:    %f s\n",time_span);
773 //      }
774 //MR CHANGE END
775 
776 //MR CHANGE
777 //  Measure the time - Start the timer
778 //      MPI_Barrier(MPI_COMM_WORLD);
779 //      timer_start = rdtsc();
780 //MR CHANGE END
781 
782       phio_closefile_write(&f_descriptor);
783 
784 //MR CHANGE
785 //    Measure the time - End of timer
786 //      MPI_Barrier(MPI_COMM_WORLD);
787 //      timer_end = rdtsc();
788 //      time_span=(double)((timer_end-timer_start)/cpu_speed);
789       if (*pid==0) {
790 //        printf("Time: 'finalizephmpiio' is:    %f s\n",time_span);
791         printf("Last field %d '%s' finished! \n",nfields, fieldtag);
792         printf("\n");
793 //        printf("*****************************\n");
794       }
795     }
796 //MR CHANGE END
797 
798     ///////////////////////////////////////////////////////////////////////////////////////////
799 
800 
801 }
802 
803 
804 void
805 Write_Displ(  int* pid,
806               int* stepno,
807               int* nshg,
808               int* numVars,
809               double* array1 ) {
810 
811 
812     char fname[255];
813     char rfile[60];
814     int irstou;
815     int magic_number = 362436;
816     int* mptr = &magic_number;
817     time_t timenow = time ( &timenow);
818     double version=0.0;
819     int isize, nitems;
820     int iarray[10];
821 
822     sprintf(rfile,"restart.%d.%d",*stepno,*pid+1);
823     openfile(rfile,"append", &irstou);
824 
825     isize = (*nshg)*(*numVars);
826     nitems = 3;
827     iarray[ 0 ] = (*nshg);
828     iarray[ 1 ] = (*numVars);
829     iarray[ 2 ] = (*stepno);
830     writeheader( &irstou, "displacement", (void*)iarray, &nitems, &isize, "double", phasta_iotype );
831 
832     nitems = (*nshg)*(*numVars);
833     writedatablock( &irstou, "displacement", (void*)(array1), &nitems, "double", phasta_iotype );
834 
835     closefile( &irstou, "append" );
836 }
837 
838 void
839 Write_Field(  int *pid,
840               char* filemode,
841               char* fieldtag,
842               int* tagsize,
843               void* array,
844               char* arraytype,
845               int* nshg,
846               int* numvars,
847               int* stepno) {
848 
849     //printf("Rank is %d, field is %s, tagsize is %d, nshg is %d, numvars is %d\n",*pid,fieldtag,*tagsize,*nshg,*numvars);
850 
851 //     char rfile[32];
852     // assuming restart.sn.(pid+1)
853 //     sprintf(rfile,"restart.%d.%d",*stepno,*pid+1);
854 
855     char *fieldlabel = (char *)malloc((*tagsize+1)*sizeof(char));
856     strncpy(fieldlabel, fieldtag, *tagsize);
857     fieldlabel[*tagsize] = '\0';
858 
859     int irstou;
860     int magic_number = 362436;
861     int* mptr = &magic_number;
862     double version=0.0;
863     int isize, nitems;
864     int iarray[10];
865 
866     char fmode[10];
867     if(!strncmp(filemode,"w",1))
868       strcpy(fmode,"write");
869     else // default is append
870       strcpy(fmode,"append");
871 
872     char datatype[10];
873     if(!strncmp(arraytype,"i",1))
874       strcpy(datatype,"int");
875     else // default is double
876       strcpy(datatype,"double");
877 
878 /*     openfile_(rfile, fmode, &irstou);
879 
880      nitems = 3; // assuming field will write 3 items in iarray
881      iarray[ 0 ] = (*nshg);
882      iarray[ 1 ] = (*numvars);
883      iarray[ 2 ] = (*stepno);
884 
885      isize = (*nshg)*(*numvars);
886      writeheader_( &irstou, fieldlabel, (void*)iarray, &nitems, &isize, datatype, phasta_iotype );
887 
888      nitems = (*nshg)*(*numvars);
889      writedatablock_( &irstou, fieldlabel, array, &nitems, datatype, phasta_iotype );
890      closefile_( &irstou, fmode);
891 */
892     /////////////////////////////// Start of writing using new-lib ////////////////////////////
893 
894     int nfiles;
895     int nfields;
896     int numparts;
897     int irank;
898     int nprocs;
899 
900 //    unsigned long long timer_start;
901 //    unsigned long long timer_end;
902 //    double time_span;
903 
904     nfiles = outpar.nsynciofiles;
905     nfields = outpar.nsynciofieldswriterestart;
906     numparts = workfc.numpe;
907     irank = *pid; //workfc.myrank;
908     nprocs = workfc.numpe;
909 
910     int nppf = numparts/nfiles;
911     int GPID;
912 
913     // Calculate number of parts each  proc deal with and where it start and end ...
914     int nppp = numparts/nprocs;// nppp : Number of parts per proc ...
915     int startpart = irank * nppp +1;// Part id from which I (myrank) start ...
916     int endpart = startpart + nppp - 1;// Part id to which I (myrank) end ...
917 
918     char filename[255],path[255],fieldtag_s[255];
919     bzero((void*)filename,255);
920     bzero((void*)fieldtag_s,255);
921 
922     strncpy(fieldlabel, fieldtag, *tagsize);
923 
924     field_flag++;
925     if(*pid==0) {
926 //      printf("\n*****************************\n");
927       printf("\n");
928       printf("The %d/%d th field to be written is '%s'\n",field_flag,nfields,fieldlabel);
929     }
930 
931     sprintf(filename,"restart-dat.%d.%d",*stepno,((int)(irank/(nprocs/nfiles))+1));
932 
933 //     MPI_Barrier(MPI_COMM_WORLD);
934 //     timer_start = rdtsc();
935 
936     int i;
937     for ( i = 0; i < nppp; i++  ) {
938         GPID = startpart + i;
939 
940         // Write solution field ...
941         sprintf(fieldtag_s,"%s@%d",fieldlabel,GPID);
942 
943         isize = (*nshg)*(*numvars);
944         nitems = 3;
945         iarray[ 0 ] = (*nshg);
946         iarray[ 1 ] = (*numvars);
947         iarray[ 2 ] = (*stepno);
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         writeheader( &f_descriptor, fieldtag_s, (void*)iarray, &nitems, &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 //        if (*pid==0) {
963 //          printf("Time: header for '%s':    %f s\n",fieldtag_s,time_span);
964 //        }
965 //MR CHANGE END
966 
967         nitems = (*nshg)*(*numvars);
968 
969 //MR CHANGE
970 //  Measure the time - Start the timer
971 //        MPI_Barrier(MPI_COMM_WORLD);
972 //        timer_start = rdtsc();
973 //MR CHANGE END
974 
975         writedatablock( &f_descriptor, fieldtag_s, array, &isize, datatype, phasta_iotype );
976 
977 //MR CHANGE
978 //  Measure the time - End of timer
979 //        MPI_Barrier(MPI_COMM_WORLD);
980 //        timer_end = rdtsc();
981 //        time_span=(double)((timer_end-timer_start)/cpu_speed);
982 
983 //        int isizemin,isizemax,isizetot;
984 //        double sizemin,sizemax,sizeavg,sizetot,rate;
985 
986 //        MPI_Allreduce(&isize,&isizemin,1,MPI_INT,MPI_MIN,MPI_COMM_WORLD);
987 //        MPI_Allreduce(&isize,&isizemax,1,MPI_INT,MPI_MAX,MPI_COMM_WORLD);
988 //        MPI_Allreduce(&isize,&isizetot,1,MPI_INT,MPI_SUM,MPI_COMM_WORLD);
989 
990 //        sizemin=(double)(8.0*isizemin/1024.0/1024.0);
991 //        sizemax=(double)(8.0*isizemax/1024.0/1024.0);
992 //        sizetot=(double)(8.0*isizetot/1024.0/1024.0);
993 //        sizeavg=sizetot/workfc.numpe;
994 //        rate=sizetot/time_span;
995 
996 //        if (*pid==0) {
997 //          printf("Time: block for '%s':    %f s\n",fieldtag_s,time_span);
998 //          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);
999 //          printf("*****************************\n");
1000 //        }
1001 //MR CHANGE END
1002 
1003     }
1004 
1005 //     MPI_Barrier(MPI_COMM_WORLD);
1006 //     timer_end = rdtsc();
1007 //     time_span=(double)(timer_end-timer_start)/cpu_speed;
1008 
1009 //     if (*pid==0) {
1010 //         printf("Field '%s' written in:     %f s\n",fieldtag,time_span);
1011 //         printf("Write field '%s' finished! \n",fieldtag_s);
1012 //     }
1013 
1014 //     if (field_flag==nfields){
1015 //       closefile_(&f_descriptor, "write");
1016 //       finalizephmpiio_(&f_descriptor);
1017 //       if(*pid==0) {
1018 //         printf("Last field %d '%s' finished! \n",nfields, fieldtag_s);
1019 //         printf("\n*****************************\n");
1020 //       }
1021 //     }
1022 
1023     if (field_flag==nfields){
1024 
1025 //MR CHANGE
1026 //  Measure the time - Start the timer
1027 //      MPI_Barrier(MPI_COMM_WORLD);
1028 //      timer_start = rdtsc();
1029 //MR CHANGE END
1030 
1031 //MR CHANGE
1032 //    Measure the time - End of timer
1033 //      MPI_Barrier(MPI_COMM_WORLD);
1034 //      timer_end = rdtsc();
1035 //      time_span=(double)((timer_end-timer_start)/cpu_speed);
1036 //      if (*pid==0) {
1037 //        printf("\n*****************************\n");
1038 //        printf("Time: 'closefile' is:    %f s\n",time_span);
1039 //      }
1040 //MR CHANGE END
1041 
1042 //MR CHANGE
1043 //  Measure the time - Start the timer
1044 //      MPI_Barrier(MPI_COMM_WORLD);
1045 //      timer_start = rdtsc();
1046 //MR CHANGE END
1047 
1048       phio_closefile_write(&f_descriptor);
1049 
1050 //MR CHANGE
1051 //    Measure the time - End of timer
1052 //      MPI_Barrier(MPI_COMM_WORLD);
1053 //      timer_end = rdtsc();
1054 //      time_span=(double)((timer_end-timer_start)/cpu_speed);
1055       if (*pid==0) {
1056 //        printf("Time: 'finalizephmpiio' is:    %f s\n",time_span);
1057         printf("Last field %d '%s' finished! \n",nfields, fieldtag);
1058         printf("\n");
1059 //        printf("*****************************\n");
1060       }
1061     }
1062 //MR CHANGE END
1063 
1064     ///////////////////////////////////////////////////////////////////////////////////////////
1065 
1066     free(fieldlabel);
1067 }
1068 
1069 //MR CHANGE
1070 void
1071 Write_PhAvg(  int* pid,
1072               char* filemode,
1073               char* fieldtag,
1074               int* tagsize,
1075               int* iphase,
1076               void* array,
1077               char* arraytype,
1078               int* nshg,
1079               int* numvars,
1080               int* stepno) {
1081 
1082     char rfile[32];
1083     // assuming restart_phase_avg_<sn>.<iphase>.<pid+1>
1084     sprintf(rfile,"restart_phase_avg_%d.%d.%d",*stepno,*iphase,*pid+1);
1085 
1086     char *fieldlabel = (char *)malloc((*tagsize+1)*sizeof(char));
1087     strncpy(fieldlabel, fieldtag, *tagsize);
1088     fieldlabel[*tagsize] = '\0';
1089 
1090     int irstou;
1091     int isize, nitems;
1092     int iarray[10];
1093 
1094     char fmode[10];
1095     if(!strncmp(filemode,"w",1))
1096       strcpy(fmode,"write");
1097     else // default is append
1098       strcpy(fmode,"append");
1099 
1100     char datatype[10];
1101     if(!strncmp(arraytype,"i",1))
1102       strcpy(datatype,"int");
1103     else // default is double
1104       strcpy(datatype,"double");
1105 
1106     openfile(rfile, fmode, &irstou);
1107 
1108     if(!strcmp(fmode,"write")) {
1109       // may be create a routine for 'top' portion under write mode
1110       int magic_number = 362436;
1111       int* mptr = &magic_number;
1112       time_t timenow = time ( &timenow);
1113       double version=0.0;
1114 
1115       /* writing the top ascii header for the restart file */
1116 
1117       writestring( &irstou,"# PHASTA Input File Version 2.0\n");
1118       writestring( &irstou,
1119                     "# format \"keyphrase : sizeofnextblock usual headers\"\n");
1120 
1121       char fname[255];
1122       bzero( (void*)fname, 255 );
1123       sprintf(fname,"# Output generated by phasta version (NOT YET CURRENT): %lf \n", version);
1124       writestring( &irstou, fname );
1125 
1126       bzero( (void*)fname, 255 );
1127       gethostname(fname,255);
1128       writestring( &irstou,"# This result was produced on: ");
1129       writestring( &irstou, fname );
1130       writestring( &irstou,"\n");
1131 
1132       bzero( (void*)fname, 255 );
1133       sprintf(fname,"# %s\n", ctime( &timenow ));
1134       writestring( &irstou, fname );
1135 
1136       isize = 1;
1137       nitems = 1;
1138       iarray[ 0 ] = 1;
1139       writeheader( &irstou, "byteorder magic number ",
1140                     (void*)iarray, &nitems, &isize, "integer", phasta_iotype );
1141       nitems = 1;
1142       writedatablock( &irstou, "byteorder magic number ",
1143                        (void*)mptr, &nitems, "integer", phasta_iotype );
1144     }
1145 
1146     nitems = 3; // assuming field will write 3 items in iarray
1147     iarray[ 0 ] = (*nshg);
1148     iarray[ 1 ] = (*numvars);
1149     iarray[ 2 ] = (*stepno);
1150 
1151     isize = (*nshg)*(*numvars);
1152     writeheader( &irstou, fieldlabel, (void*)iarray, &nitems, &isize, datatype, phasta_iotype );
1153 
1154     nitems = (*nshg)*(*numvars);
1155     writedatablock( &irstou, fieldlabel, array, &nitems, datatype, phasta_iotype );
1156 
1157     closefile( &irstou, fmode);
1158 
1159     free(fieldlabel);
1160 }
1161 
1162 //MR CHANGE
1163 void
1164 Write_PhAvg2( int* pid,
1165               char* filemode,
1166               char* fieldtag,
1167               int* tagsize,
1168               int* iphase,
1169               int* nphasesincycle,
1170               void* array,
1171               char* arraytype,
1172               int* nshg,
1173               int* numvars,
1174               int* stepno) {
1175 
1176 //     char rfile[32];
1177     // assuming restart.sn.(pid+1)
1178 //     sprintf(rfile,"restart.%d.%d",*stepno,*pid+1);
1179 
1180     int addtagsize=0; // phase number is added to the name of the field
1181     if(*iphase<10)
1182       addtagsize=1;
1183     else if(*iphase<100)
1184       addtagsize=2;
1185     else if(*iphase<1000)
1186       addtagsize=3;
1187 
1188     int tagsize2;
1189     tagsize2=*tagsize+addtagsize;
1190 
1191 //     char *fieldlabel = (char *)malloc((*tagsize+1)*sizeof(char));
1192 //     strncpy(fieldlabel, fieldtag, *tagsize);
1193 //     fieldlabel[*tagsize] = '\0';
1194 
1195     char *fieldlabel = (char *)malloc((tagsize2+1)*sizeof(char));
1196     strncpy(fieldlabel, fieldtag, *tagsize);
1197     fieldlabel[tagsize2] = '\0';
1198 
1199     char straddtagsize[10];
1200     sprintf(straddtagsize,"%d",*iphase);
1201 
1202     if(*iphase<10) {
1203       fieldlabel[tagsize2-1]=straddtagsize[0];
1204     }
1205     else if(*iphase<100) {
1206       fieldlabel[tagsize2-2]=straddtagsize[0];
1207       fieldlabel[tagsize2-1]=straddtagsize[1];
1208     }
1209     else if(*iphase<1000) {
1210       fieldlabel[tagsize2-3]=straddtagsize[0];
1211       fieldlabel[tagsize2-2]=straddtagsize[1];
1212       fieldlabel[tagsize2-1]=straddtagsize[2];
1213     }
1214 
1215     int irstou;
1216     int magic_number = 362436;
1217     int* mptr = &magic_number;
1218     double version=0.0;
1219     int isize, nitems;
1220     int iarray[10];
1221 
1222     char fmode[10];
1223     if(!strncmp(filemode,"w",1))
1224       strcpy(fmode,"write");
1225     else // default is append
1226       strcpy(fmode,"append");
1227 
1228     char datatype[10];
1229     if(!strncmp(arraytype,"i",1))
1230       strcpy(datatype,"int");
1231     else // default is double
1232       strcpy(datatype,"double");
1233 
1234 //
1235 // //     if(*iphase==1) //open the file but then keep it open for the remaining cycles
1236 //     openfile_(rfile, fmode, &irstou);
1237 //
1238 // //     printf("iphase: %d - pid: %d - irstou %d\n",*iphase,*pid,irstou);
1239 //
1240 //
1241 //     nitems = 3; // assuming field will write 3 items in iarray
1242 //     iarray[ 0 ] = (*nshg);
1243 //     iarray[ 1 ] = (*numvars);
1244 //     iarray[ 2 ] = (*stepno);
1245 //
1246 //     isize = (*nshg)*(*numvars);
1247 //     writeheader_( &irstou, fieldlabel, (void*)iarray, &nitems, &isize, datatype, phasta_iotype );
1248 //
1249 //     nitems = (*nshg)*(*numvars);
1250 //     writedatablock_( &irstou, fieldlabel, array, &nitems, datatype, phasta_iotype );
1251 //
1252 // //     if(*iphase==*nphasesincycle) //close the file after nphasesincycle
1253 //       closefile_( &irstou, fmode);
1254 //
1255 
1256     /////////////////////////////// Start of writing using new-lib ////////////////////////////
1257 
1258     int nfiles;
1259     int nfields;
1260     int numparts;
1261     int irank;
1262     int nprocs;
1263 //    unsigned long long timer_start;
1264 //    unsigned long long timer_end;
1265 //    double time_span;
1266 
1267     nfiles = outpar.nsynciofiles;
1268     nfields = outpar.nsynciofieldswriterestart;
1269     numparts = workfc.numpe;
1270     irank = *pid; //workfc.myrank;
1271     nprocs = workfc.numpe;
1272 
1273     int nppf = numparts/nfiles;
1274     int GPID;
1275 
1276     // Calculate number of parts each  proc deal with and where it start and end ...
1277     int nppp = numparts/nprocs;// nppp : Number of parts per proc ...
1278     int startpart = irank * nppp +1;// Part id from which I (myrank) start ...
1279     int endpart = startpart + nppp - 1;// Part id to which I (myrank) end ...
1280 
1281     //int descriptor;
1282     char filename[255],path[255],fieldtag_s[255];
1283     bzero((void*)filename,255);
1284     bzero((void*)fieldtag_s,255);
1285 
1286 //     char * namer;
1287 //     namer = strtok(fieldlabel," ");
1288 //     strncpy(fieldlabel, fieldtag, *tagsize);
1289 
1290     field_flag++;
1291     if(*pid==0) {
1292 //      printf("\n*****************************\n");
1293       printf("\n");
1294       printf("The %d/%d th field to be written is '%s'\n",field_flag,nfields,fieldlabel);
1295     }
1296 
1297     sprintf(filename,"restart-dat.%d.%d",*stepno,((int)(irank/(nprocs/nfiles))+1));
1298 
1299     int i;
1300     for ( i = 0; i < nppp; i++  ) {
1301         GPID = startpart + i;
1302 
1303         // Write solution field ...
1304         sprintf(fieldtag_s,"%s@%d",fieldlabel,GPID);
1305 
1306         //printf("This is %d and fieldtag_s is %s \n",myrank,fieldtag_s);
1307 
1308         isize = (*nshg)*(*numvars);
1309         nitems = 3;
1310         iarray[ 0 ] = (*nshg);
1311         iarray[ 1 ] = (*numvars);
1312         iarray[ 2 ] = (*stepno);
1313 
1314 //MR CHANGE
1315 //  Measure the time - Start the timer
1316 //        MPI_Barrier(MPI_COMM_WORLD);
1317 //        timer_start = rdtsc();
1318 //MR CHANGE END
1319 
1320         writeheader( &f_descriptor, fieldtag_s, (void*)iarray, &nitems, &isize, "double", phasta_iotype);
1321 
1322 //MR CHANGE
1323 //  Measure the time - End of timer
1324 //        MPI_Barrier(MPI_COMM_WORLD);
1325 //        timer_end = rdtsc();
1326 //        time_span=(double)((timer_end-timer_start)/cpu_speed);
1327 //        if (*pid==0) {
1328 //          printf("Time: header for '%s':    %f s\n",fieldtag_s,time_span);
1329 //        }
1330 //MR CHANGE END
1331 
1332         nitems = (*nshg)*(*numvars);
1333 
1334 //MR CHANGE
1335 //  Measure the time - Start the timer
1336 //        MPI_Barrier(MPI_COMM_WORLD);
1337 //        timer_start = rdtsc();
1338 //MR CHANGE END
1339 
1340         writedatablock( &f_descriptor, fieldtag_s, array, &isize, "double", phasta_iotype );
1341 
1342 //MR CHANGE
1343 //  Measure the time - End of timer
1344 //        MPI_Barrier(MPI_COMM_WORLD);
1345 //        timer_end = rdtsc();
1346 //        time_span=(double)((timer_end-timer_start)/cpu_speed);
1347 
1348 //        int isizemin,isizemax,isizetot;
1349 //        double sizemin,sizemax,sizeavg,sizetot,rate;
1350 
1351 //        MPI_Allreduce(&isize,&isizemin,1,MPI_INT,MPI_MIN,MPI_COMM_WORLD);
1352 //        MPI_Allreduce(&isize,&isizemax,1,MPI_INT,MPI_MAX,MPI_COMM_WORLD);
1353 //        MPI_Allreduce(&isize,&isizetot,1,MPI_INT,MPI_SUM,MPI_COMM_WORLD);
1354 
1355 //        sizemin=(double)(8.0*isizemin/1024.0/1024.0);
1356 //        sizemax=(double)(8.0*isizemax/1024.0/1024.0);
1357 //        sizetot=(double)(8.0*isizetot/1024.0/1024.0);
1358 //        sizeavg=sizetot/workfc.numpe;
1359 //        rate=sizetot/time_span;
1360 
1361 //        if (*pid==0) {
1362 //          printf("Time: block for '%s':    %f s\n",fieldtag_s,time_span);
1363 //          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);
1364 //          printf("*****************************\n");
1365 //        }
1366 //MR CHANGE END
1367 
1368     }
1369 
1370 //     if (*pid==0) {
1371 //         printf("Field '%s' written in:     %f s\n",fieldtag,time_span);
1372 //         printf("Write field '%s' finished! \n",fieldtag_s);
1373 //     }
1374 
1375 //
1376 //     if (field_flag==nfields){
1377 //       closefile_(&f_descriptor, "write");
1378 //       finalizephmpiio_(&f_descriptor);
1379 //       if(*pid==0) {
1380 //         printf("Last field %d '%s' finished! \n",nfields, fieldtag_s);
1381 //         printf("\n*****************************\n");
1382 //       }
1383 //     }
1384 
1385     if (field_flag==nfields){
1386 
1387 //MR CHANGE
1388 //  Measure the time - Start the timer
1389 //      MPI_Barrier(MPI_COMM_WORLD);
1390 //      timer_start = rdtsc();
1391 //MR CHANGE END
1392 
1393 //MR CHANGE
1394 //    Measure the time - End of timer
1395 //      MPI_Barrier(MPI_COMM_WORLD);
1396 //      timer_end = rdtsc();
1397 //      time_span=(double)((timer_end-timer_start)/cpu_speed);
1398 //     if (*pid==0) {
1399 //        printf("\n*****************************\n");
1400 //        printf("Time: 'closefile' is:    %f s\n",time_span);
1401 //      }
1402 //MR CHANGE END
1403 
1404 //MR CHANGE
1405 //  Measure the time - Start the timer
1406 //      MPI_Barrier(MPI_COMM_WORLD);
1407 //      timer_start = rdtsc();
1408 //MR CHANGE END
1409 
1410       phio_closefile_write(&f_descriptor);
1411 
1412 //MR CHANGE
1413 //    Measure the time - End of timer
1414 //      MPI_Barrier(MPI_COMM_WORLD);
1415 //      timer_end = rdtsc();
1416 //      time_span=(double)((timer_end-timer_start)/cpu_speed);
1417       if (*pid==0) {
1418 //        printf("Time: 'finalizephmpiio' is:    %f s\n",time_span);
1419 //        printf("Last field %d '%s' finished! \n",nfields, fieldtag);
1420         printf("\n");
1421 //        printf("*****************************\n");
1422       }
1423     }
1424 //MR CHANGE END
1425 
1426     ///////////////////////////////////////////////////////////////////////////////////////////
1427 
1428     free(fieldlabel);
1429 }
1430 
1431 
1432 void
1433 Write_d2wall(   int* pid,
1434                 int* numnp,
1435                 double* array1 ) {
1436 
1437 //    time_t timenow = time ( &timenow);
1438     int isize, nitems;
1439     int iarray[10];
1440 
1441 //    MPI_Barrier(MPI_COMM_WORLD);
1442 
1443     /////////////////////////////// Start of writing using new-lib ////////////////////////////
1444 
1445     int nfiles;
1446     int nfields;
1447     int numparts;
1448     int irank;
1449     int nprocs;
1450 
1451     //  First, count the number of fields to write and store the result in
1452     //countfieldstowriterestart();
1453 
1454     //  Retrieve and compute the parameters required for SyncIO
1455     nfiles = outpar.nsynciofiles;
1456     nfields = 1; //outpar.nsynciofieldswriterestart;  // Only the distance to the walls in d2wall
1457     numparts = workfc.numpe;
1458     irank = *pid; //workfc.myrank;
1459     nprocs = workfc.numpe;
1460     int nppf = numparts/nfiles;
1461     int GPID;
1462 
1463     // Calculate number of parts each proc deal with and where it start and end ...
1464     int nppp = numparts/nprocs;// nppp : Number of parts per proc ...
1465     int startpart = irank * nppp +1;// Part id from which I (myrank) start ...
1466     int endpart = startpart + nppp - 1;// Part id to which I (myrank) end ...
1467 
1468     int descriptor;
1469     char filename[255],path[255],fieldtag_s[255];
1470     bzero((void*)filename,255);
1471     bzero((void*)fieldtag_s,255);
1472 
1473     phio_openfile_write("d2wall.", &nfiles, &nfields, &nppf, &f_descriptor);
1474 
1475     field_flag=0;
1476 
1477      int i;
1478      for ( i = 0; i < nppp; i++) { //This loop is useful only if several parts per processor
1479      // GPID : global part id, corresponds to rank ...
1480         // e.g : (in this example)
1481         // proc 0 : 1--4
1482         // proc 1 : 5--8 ...
1483         GPID = startpart + i;
1484 
1485         // Write solution field ...
1486         sprintf(fieldtag_s,"d2wall@%d",GPID);
1487 
1488         isize = (*numnp);
1489         nitems = 2;
1490         iarray[ 0 ] = (*numnp);
1491         iarray[ 1 ] = 1; //numVars = 1
1492 
1493         writeheader( &f_descriptor, fieldtag_s, (void*)iarray, &nitems, &isize, "double", phasta_iotype);
1494 
1495         //nitems = (*nshg)*(*numVars);
1496         //nitems = (*numnp);
1497 
1498         writedatablock( &f_descriptor, fieldtag_s, (void*)(array1), &isize, "double", phasta_iotype );
1499 
1500 
1501     }
1502     field_flag++;
1503 
1504     if (field_flag==nfields){
1505       phio_closefile_write(&f_descriptor);
1506       if (irank==0) {
1507         printf("\n");
1508       }
1509     }
1510 }
1511 
1512