xref: /phasta/phSolver/common/new_interface.c (revision 4ead34e5d5e094065da4d96ce571099cb42db973)
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 #include "common_c.h"
19 
20 #ifdef intel
21 #include <winsock2.h>
22 #else
23 #include <unistd.h>
24 #include <strings.h>
25 #endif
26 
27 void igetMinMaxAvg(int *ivalue, double *stats, int *statRanks) {
28   int isThisRank;
29   double *value = (double*)malloc(sizeof(double));
30   *value = 1.0*(*ivalue);
31   rgetMinMaxAvg(value,stats,statRanks);
32   free(value);
33 }
34 
35 void rgetMinMaxAvg(double *value, double *stats, int *statRanks) {
36   int isThisRank;
37 
38   MPI_Allreduce(value,&stats[0],1,MPI_DOUBLE,MPI_MIN,MPI_COMM_WORLD);
39   isThisRank=workfc.numpe+1;
40   if(*value==stats[0])
41     isThisRank=workfc.myrank;
42   MPI_Allreduce(&isThisRank,&statRanks[0],1,MPI_INT,MPI_MIN,MPI_COMM_WORLD);
43 
44   MPI_Allreduce(value,&stats[1],1,MPI_DOUBLE,MPI_MAX,MPI_COMM_WORLD);
45   isThisRank=workfc.numpe+1;
46   if(*value==stats[1])
47     isThisRank=workfc.myrank;
48   MPI_Allreduce(&isThisRank,&statRanks[1],1,MPI_INT,MPI_MIN,MPI_COMM_WORLD);
49 
50   MPI_Allreduce(value,&stats[2],1,MPI_DOUBLE,MPI_SUM,MPI_COMM_WORLD);
51   stats[2] /= workfc.numpe;
52 
53   double sqValue = (*value)*(*value), sqValueAvg = 0.;
54   MPI_Allreduce(&sqValue,&sqValueAvg,1,MPI_DOUBLE,MPI_SUM,MPI_COMM_WORLD);
55   sqValueAvg /= workfc.numpe;
56 
57   stats[3] = sqrt(sqValueAvg-stats[2]*stats[2]);
58 }
59 
60 void print_mesh_stats(void) {
61   int statRanks[2];
62   double iStats[4], rStats[4];
63 
64   igetMinMaxAvg(&conpar.nshg,iStats,statRanks);
65   if(workfc.myrank==workfc.master)
66     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]);
67   igetMinMaxAvg(&conpar.numel,iStats,statRanks);
68   if(workfc.myrank==workfc.master)
69     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]);
70   igetMinMaxAvg(&conpar.numelb,iStats,statRanks);
71   if(workfc.myrank==workfc.master)
72     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]);
73   igetMinMaxAvg(&conpar.nnz_tot,iStats,statRanks);
74   if(workfc.myrank==workfc.master) {
75     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]);
76     printf("\n");
77   }
78 }
79 
80 void print_mpi_stats(void) {
81   int statRanks[2];
82   double iStats[4], rStats[4];
83 
84 // NS equations
85   igetMinMaxAvg(&mpistats.iISend,iStats,statRanks);
86   if(workfc.myrank==workfc.master)
87     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]);
88   igetMinMaxAvg(&mpistats.iIRecv,iStats,statRanks);
89   if(workfc.myrank==workfc.master)
90     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]);
91   igetMinMaxAvg(&mpistats.iWaitAll,iStats,statRanks);
92   if(workfc.myrank==workfc.master)
93     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]);
94   igetMinMaxAvg(&mpistats.iAllR,iStats,statRanks);
95   if(workfc.myrank==workfc.master)
96     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]);
97 
98   rgetMinMaxAvg(&mpistats.rISend,rStats,statRanks);
99   if(workfc.myrank==workfc.master)
100     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]);
101   rgetMinMaxAvg(&mpistats.rIRecv,rStats,statRanks);
102   if(workfc.myrank==workfc.master)
103     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]);
104   rgetMinMaxAvg(&mpistats.rWaitAll,rStats,statRanks);
105   if(workfc.myrank==workfc.master)
106     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]);
107   rgetMinMaxAvg(&mpistats.rCommu,rStats,statRanks);
108   if(workfc.myrank==workfc.master)
109     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]);
110   rgetMinMaxAvg(&mpistats.rAllR,rStats,statRanks);
111   if(workfc.myrank==workfc.master) {
112     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]);
113     printf("\n");
114   }
115 // Scalars
116   igetMinMaxAvg(&mpistats.iISendScal,iStats,statRanks);
117   if(workfc.myrank==workfc.master)
118     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]);
119   igetMinMaxAvg(&mpistats.iIRecvScal,iStats,statRanks);
120   if(workfc.myrank==workfc.master)
121     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]);
122   igetMinMaxAvg(&mpistats.iWaitAllScal,iStats,statRanks);
123   if(workfc.myrank==workfc.master)
124     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]);
125   igetMinMaxAvg(&mpistats.iAllRScal,iStats,statRanks);
126   if(workfc.myrank==workfc.master)
127     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]);
128 
129   rgetMinMaxAvg(&mpistats.rISendScal,rStats,statRanks);
130   if(workfc.myrank==workfc.master)
131     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]);
132   rgetMinMaxAvg(&mpistats.rIRecvScal,rStats,statRanks);
133   if(workfc.myrank==workfc.master)
134     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]);
135   rgetMinMaxAvg(&mpistats.rWaitAllScal,rStats,statRanks);
136   if(workfc.myrank==workfc.master)
137     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]);
138   rgetMinMaxAvg(&mpistats.rCommuScal,rStats,statRanks);
139   if(workfc.myrank==workfc.master)
140     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]);
141   rgetMinMaxAvg(&mpistats.rAllRScal,rStats,statRanks);
142   if(workfc.myrank==workfc.master)
143     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]);
144 
145 
146 }
147 
148 void print_system_stats(double *tcorecp, double *tcorecpscal) {
149   int statRanks[2];
150   double iStats[4], rStats[4];
151   double syst_assembly, syst_solve;
152 
153 // NS equations
154   syst_assembly = tcorecp[0];
155   syst_solve = tcorecp[1];
156 
157   rgetMinMaxAvg(&syst_assembly,rStats,statRanks);
158   if(workfc.myrank==workfc.master)
159     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]);
160 
161   rgetMinMaxAvg(&syst_solve,rStats,statRanks);
162   if(workfc.myrank==workfc.master)
163     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]);
164 
165 // Scalars
166   syst_assembly = tcorecpscal[0];
167   syst_solve = tcorecpscal[1];
168 
169   rgetMinMaxAvg(&syst_assembly,rStats,statRanks);
170   if(workfc.myrank==workfc.master)
171     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]);
172 
173   rgetMinMaxAvg(&syst_solve,rStats,statRanks);
174   if(workfc.myrank==workfc.master) {
175     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]);
176     printf("\n");
177   }
178 }
179 
180 
181 
182 void countfieldstowriterestart()
183 {
184   int nfields = 2; //solution, time derivatives
185 
186   if(outpar.ivort == 1){
187     nfields++; //vorticity
188   }
189 
190   if(abs(turbvar.itwmod) != 1 && outpar.iowflux == 1) {
191     nfields++; //instantaneous wss in bflux.f
192   }
193 
194   if(timdat.istep == inpdat.nstep[timdat.itseq-1]){ //Last time step of the computation
195 
196     //projection vectors and pressure projection vectors (call saveLesRestart in itrdrv)
197     if(matdat.matflg[0][0] ==-1) {
198         nfields = nfields +2;
199     }
200 
201     //if Print Error Indicators = true (call write_error in itrdrv)
202     if(turbvar.ierrcalc == 1){
203       nfields++;
204     }
205 
206     //if Print ybar = True (call write_field(myrank,'a','ybar',4,... in itrdrv)
207     if(outpar.ioybar == 1){
208       nfields++;  //ybar
209 
210       //phase average fields
211       if(outpar.nphasesincycle >0) {
212         nfields = nfields + outpar.nphasesincycle;
213       }
214 
215       if(abs(turbvar.itwmod) != 1 && outpar.iowflux == 1) {
216         nfields++; //wssbar
217       }
218 
219     }
220 
221     if(turbvari.irans < 0) {
222       nfields++; //dwal
223     }
224 
225   }
226 
227   outpar.nsynciofieldswriterestart = nfields;
228 
229   if(workfc.myrank == 0) {
230     printf("Number of fields to write in restart files: %d\n", nfields);
231   }
232 }
233 
234 
235 void
236 Write_Restart(  int* pid,
237                 int* stepno,
238                 int* nshg,
239                 int* numVars,
240                 double* array1,
241                 double* array2 ) {
242 
243     char fname[255];
244     char rfile[60];
245     char existingfile[30], linkfile[30];
246     int irstou;
247     int magic_number = 362436;
248     int* mptr = &magic_number;
249     double version=0.0;
250     int isize, nitems;
251     int iarray[10];
252     int nfiles;
253     int nfields;
254     int numparts;
255     int irank;
256     int nprocs;
257     int nppf;
258 
259     //  First, count the number of fields to write and store the result in
260     countfieldstowriterestart();
261 
262     //  Retrieve and compute the parameters required for SyncIO
263     nfiles = outpar.nsynciofiles;
264     nfields = outpar.nsynciofieldswriterestart;
265     numparts = workfc.numpe;
266     irank = *pid; //workfc.myrank;
267     nprocs = workfc.numpe;
268     int nppp = numparts/nprocs;   // always 1 for PHASTA
269     int descriptor;
270     char filename[255];
271     bzero((void*)filename,255);
272 
273     if(nfiles == 0 ){
274       sprintf(filename,"restart.%d.", *stepno);
275       nppf=1;
276     } else {
277       nppf=numparts/nfiles;
278       sprintf(filename,"restart-dat.%d.", *stepno);
279     }
280     phio_openfile_write(filename, &nfiles, &nfields, &nppf, &f_descriptor);
281 
282     field_flag=0;
283 
284      int i;
285      for ( i = 0; i < nppp; i++) { //This loop is useful only if several parts per processor
286         // Write solution field ...
287         isize = (*nshg)*(*numVars);
288         nitems = 3;
289         iarray[ 0 ] = (*nshg);
290         iarray[ 1 ] = (*numVars);
291         iarray[ 2 ] = (*stepno);
292 
293         phio_writeheader(f_descriptor, "solution", (void*)iarray, &nitems,
294             &isize, "double", phasta_iotype);
295         nitems = (*nshg)*(*numVars);
296         phio_writedatablock(f_descriptor, "solution", (void*)(array1),
297             &isize, "double", phasta_iotype );
298     }
299     field_flag++;
300 
301     for ( i = 0; i < nppp; i++) {
302         // Write solution field ...
303         isize = (*nshg)*(*numVars);
304         nitems = 3;
305         iarray[ 0 ] = (*nshg);
306         iarray[ 1 ] = (*numVars);
307         iarray[ 2 ] = (*stepno);
308         phio_writeheader(f_descriptor, "time derivative of solution",
309             (void*)iarray, &nitems, &isize, "double", phasta_iotype);
310         nitems = (*nshg)*(*numVars);
311         phio_writedatablock(f_descriptor, "time derivative of solution",
312             (void*)(array2), &isize, "double", phasta_iotype );
313     }
314     field_flag++;
315 
316     if (field_flag==nfields){
317       phio_closefile_write(f_descriptor);
318       if (*pid==0) {
319         printf("\n");
320       }
321     }
322 }
323 
324 void
325 Write_Error(  int* pid,
326               int* stepno,
327               int* nshg,
328               int* numVars,
329               double* array1 ) {
330     char fname[255];
331     char rfile[60];
332     int irstou;
333     int magic_number = 362436;
334     int* mptr = &magic_number;
335     double version=0.0;
336     int isize, nitems;
337     int iarray[10];
338     int nfiles;
339     int nfields;
340     int numparts;
341     int irank;
342     int nprocs;
343 
344     nfiles = outpar.nsynciofiles;
345     nfields = outpar.nsynciofieldswriterestart;
346     numparts = workfc.numpe;
347     irank = *pid; //workfc.myrank;
348     nprocs = workfc.numpe;
349 
350     // Calculate number of parts each  proc deal with and where it start and end ...
351     int nppp = numparts/nprocs;// nppp : Number of parts per proc ...
352     int startpart = irank * nppp +1;// Part id from which I (myrank) start ...
353     int endpart = startpart + nppp - 1;// Part id to which I (myrank) end ...
354 
355     field_flag++;
356 
357     int i;
358     for ( i = 0; i < nppp; i++  ) {
359 
360         if(*pid==0) {
361           printf("\n");
362           printf("The %d/%d th field to be written is 'errors'\n",field_flag,nfields);
363         }
364 
365         isize = (*nshg)*(*numVars);
366         nitems = 3;
367         iarray[ 0 ] = (*nshg);
368         iarray[ 1 ] = (*numVars);
369         iarray[ 2 ] = (*stepno);
370 
371         phio_writeheader(f_descriptor, "errors", (void*)iarray, &nitems,
372             &isize, "double", phasta_iotype);
373 
374         phio_writedatablock(f_descriptor, "errors", (void*)array1, &isize,
375             "double", phasta_iotype );
376     }
377     if (field_flag==nfields){
378       phio_closefile_write(f_descriptor);
379       if (*pid==0) {
380         printf("Last field %d 'errors' finished! \n",nfields);
381         printf("\n");
382       }
383     }
384 }
385 
386 void
387 Write_Displ(  int* pid,
388               int* stepno,
389               int* nshg,
390               int* numVars,
391               double* array1 ) {
392   fprintf(stderr, "This function is dead...exiting\n");
393   exit(1);
394 }
395 
396 void
397 Write_Field(  int *pid,
398               char* filemode,
399               char* fieldtag,
400               int* tagsize,
401               void* array,
402               char* arraytype,
403               int* nshg,
404               int* numvars,
405               int* stepno) {
406     char *fieldlabel = (char *)malloc((*tagsize+1)*sizeof(char));
407     strncpy(fieldlabel, fieldtag, *tagsize);
408     fieldlabel[*tagsize] = '\0';
409 
410     int irstou;
411     int magic_number = 362436;
412     int* mptr = &magic_number;
413     double version=0.0;
414     int isize, nitems;
415     int iarray[10];
416 
417     char fmode[10];
418     if(!strncmp(filemode,"w",1))
419       strcpy(fmode,"write");
420     else // default is append
421       strcpy(fmode,"append");
422 
423     char datatype[10];
424     if(!strncmp(arraytype,"i",1))
425       strcpy(datatype,"int");
426     else // default is double
427       strcpy(datatype,"double");
428 
429     int nfiles;
430     int nfields;
431     int numparts;
432     int irank;
433     int nprocs;
434 
435     nfiles = outpar.nsynciofiles;
436     nfields = outpar.nsynciofieldswriterestart;
437     numparts = workfc.numpe;
438     irank = *pid; //workfc.myrank;
439     nprocs = workfc.numpe;
440 
441     // Calculate number of parts each  proc deal with and where it start and end ...
442     int nppp = numparts/nprocs;// nppp : Number of parts per proc ...
443     int startpart = irank * nppp +1;// Part id from which I (myrank) start ...
444     int endpart = startpart + nppp - 1;// Part id to which I (myrank) end ...
445 
446     strncpy(fieldlabel, fieldtag, *tagsize);
447 
448     field_flag++;
449     if(*pid==0) {
450       printf("\n");
451       printf("The %d/%d th field to be written is '%s'\n",field_flag,nfields,fieldlabel);
452     }
453 
454     int i;
455     for ( i = 0; i < nppp; i++  ) {
456         // Write solution field ...
457         isize = (*nshg)*(*numvars);
458         nitems = 3;
459         iarray[ 0 ] = (*nshg);
460         iarray[ 1 ] = (*numvars);
461         iarray[ 2 ] = (*stepno);
462 
463         phio_writeheader(f_descriptor, fieldlabel, (void*)iarray, &nitems,
464             &isize, datatype, phasta_iotype);
465         nitems = (*nshg)*(*numvars);
466         phio_writedatablock(f_descriptor, fieldlabel, array, &isize,
467             datatype, phasta_iotype );
468     }
469     if (field_flag==nfields){
470       phio_closefile_write(f_descriptor);
471       if (*pid==0) {
472         printf("Last field %d '%s' finished! \n",nfields, fieldtag);
473         printf("\n");
474       }
475     }
476     free(fieldlabel);
477 }
478 
479 void
480 Write_PhAvg2( int* pid,
481               char* filemode,
482               char* fieldtag,
483               int* tagsize,
484               int* iphase,
485               int* nphasesincycle,
486               void* array,
487               char* arraytype,
488               int* nshg,
489               int* numvars,
490               int* stepno) {
491     int addtagsize=0; // phase number is added to the name of the field
492     if(*iphase<10)
493       addtagsize=1;
494     else if(*iphase<100)
495       addtagsize=2;
496     else if(*iphase<1000)
497       addtagsize=3;
498 
499     int tagsize2;
500     tagsize2=*tagsize+addtagsize;
501 
502     char *fieldlabel = (char *)malloc((tagsize2+1)*sizeof(char));
503     strncpy(fieldlabel, fieldtag, *tagsize);
504     fieldlabel[tagsize2] = '\0';
505 
506     char straddtagsize[10];
507     sprintf(straddtagsize,"%d",*iphase);
508 
509     if(*iphase<10) {
510       fieldlabel[tagsize2-1]=straddtagsize[0];
511     }
512     else if(*iphase<100) {
513       fieldlabel[tagsize2-2]=straddtagsize[0];
514       fieldlabel[tagsize2-1]=straddtagsize[1];
515     }
516     else if(*iphase<1000) {
517       fieldlabel[tagsize2-3]=straddtagsize[0];
518       fieldlabel[tagsize2-2]=straddtagsize[1];
519       fieldlabel[tagsize2-1]=straddtagsize[2];
520     }
521 
522     int irstou;
523     int magic_number = 362436;
524     int* mptr = &magic_number;
525     double version=0.0;
526     int isize, nitems;
527     int iarray[10];
528 
529     char fmode[10];
530     if(!strncmp(filemode,"w",1))
531       strcpy(fmode,"write");
532     else // default is append
533       strcpy(fmode,"append");
534 
535     char datatype[10];
536     if(!strncmp(arraytype,"i",1))
537       strcpy(datatype,"int");
538     else // default is double
539       strcpy(datatype,"double");
540 
541     int nfiles;
542     int nfields;
543     int numparts;
544     int irank;
545     int nprocs;
546 
547     nfiles = outpar.nsynciofiles;
548     nfields = outpar.nsynciofieldswriterestart;
549     numparts = workfc.numpe;
550     irank = *pid; //workfc.myrank;
551     nprocs = workfc.numpe;
552 
553     // Calculate number of parts each  proc deal with and where it start and end ...
554     int nppp = numparts/nprocs;// nppp : Number of parts per proc ...
555     int startpart = irank * nppp +1;// Part id from which I (myrank) start ...
556     int endpart = startpart + nppp - 1;// Part id to which I (myrank) end ...
557 
558     field_flag++;
559     if(*pid==0) {
560       printf("\n");
561       printf("The %d/%d th field to be written is '%s'\n",field_flag,nfields,fieldlabel);
562     }
563 
564     int i;
565     for ( i = 0; i < nppp; i++  ) {
566         // Write solution field ...
567         isize = (*nshg)*(*numvars);
568         nitems = 3;
569         iarray[ 0 ] = (*nshg);
570         iarray[ 1 ] = (*numvars);
571         iarray[ 2 ] = (*stepno);
572         phio_writeheader(f_descriptor, fieldlabel, (void*)iarray, &nitems,
573             &isize, "double", phasta_iotype);
574         nitems = (*nshg)*(*numvars);
575         phio_writedatablock(f_descriptor, fieldlabel, array, &isize,
576             "double", phasta_iotype );
577     }
578     if (field_flag==nfields){
579       phio_closefile_write(f_descriptor);
580       if (*pid==0) {
581         printf("\n");
582       }
583     }
584     free(fieldlabel);
585 }
586 
587 
588 void
589 Write_d2wall(   int* pid,
590                 int* numnp,
591                 double* array1 ) {
592     int isize, nitems;
593     int iarray[10];
594     int nfiles;
595     int nfields;
596     int numparts;
597     int irank;
598     int nprocs;
599 
600     //  Retrieve and compute the parameters required for SyncIO
601     nfiles = outpar.nsynciofiles;
602     nfields = 1; //outpar.nsynciofieldswriterestart;  // Only the distance to the walls in d2wall
603     numparts = workfc.numpe;
604     irank = *pid; //workfc.myrank;
605     nprocs = workfc.numpe;
606 
607     // Calculate number of parts each proc deal with and where it start and end ...
608     int nppp = numparts/nprocs;// nppp : Number of parts per proc ...
609     int startpart = irank * nppp +1;// Part id from which I (myrank) start ...
610     int endpart = startpart + nppp - 1;// Part id to which I (myrank) end ...
611 
612     int descriptor;
613     int nppf=numparts/nfiles;
614 
615     phio_openfile_write("d2wall.", &nfiles, &nfields, &nppf, &f_descriptor);
616 
617     field_flag=0;
618 
619      int i;
620      for ( i = 0; i < nppp; i++) { //This loop is useful only if several parts per processor
621         // Write solution field ...
622         isize = (*numnp);
623         nitems = 2;
624         iarray[ 0 ] = (*numnp);
625         iarray[ 1 ] = 1; //numVars = 1
626 
627         phio_writeheader( &f_descriptor, "d2wall", (void*)iarray, &nitems,
628             &isize, "double", phasta_iotype);
629         phio_writedatablock( &f_descriptor, "d2wall", (void*)(array1), &isize,
630             "double", phasta_iotype );
631     }
632     field_flag++;
633 
634     if (field_flag==nfields){
635       phio_closefile_write(&f_descriptor);
636       if (irank==0) {
637         printf("\n");
638       }
639     }
640 }
641