xref: /phasta/phSolver/common/new_interface.c (revision dc9538425e8e2132ad0d36a4b7a0472f4a4e7110)
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 
285      int i;
286      for ( i = 0; i < nppp; i++) { //This loop is useful only if several parts per processor
287         // Write solution field ...
288         isize = (*nshg)*(*numVars);
289         nitems = 3;
290         iarray[ 0 ] = (*nshg);
291         iarray[ 1 ] = (*numVars);
292         iarray[ 2 ] = (*stepno);
293 
294         phio_writeheader(f_descriptor, "solution", (void*)iarray, &nitems,
295             &isize, "double", phasta_iotype);
296         nitems = (*nshg)*(*numVars);
297         phio_writedatablock(f_descriptor, "solution", (void*)(array1),
298             &isize, "double", phasta_iotype );
299     }
300     field_flag++;
301 
302     for ( i = 0; i < nppp; i++) {
303         // Write solution field ...
304         isize = (*nshg)*(*numVars);
305         nitems = 3;
306         iarray[ 0 ] = (*nshg);
307         iarray[ 1 ] = (*numVars);
308         iarray[ 2 ] = (*stepno);
309         phio_writeheader(f_descriptor, "time derivative of solution",
310             (void*)iarray, &nitems, &isize, "double", phasta_iotype);
311         nitems = (*nshg)*(*numVars);
312         phio_writedatablock(f_descriptor, "time derivative of solution",
313             (void*)(array2), &isize, "double", phasta_iotype );
314     }
315     field_flag++;
316 
317     if (field_flag==nfields){
318       phio_closefile_write(f_descriptor);
319       if (*pid==0) {
320         printf("\n");
321       }
322     }
323 }
324 
325 void
326 Write_Error(  int* pid,
327               int* stepno,
328               int* nshg,
329               int* numVars,
330               double* array1 ) {
331     char fname[255];
332     char rfile[60];
333     int irstou;
334     int magic_number = 362436;
335     int* mptr = &magic_number;
336     double version=0.0;
337     int isize, nitems;
338     int iarray[10];
339     int nfiles;
340     int nfields;
341     int numparts;
342     int irank;
343     int nprocs;
344 
345     nfiles = outpar.nsynciofiles;
346     nfields = outpar.nsynciofieldswriterestart;
347     numparts = workfc.numpe;
348     irank = *pid; //workfc.myrank;
349     nprocs = workfc.numpe;
350 
351     // Calculate number of parts each  proc deal with and where it start and end ...
352     int nppp = numparts/nprocs;// nppp : Number of parts per proc ...
353     int startpart = irank * nppp +1;// Part id from which I (myrank) start ...
354     int endpart = startpart + nppp - 1;// Part id to which I (myrank) end ...
355 
356     field_flag++;
357 
358     int i;
359     for ( i = 0; i < nppp; i++  ) {
360 
361         if(*pid==0) {
362           printf("\n");
363           printf("The %d/%d th field to be written is 'errors'\n",field_flag,nfields);
364         }
365 
366         isize = (*nshg)*(*numVars);
367         nitems = 3;
368         iarray[ 0 ] = (*nshg);
369         iarray[ 1 ] = (*numVars);
370         iarray[ 2 ] = (*stepno);
371 
372         phio_writeheader(f_descriptor, "errors", (void*)iarray, &nitems,
373             &isize, "double", phasta_iotype);
374 
375         phio_writedatablock(f_descriptor, "errors", (void*)array1, &isize,
376             "double", phasta_iotype );
377     }
378     if (field_flag==nfields){
379       phio_closefile_write(f_descriptor);
380       if (*pid==0) {
381         printf("Last field %d 'errors' finished! \n",nfields);
382         printf("\n");
383       }
384     }
385 }
386 
387 void
388 Write_Displ(  int* pid,
389               int* stepno,
390               int* nshg,
391               int* numVars,
392               double* array1 ) {
393   fprintf(stderr, "This function is dead...exiting\n");
394   exit(1);
395 }
396 
397 void
398 Write_Field(  int *pid,
399               char* filemode,
400               char* fieldtag,
401               int* tagsize,
402               void* array,
403               char* arraytype,
404               int* nshg,
405               int* numvars,
406               int* stepno) {
407     char *fieldlabel = (char *)malloc((*tagsize+1)*sizeof(char));
408     strncpy(fieldlabel, fieldtag, *tagsize);
409     fieldlabel[*tagsize] = '\0';
410 
411     int irstou;
412     int magic_number = 362436;
413     int* mptr = &magic_number;
414     double version=0.0;
415     int isize, nitems;
416     int iarray[10];
417 
418     char fmode[10];
419     if(!strncmp(filemode,"w",1))
420       strcpy(fmode,"write");
421     else // default is append
422       strcpy(fmode,"append");
423 
424     char datatype[10];
425     if(!strncmp(arraytype,"i",1))
426       strcpy(datatype,"int");
427     else // default is double
428       strcpy(datatype,"double");
429 
430     int nfiles;
431     int nfields;
432     int numparts;
433     int irank;
434     int nprocs;
435 
436     nfiles = outpar.nsynciofiles;
437     nfields = outpar.nsynciofieldswriterestart;
438     numparts = workfc.numpe;
439     irank = *pid; //workfc.myrank;
440     nprocs = workfc.numpe;
441 
442     // Calculate number of parts each  proc deal with and where it start and end ...
443     int nppp = numparts/nprocs;// nppp : Number of parts per proc ...
444     int startpart = irank * nppp +1;// Part id from which I (myrank) start ...
445     int endpart = startpart + nppp - 1;// Part id to which I (myrank) end ...
446 
447     strncpy(fieldlabel, fieldtag, *tagsize);
448 
449     field_flag++;
450     if(*pid==0) {
451       printf("\n");
452       printf("The %d/%d th field to be written is '%s'\n",field_flag,nfields,fieldlabel);
453     }
454 
455     int i;
456     for ( i = 0; i < nppp; i++  ) {
457         // Write solution field ...
458         isize = (*nshg)*(*numvars);
459         nitems = 3;
460         iarray[ 0 ] = (*nshg);
461         iarray[ 1 ] = (*numvars);
462         iarray[ 2 ] = (*stepno);
463 
464         phio_writeheader(f_descriptor, fieldlabel, (void*)iarray, &nitems,
465             &isize, datatype, phasta_iotype);
466         nitems = (*nshg)*(*numvars);
467         phio_writedatablock(f_descriptor, fieldlabel, array, &isize,
468             datatype, phasta_iotype );
469     }
470     if (field_flag==nfields){
471       phio_closefile_write(f_descriptor);
472       if (*pid==0) {
473         printf("Last field %d '%s' finished! \n",nfields, fieldtag);
474         printf("\n");
475       }
476     }
477     free(fieldlabel);
478 }
479 
480 void
481 Write_PhAvg2( int* pid,
482               char* filemode,
483               char* fieldtag,
484               int* tagsize,
485               int* iphase,
486               int* nphasesincycle,
487               void* array,
488               char* arraytype,
489               int* nshg,
490               int* numvars,
491               int* stepno) {
492     int addtagsize=0; // phase number is added to the name of the field
493     if(*iphase<10)
494       addtagsize=1;
495     else if(*iphase<100)
496       addtagsize=2;
497     else if(*iphase<1000)
498       addtagsize=3;
499 
500     int tagsize2;
501     tagsize2=*tagsize+addtagsize;
502 
503     char *fieldlabel = (char *)malloc((tagsize2+1)*sizeof(char));
504     strncpy(fieldlabel, fieldtag, *tagsize);
505     fieldlabel[tagsize2] = '\0';
506 
507     char straddtagsize[10];
508     sprintf(straddtagsize,"%d",*iphase);
509 
510     if(*iphase<10) {
511       fieldlabel[tagsize2-1]=straddtagsize[0];
512     }
513     else if(*iphase<100) {
514       fieldlabel[tagsize2-2]=straddtagsize[0];
515       fieldlabel[tagsize2-1]=straddtagsize[1];
516     }
517     else if(*iphase<1000) {
518       fieldlabel[tagsize2-3]=straddtagsize[0];
519       fieldlabel[tagsize2-2]=straddtagsize[1];
520       fieldlabel[tagsize2-1]=straddtagsize[2];
521     }
522 
523     int irstou;
524     int magic_number = 362436;
525     int* mptr = &magic_number;
526     double version=0.0;
527     int isize, nitems;
528     int iarray[10];
529 
530     char fmode[10];
531     if(!strncmp(filemode,"w",1))
532       strcpy(fmode,"write");
533     else // default is append
534       strcpy(fmode,"append");
535 
536     char datatype[10];
537     if(!strncmp(arraytype,"i",1))
538       strcpy(datatype,"int");
539     else // default is double
540       strcpy(datatype,"double");
541 
542     int nfiles;
543     int nfields;
544     int numparts;
545     int irank;
546     int nprocs;
547 
548     nfiles = outpar.nsynciofiles;
549     nfields = outpar.nsynciofieldswriterestart;
550     numparts = workfc.numpe;
551     irank = *pid; //workfc.myrank;
552     nprocs = workfc.numpe;
553 
554     // Calculate number of parts each  proc deal with and where it start and end ...
555     int nppp = numparts/nprocs;// nppp : Number of parts per proc ...
556     int startpart = irank * nppp +1;// Part id from which I (myrank) start ...
557     int endpart = startpart + nppp - 1;// Part id to which I (myrank) end ...
558 
559     field_flag++;
560     if(*pid==0) {
561       printf("\n");
562       printf("The %d/%d th field to be written is '%s'\n",field_flag,nfields,fieldlabel);
563     }
564 
565     int i;
566     for ( i = 0; i < nppp; i++  ) {
567         // Write solution field ...
568         isize = (*nshg)*(*numvars);
569         nitems = 3;
570         iarray[ 0 ] = (*nshg);
571         iarray[ 1 ] = (*numvars);
572         iarray[ 2 ] = (*stepno);
573         phio_writeheader(f_descriptor, fieldlabel, (void*)iarray, &nitems,
574             &isize, "double", phasta_iotype);
575         nitems = (*nshg)*(*numvars);
576         phio_writedatablock(f_descriptor, fieldlabel, array, &isize,
577             "double", phasta_iotype );
578     }
579     if (field_flag==nfields){
580       phio_closefile_write(f_descriptor);
581       if (*pid==0) {
582         printf("\n");
583       }
584     }
585     free(fieldlabel);
586 }
587 
588 
589 void
590 Write_d2wall(   int* pid,
591                 int* numnp,
592                 double* array1 ) {
593     int isize, nitems;
594     int iarray[10];
595     int nfiles;
596     int nfields;
597     int numparts;
598     int irank;
599     int nprocs;
600 
601     //  Retrieve and compute the parameters required for SyncIO
602     nfiles = outpar.nsynciofiles;
603     nfields = 1; //outpar.nsynciofieldswriterestart;  // Only the distance to the walls in d2wall
604     numparts = workfc.numpe;
605     irank = *pid; //workfc.myrank;
606     nprocs = workfc.numpe;
607 
608     // Calculate number of parts each proc deal with and where it start and end ...
609     int nppp = numparts/nprocs;// nppp : Number of parts per proc ...
610     int startpart = irank * nppp +1;// Part id from which I (myrank) start ...
611     int endpart = startpart + nppp - 1;// Part id to which I (myrank) end ...
612 
613     int descriptor;
614     int nppf=numparts/nfiles;
615 
616     phio_openfile_write("d2wall.", &nfiles, &nfields, &nppf, &f_descriptor);
617 
618     field_flag=0;
619 
620      int i;
621      for ( i = 0; i < nppp; i++) { //This loop is useful only if several parts per processor
622         // Write solution field ...
623         isize = (*numnp);
624         nitems = 2;
625         iarray[ 0 ] = (*numnp);
626         iarray[ 1 ] = 1; //numVars = 1
627 
628         phio_writeheader(f_descriptor, "d2wall", (void*)iarray, &nitems,
629             &isize, "double", phasta_iotype);
630         phio_writedatablock(f_descriptor, "d2wall", (void*)(array1), &isize,
631             "double", phasta_iotype );
632     }
633     field_flag++;
634 
635     if (field_flag==nfields){
636       phio_closefile_write(f_descriptor);
637       if (irank==0) {
638         printf("\n");
639       }
640     }
641 }
642