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