xref: /phasta/phastaIO/phastaIO.cc (revision fe03b87fb3879205ff37d360f6950bfcc4ac4941)
1 /*
2  *
3  * This is the SyncIO library that uses MPI-IO collective functions to
4  * implement a flexible I/O checkpoint solution for a large number of
5  * processors.
6  *
7  * Previous developer: Ning Liu         (liun2@cs.rpi.edu)
8  *                     Jing Fu          (fuj@cs.rpi.edu)
9  * Current developers: Michel Rasquin   (Michel.Rasquin@colorado.edu),
10  *                     Ben Matthews     (benjamin.a.matthews@colorado.edu)
11  *
12  */
13 
14 #include <map>
15 #include <vector>
16 #include <string>
17 #include <string.h>
18 #include <ctype.h>
19 #include <stdlib.h>
20 #include <stdio.h>
21 #include <math.h>
22 #include <sstream>
23 #include "phastaIO.h"
24 #include "mpi.h"
25 #include "phiotmrc.h"
26 
27 #define VERSION_INFO_HEADER_SIZE 8192
28 #define DB_HEADER_SIZE 1024
29 #define ONE_MEGABYTE 1048576
30 #define TWO_MEGABYTE 2097152
31 #define ENDIAN_TEST_NUMBER 12180 // Troy's Zip Code!!
32 #define MAX_PHASTA_FILES 64
33 #define MAX_PHASTA_FILE_NAME_LENGTH 1024
34 #define MAX_FIELDS_NUMBER ((VERSION_INFO_HEADER_SIZE/MAX_FIELDS_NAME_LENGTH)-4) // The meta data include - MPI_IO_Tag, nFields, nFields*names of the fields, nppf
35 // -3 for MPI_IO_Tag, nFields and nppf, -4 for extra security (former nFiles)
36 #define MAX_FIELDS_NAME_LENGTH 128
37 #define DefaultMHSize (4*ONE_MEGABYTE)
38 //#define DefaultMHSize (8350) //For test
39 #define LATEST_WRITE_VERSION 1
40 #define inv1024sq 953.674316406e-9 // = 1/1024/1024
41 int MasterHeaderSize = -1;
42 
43 bool PRINT_PERF = false; // default print no perf results
44 int irank = -1; // global rank, should never be manually manipulated
45 int mysize = -1;
46 
47 // Static variables are bad but used here to store the subcommunicator and associated variables
48 // Prevent MPI_Comm_split to be called more than once, especially on BGQ with the V1R2M1 driver (leak detected in MPI_Comm_split - IBM working on it)
49 static int s_assign_local_comm = 0;
50 static MPI_Comm s_local_comm;
51 static int s_local_size = -1;
52 static int s_local_rank = -1;
53 
54 // the following defines are for debug printf
55 #define PHASTAIO_DEBUG 0 //default to not print any debugging info
56 
57 void phprintf(const char* fmt, ...) {
58 #if PHASTAIO_DEBUG
59   char format[1024];
60   snprintf(format, sizeof(format), "phastaIO debug: %s", fmt);
61   va_list ap;
62   va_start(ap,fmt);
63   vprintf(format,ap)
64   va_end(ap);
65 #endif
66 }
67 
68 void phprintf_0(const char* fmt, ...) {
69 #if PHASTAIO_DEBUG
70   int rank = 0;
71   MPI_Comm_rank(MPI_COMM_WORLD, &rank);
72   if(rank == 0){
73     char format[1024];
74     snprintf(format, sizeof(format), "phastaIO debug: irank=0 %s", fmt);
75     va_list ap;
76     va_start(ap,s);
77     vprintf(format, ap);
78     va_end(ap);
79   }
80 #endif
81 }
82 
83 enum PhastaIO_Errors
84 {
85 	MAX_PHASTA_FILES_EXCEEDED = -1,
86 	UNABLE_TO_OPEN_FILE = -2,
87 	NOT_A_MPI_FILE = -3,
88 	GPID_EXCEEDED = -4,
89 	DATA_TYPE_ILLEGAL = -5,
90 };
91 
92 using namespace std;
93 
94 namespace{
95 
96         map<int, std::string> LastHeaderKey;
97 	vector< FILE* > fileArray;
98 	vector< bool > byte_order;
99 	vector< int > header_type;
100 	int DataSize=0;
101 	bool LastHeaderNotFound = false;
102 	bool Wrong_Endian = false ;
103 	bool Strict_Error = false ;
104 	bool binary_format = true;
105 
106 	/***********************************************************************/
107 	/***************** NEW PHASTA IO CODE STARTS HERE **********************/
108 	/***********************************************************************/
109 
110 	typedef struct
111 	{
112 		char filename[MAX_PHASTA_FILE_NAME_LENGTH];   /* defafults to 1024 */
113 		unsigned long my_offset;
114 		unsigned long next_start_address;
115 		unsigned long **my_offset_table;
116 		unsigned long **my_read_table;
117 
118 		double * double_chunk;
119 		double * read_double_chunk;
120 
121 		int field_count;
122 		int part_count;
123 		int read_field_count;
124 		int read_part_count;
125 		int GPid;
126 		int start_id;
127 
128 		int mhsize;
129 
130 		int myrank;
131 		int numprocs;
132 		int local_myrank;
133 		int local_numprocs;
134 
135 		int nppp;
136 		int nPPF;
137 		int nFiles;
138 		int nFields;
139 
140 		int * int_chunk;
141 		int * read_int_chunk;
142 
143 		int Wrong_Endian;                            /* default to false */
144 		char * master_header;
145 		MPI_File file_handle;
146 		MPI_Comm local_comm;
147 	} phastaio_file_t;
148 
149 	typedef struct
150 	{
151 		int nppf, nfields;
152 		char * masterHeader;
153 	}serial_file;
154 
155 	serial_file *SerialFile;
156 	phastaio_file_t *PhastaIOActiveFiles[MAX_PHASTA_FILES];
157 	int PhastaIONextActiveIndex = 0; /* indicates next index to allocate */
158 
159 	// the caller has the responsibility to delete the returned string
160 	// TODO: StringStipper("nbc value? ") returns NULL?
161 	char*
162 		StringStripper( const char  istring[] ) {
163 
164 			int length = strlen( istring );
165 
166 			char* dest = (char *)malloc( length + 1 );
167 
168 			strcpy( dest, istring );
169 			dest[ length ] = '\0';
170 
171 			if ( char* p = strpbrk( dest, " ") )
172 				*p = '\0';
173 
174 			return dest;
175 		}
176 
177 	inline int
178 		cscompare( const char teststring[],
179 				const char targetstring[] ) {
180 
181 			char* s1 = const_cast<char*>(teststring);
182 			char* s2 = const_cast<char*>(targetstring);
183 
184 			while( *s1 == ' ') s1++;
185 			while( *s2 == ' ') s2++;
186 			while( ( *s1 )
187 					&& ( *s2 )
188 					&& ( *s2 != '?')
189 					&& ( tolower( *s1 )==tolower( *s2 ) ) ) {
190 				s1++;
191 				s2++;
192 				while( *s1 == ' ') s1++;
193 				while( *s2 == ' ') s2++;
194 			}
195 			if ( !( *s1 ) || ( *s1 == '?') ) return 1;
196 			else return 0;
197 		}
198 
199 	inline void
200 		isBinary( const char iotype[] ) {
201 
202 			char* fname = StringStripper( iotype );
203 			if ( cscompare( fname, "binary" ) ) binary_format = true;
204 			else binary_format = false;
205 			free (fname);
206 
207 		}
208 
209 	inline size_t
210 		typeSize( const char typestring[] ) {
211 
212 			char* ts1 = StringStripper( typestring );
213 
214 			if ( cscompare( "integer", ts1 ) ) {
215 				free (ts1);
216 				return sizeof(int);
217 			} else if ( cscompare( "double", ts1 ) ) {
218 				free (ts1);
219 				return sizeof( double );
220 			} else {
221 				free (ts1);
222 				fprintf(stderr,"unknown type : %s\n",ts1);
223 				return 0;
224 			}
225 		}
226 
227 	int
228 		readHeader( FILE*       fileObject,
229 				const char  phrase[],
230 				int*        params,
231 				int         expect ) {
232 
233 			char* text_header;
234 			char* token;
235 			char Line[1024];
236 			char junk;
237 			bool FOUND = false ;
238 			int real_length;
239 			int skip_size, integer_value;
240 			int rewind_count=0;
241 
242 			if( !fgets( Line, 1024, fileObject ) && feof( fileObject ) ) {
243 				rewind( fileObject );
244 				clearerr( fileObject );
245 				rewind_count++;
246 				fgets( Line, 1024, fileObject );
247 			}
248 
249 			while( !FOUND  && ( rewind_count < 2 ) )  {
250 				if ( ( Line[0] != '\n' ) && ( real_length = strcspn( Line, "#" )) ){
251 					text_header = (char *)malloc( real_length + 1 );
252 
253 					strncpy( text_header, Line, real_length );
254 					text_header[ real_length ] =static_cast<char>(NULL);
255 					token = strtok ( text_header, ":" );
256 					//if( cscompare( phrase , token ) ) {
257                                         // Double comparison required because different fields can still start
258                                         // with the same name for mixed meshes (nbc code, nbc values, etc).
259                                         // Would be easy to fix cscompare instead but would it break sth else?
260 					if( cscompare( phrase , token ) && cscompare( token , phrase ) ) {
261 						FOUND = true ;
262 						token = strtok( NULL, " ,;<>" );
263 						skip_size = atoi( token );
264 						int i;
265 						for( i=0; i < expect && ( token = strtok( NULL," ,;<>") ); i++) {
266 							params[i] = atoi( token );
267 						}
268 						if ( i < expect ) {
269 							fprintf(stderr,"Aloha Expected # of ints not found for: %s\n",phrase );
270 						}
271 					} else if ( cscompare(token,"byteorder magic number") ) {
272 						if ( binary_format ) {
273 							fread((void*)&integer_value,sizeof(int),1,fileObject);
274 							fread( &junk, sizeof(char), 1 , fileObject );
275 							if ( 362436 != integer_value ) Wrong_Endian = true;
276 						} else{
277 							fscanf(fileObject, "%d\n", &integer_value );
278 						}
279 					} else {
280 						/* some other header, so just skip over */
281 						token = strtok( NULL, " ,;<>" );
282 						skip_size = atoi( token );
283 						if ( binary_format)
284 							fseek( fileObject, skip_size, SEEK_CUR );
285 						else
286 							for( int gama=0; gama < skip_size; gama++ )
287 								fgets( Line, 1024, fileObject );
288 					}
289 					free (text_header);
290 				} // end of if before while loop
291 
292 				if ( !FOUND )
293 					if( !fgets( Line, 1024, fileObject ) && feof( fileObject ) ) {
294 						rewind( fileObject );
295 						clearerr( fileObject );
296 						rewind_count++;
297 						fgets( Line, 1024, fileObject );
298 					}
299 			}
300 
301 			if ( !FOUND ) {
302 				//fprintf(stderr, "Error: Could not find: %s\n", phrase);
303 				if(irank == 0) printf("WARNING: Could not find: %s\n", phrase);
304 				return 1;
305 			}
306 			return 0;
307 		}
308 
309 } // end unnamed namespace
310 
311 
312 // begin of publicly visible functions
313 
314 /**
315  * This function takes a long long pointer and assign (start) phiotmrc value to it
316  */
317 void startTimer(double* start) {
318 
319         if( !PRINT_PERF ) return;
320 
321 	MPI_Barrier(MPI_COMM_WORLD);
322 	*start =  phiotmrc();
323 }
324 
325 /**
326  * This function takes a long long pointer and assign (end) phiotmrc value to it
327  */
328 void endTimer(double* end) {
329 
330         if( !PRINT_PERF ) return;
331 
332 	*end = phiotmrc();
333 	MPI_Barrier(MPI_COMM_WORLD);
334 }
335 
336 /**
337  * choose to print some performance results (or not) according to
338  * the PRINT_PERF macro
339  */
340 void printPerf(
341 		const char* func_name,
342 		double start,
343 		double end,
344 		unsigned long datasize,
345 		int printdatainfo,
346 		const char* extra_msg) {
347 
348 	if( !PRINT_PERF ) return;
349 
350 	unsigned long data_size = datasize;
351 
352 	double time = end - start;
353 
354 	unsigned long isizemin,isizemax,isizetot;
355 	double sizemin,sizemax,sizeavg,sizetot,rate;
356 	double tmin, tmax, tavg, ttot;
357 
358 	MPI_Allreduce(&time, &tmin,1, MPI_DOUBLE, MPI_MIN, MPI_COMM_WORLD);
359 	MPI_Allreduce(&time, &tmax,1, MPI_DOUBLE, MPI_MAX, MPI_COMM_WORLD);
360 	MPI_Allreduce(&time, &ttot,1, MPI_DOUBLE, MPI_SUM, MPI_COMM_WORLD);
361 	tavg = ttot/mysize;
362 
363 	if(irank == 0) {
364 		if ( PhastaIONextActiveIndex == 0 ) printf("** 1PFPP ");
365 		else  printf("** syncIO ");
366 		printf("%s(): Tmax = %f sec, Tmin = %f sec, Tavg = %f sec", func_name, tmax, tmin, tavg);
367 	}
368 
369 	if(printdatainfo == 1) { // if printdatainfo ==1, compute I/O rate and block size
370 		MPI_Allreduce(&data_size,&isizemin,1,MPI_LONG_LONG_INT,MPI_MIN,MPI_COMM_WORLD);
371 		MPI_Allreduce(&data_size,&isizemax,1,MPI_LONG_LONG_INT,MPI_MAX,MPI_COMM_WORLD);
372 		MPI_Allreduce(&data_size,&isizetot,1,MPI_LONG_LONG_INT,MPI_SUM,MPI_COMM_WORLD);
373 
374 		sizemin=(double)(isizemin*inv1024sq);
375 		sizemax=(double)(isizemax*inv1024sq);
376 		sizetot=(double)(isizetot*inv1024sq);
377 		sizeavg=(double)(1.0*sizetot/mysize);
378 		rate=(double)(1.0*sizetot/tmax);
379 
380 		if( irank == 0) {
381 			printf(", Rate = %f MB/s [%s] \n \t\t\t block size: Min= %f MB; Max= %f MB; Avg= %f MB; Tot= %f MB\n", rate, extra_msg, sizemin,sizemax,sizeavg,sizetot);
382 		}
383 	}
384 	else {
385 		if(irank == 0) {
386 			printf(" \n");
387 			//printf(" (%s) \n", extra_msg);
388 		}
389 	}
390 }
391 
392 /**
393  * This function is normally called at the beginning of a read operation, before
394  * init function.
395  * This function (uses rank 0) reads out nfields, nppf, master header size,
396  * endianess and allocates for masterHeader string.
397  * These values are essential for following read operations. Rank 0 will bcast
398  * these values to other ranks in the commm world
399  *
400  * If the file set is of old POSIX format, it would throw error and exit
401  */
402 void queryphmpiio(const char filename[],int *nfields, int *nppf)
403 {
404 	MPI_Comm_rank(MPI_COMM_WORLD, &irank);
405 	MPI_Comm_size(MPI_COMM_WORLD, &mysize);
406 
407 	if(irank == 0) {
408 		FILE * fileHandle;
409 		char* fname = StringStripper( filename );
410 
411 		fileHandle = fopen (fname,"rb");
412 		if (fileHandle == NULL ) {
413 			printf("\nError: File %s doesn't exist! Please check!\n",fname);
414 		}
415 		else {
416 			SerialFile =(serial_file *)calloc( 1,  sizeof( serial_file) );
417 			int meta_size_limit = VERSION_INFO_HEADER_SIZE;
418 			SerialFile->masterHeader = (char *)malloc( meta_size_limit );
419 			fread(SerialFile->masterHeader, 1, meta_size_limit, fileHandle);
420 
421 			char read_out_tag[MAX_FIELDS_NAME_LENGTH];
422 			char version[MAX_FIELDS_NAME_LENGTH/4];
423 			int mhsize;
424 			char * token;
425 			int magic_number;
426 
427 			memcpy( read_out_tag,
428 					SerialFile->masterHeader,
429 					MAX_FIELDS_NAME_LENGTH-1 );
430 
431 			if ( cscompare ("MPI_IO_Tag",read_out_tag) ) {
432 				// Test endianess ...
433 				memcpy (&magic_number,
434 						SerialFile->masterHeader + sizeof("MPI_IO_Tag : ")-1, //-1 sizeof returns the size of the string+1 for "\0"
435 						sizeof(int) );                                        // masterheader should look like "MPI_IO_Tag : 12180 " with 12180 in binary format
436 
437 				if ( magic_number != ENDIAN_TEST_NUMBER ) {
438 					printf("Endian is different!\n");
439 					// Will do swap later
440 				}
441 
442 				// test version, old version, default masterheader size is 4M
443 				// newer version masterheader size is read from first line
444 				memcpy(version,
445 						SerialFile->masterHeader + MAX_FIELDS_NAME_LENGTH/2,
446 						MAX_FIELDS_NAME_LENGTH/4 - 1); //TODO: why -1?
447 
448 				if( cscompare ("version",version) ) {
449 					// if there is "version" tag in the file, then it is newer format
450 					// read master header size from here, otherwise use default
451 					// Note: if version is "1", we know mhsize is at 3/4 place...
452 
453 					token = strtok(version, ":");
454 					token = strtok(NULL, " ,;<>" );
455 					int iversion = atoi(token);
456 
457 					if( iversion == 1) {
458 						memcpy( &mhsize,
459 								SerialFile->masterHeader + MAX_FIELDS_NAME_LENGTH/4*3 + sizeof("mhsize : ")-1,
460 								sizeof(int));
461 						if ( magic_number != ENDIAN_TEST_NUMBER )
462 							SwapArrayByteOrder(&mhsize, sizeof(int), 1);
463 
464 						if( mhsize > DefaultMHSize ) {
465 							//if actual headersize is larger than default, let's re-read
466 							free(SerialFile->masterHeader);
467 							SerialFile->masterHeader = (char *)malloc(mhsize);
468 							fseek(fileHandle, 0, SEEK_SET); // reset the file stream position
469 							fread(SerialFile->masterHeader,1,mhsize,fileHandle);
470 						}
471 					}
472 					//TODO: check if this is a valid int??
473 					MasterHeaderSize = mhsize;
474 				}
475 				else { // else it's version 0's format w/o version tag, implicating MHSize=4M
476 					MasterHeaderSize = DefaultMHSize;
477 				}
478 
479 				memcpy( read_out_tag,
480 						SerialFile->masterHeader+MAX_FIELDS_NAME_LENGTH+1,
481 						MAX_FIELDS_NAME_LENGTH ); //TODO: why +1
482 
483 				// Read in # fields ...
484 				token = strtok( read_out_tag, ":" );
485 				token = strtok( NULL," ,;<>" );
486 				*nfields = atoi( token );
487 				if ( *nfields > MAX_FIELDS_NUMBER) {
488 					printf("Error queryphmpiio: nfields is larger than MAX_FIELDS_NUMBER!\n");
489 				}
490 				SerialFile->nfields=*nfields; //TODO: sanity check of this int?
491 
492 				memcpy( read_out_tag,
493 						SerialFile->masterHeader + MAX_FIELDS_NAME_LENGTH * 2
494 						+ *nfields * MAX_FIELDS_NAME_LENGTH,
495 						MAX_FIELDS_NAME_LENGTH);
496 
497 				token = strtok( read_out_tag, ":" );
498 				token = strtok( NULL," ,;<>" );
499 				*nppf = atoi( token );
500 				SerialFile->nppf=*nppf; //TODO: sanity check of int
501 			} // end of if("MPI_IO_TAG")
502 			else {
503 				printf("Error queryphmpiio: The file you opened is not of syncIO new format, please check! read_out_tag = %s\n",read_out_tag);
504                                 exit(1);
505 			}
506 			fclose(fileHandle);
507                         free(SerialFile->masterHeader);
508                         free(SerialFile);
509 		} //end of else
510                 free(fname);
511 	}
512 
513 	// Bcast value to every one
514 	MPI_Bcast( nfields, 1, MPI_INT, 0, MPI_COMM_WORLD );
515 	MPI_Bcast( nppf, 1, MPI_INT, 0, MPI_COMM_WORLD );
516 	MPI_Bcast( &MasterHeaderSize, 1, MPI_INT, 0, MPI_COMM_WORLD );
517 	phprintf("Info queryphmpiio: myrank = %d, MasterHeaderSize = %d", irank, MasterHeaderSize);
518 }
519 
520 /**
521  * This function computes the right master header size (round to size of 2^n).
522  * This is only needed for file format version 1 in "write" mode.
523  */
524 int computeMHSize(int nfields, int nppf, int version) {
525 	int mhsize;
526 	if(version == 1) {
527 		//int meta_info_size = (2+nfields+1) * MAX_FIELDS_NAME_LENGTH; // 2 is MPI_IO_TAG and nFields, the others 1 is nppf
528 		int meta_info_size = VERSION_INFO_HEADER_SIZE;
529 		int actual_size =  nfields * nppf * sizeof(unsigned long) + meta_info_size;
530 		//printf("actual_size = %d, offset table size = %d\n", actual_size,  nfields * nppf * sizeof(long long));
531 		if (actual_size > DefaultMHSize) {
532 			mhsize = (int) ceil( (double) actual_size/DefaultMHSize); // it's rounded to ceiling of this value
533 			mhsize *= DefaultMHSize;
534 		}
535 		else {
536 			mhsize = DefaultMHSize;
537 		}
538 	}
539 	return mhsize;
540 }
541 
542 /**
543  * Computes correct color of a rank according to number of files.
544  */
545 extern "C" int computeColor( int myrank, int numprocs, int nfiles) {
546 	int color =
547 		(int)(myrank / (numprocs / nfiles));
548 	return color;
549 }
550 
551 
552 /**
553  * Check the file descriptor.
554  */
555 void checkFileDescriptor(const char fctname[],
556                          int*  fileDescriptor ) {
557 	if ( *fileDescriptor < 0 ) {
558 		printf("Error: File descriptor = %d in %s\n",*fileDescriptor,fctname);
559 		exit(1);
560 	}
561 }
562 
563 /**
564  * Initialize the file struct members and allocate space for file struct
565  * buffers.
566  *
567  * Note: this function is only called when we are using new format. Old POSIX
568  * format should skip this routine and call openfile() directly instead.
569  */
570 int initphmpiio( int *nfields, int *nppf, int *nfiles, int *filehandle, const char mode[])
571 {
572 	// we init irank again in case query not called (e.g. syncIO write case)
573 	MPI_Comm_rank(MPI_COMM_WORLD, &irank);
574 	MPI_Comm_size(MPI_COMM_WORLD, &mysize);
575 
576 	phprintf("Info initphmpiio: entering function, myrank = %d, MasterHeaderSize = %d", irank, MasterHeaderSize);
577 
578 	double timer_start, timer_end;
579 	startTimer(&timer_start);
580 
581 	char* imode = StringStripper( mode );
582 
583 	// Note: if it's read, we presume query was called prior to init and
584 	// MasterHeaderSize is already set to correct value from parsing header
585 	// otherwise it's write then it needs some computation to be set
586 	if ( cscompare( "read", imode ) ) {
587 		// do nothing
588 	}
589 	else if( cscompare( "write", imode ) ) {
590 		MasterHeaderSize =  computeMHSize(*nfields, *nppf, LATEST_WRITE_VERSION);
591 	}
592 	else {
593 		printf("Error initphmpiio: can't recognize the mode %s", imode);
594                 exit(1);
595 	}
596         free ( imode );
597 
598 	phprintf("Info initphmpiio: myrank = %d, MasterHeaderSize = %d", irank, MasterHeaderSize);
599 
600 	int i, j;
601 
602 	if( PhastaIONextActiveIndex == MAX_PHASTA_FILES ) {
603 		printf("Error initphmpiio: PhastaIONextActiveIndex = MAX_PHASTA_FILES");
604                 endTimer(&timer_end);
605 	        printPerf("initphmpiio", timer_start, timer_end, 0, 0, "");
606 		return MAX_PHASTA_FILES_EXCEEDED;
607 	}
608 	//		else if( PhastaIONextActiveIndex == 0 )  //Hang in debug mode on Intrepid
609 	//		{
610 	//			for( i = 0; i < MAX_PHASTA_FILES; i++ );
611 	//			{
612 	//				PhastaIOActiveFiles[i] = NULL;
613 	//			}
614 	//		}
615 
616 
617 	PhastaIOActiveFiles[PhastaIONextActiveIndex] = (phastaio_file_t *)calloc( 1,  sizeof( phastaio_file_t) );
618 
619 	i = PhastaIONextActiveIndex;
620 	PhastaIONextActiveIndex++;
621 
622 	//PhastaIOActiveFiles[i]->next_start_address = 2*TWO_MEGABYTE;
623 
624 	PhastaIOActiveFiles[i]->next_start_address = MasterHeaderSize;  // what does this mean??? TODO
625 
626 	PhastaIOActiveFiles[i]->Wrong_Endian = false;
627 
628 	PhastaIOActiveFiles[i]->nFields = *nfields;
629 	PhastaIOActiveFiles[i]->nPPF = *nppf;
630 	PhastaIOActiveFiles[i]->nFiles = *nfiles;
631 	MPI_Comm_rank(MPI_COMM_WORLD, &(PhastaIOActiveFiles[i]->myrank));
632 	MPI_Comm_size(MPI_COMM_WORLD, &(PhastaIOActiveFiles[i]->numprocs));
633 
634 
635 	if( *nfiles > 1 ) { // split the ranks according to each mpiio file
636 
637 		if ( s_assign_local_comm == 0) { // call mpi_comm_split for the first (and only) time
638 
639 			if (PhastaIOActiveFiles[i]->myrank == 0) printf("Building subcommunicator\n");
640 
641 			int color = computeColor(PhastaIOActiveFiles[i]->myrank, PhastaIOActiveFiles[i]->numprocs, PhastaIOActiveFiles[i]->nFiles);
642 			MPI_Comm_split(MPI_COMM_WORLD,
643 					color,
644 					PhastaIOActiveFiles[i]->myrank,
645 					&(PhastaIOActiveFiles[i]->local_comm));
646 			MPI_Comm_size(PhastaIOActiveFiles[i]->local_comm,
647 					&(PhastaIOActiveFiles[i]->local_numprocs));
648 			MPI_Comm_rank(PhastaIOActiveFiles[i]->local_comm,
649 					&(PhastaIOActiveFiles[i]->local_myrank));
650 
651 			// back up now these variables so that we do not need to call comm_split again
652 			s_local_comm = PhastaIOActiveFiles[i]->local_comm;
653 			s_local_size = PhastaIOActiveFiles[i]->local_numprocs;
654 			s_local_rank = PhastaIOActiveFiles[i]->local_myrank;
655 			s_assign_local_comm = 1;
656 		}
657 		else { // recycle the subcommunicator
658 			if (PhastaIOActiveFiles[i]->myrank == 0) printf("Recycling subcommunicator\n");
659 			PhastaIOActiveFiles[i]->local_comm = s_local_comm;
660 			PhastaIOActiveFiles[i]->local_numprocs = s_local_size;
661 			PhastaIOActiveFiles[i]->local_myrank = s_local_rank;
662 		}
663 	}
664 	else { // *nfiles == 1 here - no need to call mpi_comm_split here
665 
666 		if (PhastaIOActiveFiles[i]->myrank == 0) printf("Bypassing subcommunicator\n");
667 		PhastaIOActiveFiles[i]->local_comm = MPI_COMM_WORLD;
668 		PhastaIOActiveFiles[i]->local_numprocs = PhastaIOActiveFiles[i]->numprocs;
669 		PhastaIOActiveFiles[i]->local_myrank = PhastaIOActiveFiles[i]->myrank;
670 
671 	}
672 
673 	PhastaIOActiveFiles[i]->nppp =
674 		PhastaIOActiveFiles[i]->nPPF/PhastaIOActiveFiles[i]->local_numprocs;
675 
676 	PhastaIOActiveFiles[i]->start_id = PhastaIOActiveFiles[i]->nPPF *
677 		(int)(PhastaIOActiveFiles[i]->myrank/PhastaIOActiveFiles[i]->local_numprocs) +
678 		(PhastaIOActiveFiles[i]->local_myrank * PhastaIOActiveFiles[i]->nppp);
679 
680 	PhastaIOActiveFiles[i]->my_offset_table =
681 		( unsigned long ** ) calloc( MAX_FIELDS_NUMBER , sizeof( unsigned long *) );
682 
683 	PhastaIOActiveFiles[i]->my_read_table =
684 		( unsigned long ** ) calloc( MAX_FIELDS_NUMBER , sizeof( unsigned long *) );
685 
686 	for (j=0; j<*nfields; j++)
687 	{
688 		PhastaIOActiveFiles[i]->my_offset_table[j] =
689 			( unsigned long * ) calloc( PhastaIOActiveFiles[i]->nppp , sizeof( unsigned long) );
690 
691 		PhastaIOActiveFiles[i]->my_read_table[j] =
692 			( unsigned long * ) calloc( PhastaIOActiveFiles[i]->nppp , sizeof( unsigned long) );
693 	}
694 	*filehandle = i;
695 
696 	PhastaIOActiveFiles[i]->master_header = (char *)calloc(MasterHeaderSize,sizeof( char ));
697 	PhastaIOActiveFiles[i]->double_chunk = (double *)calloc(1,sizeof( double ));
698 	PhastaIOActiveFiles[i]->int_chunk = (int *)calloc(1,sizeof( int ));
699 	PhastaIOActiveFiles[i]->read_double_chunk = (double *)calloc(1,sizeof( double ));
700 	PhastaIOActiveFiles[i]->read_int_chunk = (int *)calloc(1,sizeof( int ));
701 
702 	// Time monitoring
703 	endTimer(&timer_end);
704 	printPerf("initphmpiio", timer_start, timer_end, 0, 0, "");
705 
706 	phprintf_0("Info initphmpiio: quiting function");
707 
708 	return i;
709 }
710 
711 /**
712  * Destruct the file struct and free buffers allocated in init function.
713  */
714 void finalizephmpiio( int *fileDescriptor )
715 {
716 	double timer_start, timer_end;
717 	startTimer(&timer_start);
718 
719 	int i, j;
720 	i = *fileDescriptor;
721 	//PhastaIONextActiveIndex--;
722 
723 	/* //free the offset table for this phasta file */
724 	//for(j=0; j<MAX_FIELDS_NUMBER; j++) //Danger: undefined behavior for my_*_table.[j] not allocated or not initialized to NULL
725 	for(j=0; j<PhastaIOActiveFiles[i]->nFields; j++)
726 	{
727 		free( PhastaIOActiveFiles[i]->my_offset_table[j]);
728 		free( PhastaIOActiveFiles[i]->my_read_table[j]);
729 	}
730 	free ( PhastaIOActiveFiles[i]->my_offset_table );
731 	free ( PhastaIOActiveFiles[i]->my_read_table );
732 	free ( PhastaIOActiveFiles[i]->master_header );
733 	free ( PhastaIOActiveFiles[i]->double_chunk );
734 	free ( PhastaIOActiveFiles[i]->int_chunk );
735 	free ( PhastaIOActiveFiles[i]->read_double_chunk );
736 	free ( PhastaIOActiveFiles[i]->read_int_chunk );
737 
738         if( PhastaIOActiveFiles[i]->nFiles > 1 && s_assign_local_comm ) { // the comm was split
739           if (PhastaIOActiveFiles[i]->myrank == 0) printf("Freeing subcommunicator\n");
740           s_assign_local_comm = 0;
741           MPI_Comm_free(&(PhastaIOActiveFiles[i]->local_comm));
742         }
743 
744 	free( PhastaIOActiveFiles[i]);
745 
746 	endTimer(&timer_end);
747 	printPerf("finalizempiio", timer_start, timer_end, 0, 0, "");
748 
749 	PhastaIONextActiveIndex--;
750 }
751 
752 
753 /**
754  * Special init for M2N in order to create a subcommunicator for the reduced solution (requires PRINT_PERF to be false for now)
755  * Initialize the file struct members and allocate space for file struct buffers.
756  *
757  * Note: this function is only called when we are using new format. Old POSIX
758  * format should skip this routine and call openfile() directly instead.
759  */
760 int initphmpiiosub( int *nfields, int *nppf, int *nfiles, int *filehandle, const char mode[],MPI_Comm my_local_comm)
761 {
762 	// we init irank again in case query not called (e.g. syncIO write case)
763 
764 	MPI_Comm_rank(my_local_comm, &irank);
765 	MPI_Comm_size(my_local_comm, &mysize);
766 
767 	phprintf("Info initphmpiio: entering function, myrank = %d, MasterHeaderSize = %d", irank, MasterHeaderSize);
768 
769 	double timer_start, timer_end;
770 	startTimer(&timer_start);
771 
772 	char* imode = StringStripper( mode );
773 
774 	// Note: if it's read, we presume query was called prior to init and
775 	// MasterHeaderSize is already set to correct value from parsing header
776 	// otherwise it's write then it needs some computation to be set
777 	if ( cscompare( "read", imode ) ) {
778 		// do nothing
779 	}
780 	else if( cscompare( "write", imode ) ) {
781 		MasterHeaderSize =  computeMHSize(*nfields, *nppf, LATEST_WRITE_VERSION);
782 	}
783 	else {
784 		printf("Error initphmpiio: can't recognize the mode %s", imode);
785                 exit(1);
786 	}
787         free ( imode );
788 
789 	phprintf("Info initphmpiio: myrank = %d, MasterHeaderSize = %d", irank, MasterHeaderSize);
790 
791 	int i, j;
792 
793 	if( PhastaIONextActiveIndex == MAX_PHASTA_FILES ) {
794 		printf("Error initphmpiio: PhastaIONextActiveIndex = MAX_PHASTA_FILES");
795                 endTimer(&timer_end);
796 	        printPerf("initphmpiio", timer_start, timer_end, 0, 0, "");
797 		return MAX_PHASTA_FILES_EXCEEDED;
798 	}
799 	//		else if( PhastaIONextActiveIndex == 0 )  //Hang in debug mode on Intrepid
800 	//		{
801 	//			for( i = 0; i < MAX_PHASTA_FILES; i++ );
802 	//			{
803 	//				PhastaIOActiveFiles[i] = NULL;
804 	//			}
805 	//		}
806 
807 
808 	PhastaIOActiveFiles[PhastaIONextActiveIndex] = (phastaio_file_t *)calloc( 1,  sizeof( phastaio_file_t) );
809 
810 	i = PhastaIONextActiveIndex;
811 	PhastaIONextActiveIndex++;
812 
813 	//PhastaIOActiveFiles[i]->next_start_address = 2*TWO_MEGABYTE;
814 
815 	PhastaIOActiveFiles[i]->next_start_address = MasterHeaderSize;  // what does this mean??? TODO
816 
817 	PhastaIOActiveFiles[i]->Wrong_Endian = false;
818 
819 	PhastaIOActiveFiles[i]->nFields = *nfields;
820 	PhastaIOActiveFiles[i]->nPPF = *nppf;
821 	PhastaIOActiveFiles[i]->nFiles = *nfiles;
822 	MPI_Comm_rank(my_local_comm, &(PhastaIOActiveFiles[i]->myrank));
823 	MPI_Comm_size(my_local_comm, &(PhastaIOActiveFiles[i]->numprocs));
824 
825 	int color = computeColor(PhastaIOActiveFiles[i]->myrank, PhastaIOActiveFiles[i]->numprocs, PhastaIOActiveFiles[i]->nFiles);
826 	MPI_Comm_split(my_local_comm,
827 			color,
828 			PhastaIOActiveFiles[i]->myrank,
829 			&(PhastaIOActiveFiles[i]->local_comm));
830 	MPI_Comm_size(PhastaIOActiveFiles[i]->local_comm,
831 			&(PhastaIOActiveFiles[i]->local_numprocs));
832 	MPI_Comm_rank(PhastaIOActiveFiles[i]->local_comm,
833 			&(PhastaIOActiveFiles[i]->local_myrank));
834 	PhastaIOActiveFiles[i]->nppp =
835 		PhastaIOActiveFiles[i]->nPPF/PhastaIOActiveFiles[i]->local_numprocs;
836 
837 	PhastaIOActiveFiles[i]->start_id = PhastaIOActiveFiles[i]->nPPF *
838 		(int)(PhastaIOActiveFiles[i]->myrank/PhastaIOActiveFiles[i]->local_numprocs) +
839 		(PhastaIOActiveFiles[i]->local_myrank * PhastaIOActiveFiles[i]->nppp);
840 
841 	PhastaIOActiveFiles[i]->my_offset_table =
842 		( unsigned long ** ) calloc( MAX_FIELDS_NUMBER , sizeof( unsigned long *) );
843 
844 	PhastaIOActiveFiles[i]->my_read_table =
845 		( unsigned long ** ) calloc( MAX_FIELDS_NUMBER , sizeof( unsigned long *) );
846 
847 	for (j=0; j<*nfields; j++)
848 	{
849 		PhastaIOActiveFiles[i]->my_offset_table[j] =
850 			( unsigned long * ) calloc( PhastaIOActiveFiles[i]->nppp , sizeof( unsigned long) );
851 
852 		PhastaIOActiveFiles[i]->my_read_table[j] =
853 			( unsigned long * ) calloc( PhastaIOActiveFiles[i]->nppp , sizeof( unsigned long) );
854 	}
855 	*filehandle = i;
856 
857 	PhastaIOActiveFiles[i]->master_header = (char *)calloc(MasterHeaderSize,sizeof( char ));
858 	PhastaIOActiveFiles[i]->double_chunk = (double *)calloc(1,sizeof( double ));
859 	PhastaIOActiveFiles[i]->int_chunk = (int *)calloc(1,sizeof( int ));
860 	PhastaIOActiveFiles[i]->read_double_chunk = (double *)calloc(1,sizeof( double ));
861 	PhastaIOActiveFiles[i]->read_int_chunk = (int *)calloc(1,sizeof( int ));
862 
863 	// Time monitoring
864 	endTimer(&timer_end);
865 	printPerf("initphmpiiosub", timer_start, timer_end, 0, 0, "");
866 
867 	phprintf_0("Info initphmpiiosub: quiting function");
868 
869 	return i;
870 }
871 
872 
873 
874 /** open file for both POSIX and MPI-IO syncIO format.
875  *
876  * If it's old POSIX format, simply call posix fopen().
877  *
878  * If it's MPI-IO foramt:
879  * in "read" mode, it builds the header table that points to the offset of
880  * fields for parts;
881  * in "write" mode, it opens the file with MPI-IO open routine.
882  */
883 void openfile(const char filename[],
884                const char mode[],
885                int*  fileDescriptor )
886 {
887 	phprintf_0("Info: entering openfile");
888 
889 	double timer_start, timer_end;
890 	startTimer(&timer_start);
891 
892 	if ( PhastaIONextActiveIndex == 0 )
893 	{
894 		FILE* file=NULL ;
895 		*fileDescriptor = 0;
896 		char* fname = StringStripper( filename );
897 		char* imode = StringStripper( mode );
898 
899 		if ( cscompare( "read", imode ) ) file = fopen(fname, "rb" );
900 		else if( cscompare( "write", imode ) ) file = fopen(fname, "wb" );
901 		else if( cscompare( "append", imode ) ) file = fopen(fname, "ab" );
902 
903 		if ( !file ){
904 			fprintf(stderr,"Error openfile: unable to open file %s",fname ) ;
905 		} else {
906 			fileArray.push_back( file );
907 			byte_order.push_back( false );
908 			header_type.push_back( sizeof(int) );
909 			*fileDescriptor = fileArray.size();
910 		}
911 		free (fname);
912 		free (imode);
913 	}
914 	else // else it would be parallel I/O, opposed to posix io
915 	{
916 		char* fname = StringStripper( filename );
917 		char* imode = StringStripper( mode );
918 		int rc;
919                 int i = *fileDescriptor;
920 		checkFileDescriptor("openfile",&i);
921 		char* token;
922 
923 		if ( cscompare( "read", imode ) )
924 		{
925 			//	      if (PhastaIOActiveFiles[i]->myrank == 0)
926 			//                printf("\n **********\nRead open ... ... regular version\n");
927 
928 			rc = MPI_File_open( PhastaIOActiveFiles[i]->local_comm,
929 					fname,
930 					MPI_MODE_RDONLY,
931 					MPI_INFO_NULL,
932 					&(PhastaIOActiveFiles[i]->file_handle) );
933 
934 			if(rc)
935 			{
936 				*fileDescriptor = UNABLE_TO_OPEN_FILE;
937 				printf("Error openfile: Unable to open file %s! File descriptor = %d\n",fname,*fileDescriptor);
938                                 endTimer(&timer_end);
939                                 printPerf("openfile", timer_start, timer_end, 0, 0, "");
940 				return;
941 			}
942 
943 			MPI_Status read_tag_status;
944 			char read_out_tag[MAX_FIELDS_NAME_LENGTH];
945 			int j;
946 			int magic_number;
947 
948 			if ( PhastaIOActiveFiles[i]->local_myrank == 0 ) {
949 				MPI_File_read_at( PhastaIOActiveFiles[i]->file_handle,
950 						0,
951 						PhastaIOActiveFiles[i]->master_header,
952 						MasterHeaderSize,
953 						MPI_CHAR,
954 						&read_tag_status );
955 			}
956 
957 			MPI_Bcast( PhastaIOActiveFiles[i]->master_header,
958 					MasterHeaderSize,
959 					MPI_CHAR,
960 					0,
961 					PhastaIOActiveFiles[i]->local_comm );
962 
963 			memcpy( read_out_tag,
964 					PhastaIOActiveFiles[i]->master_header,
965 					MAX_FIELDS_NAME_LENGTH-1 );
966 
967 			if ( cscompare ("MPI_IO_Tag",read_out_tag) )
968 			{
969 				// Test endianess ...
970 				memcpy ( &magic_number,
971 						PhastaIOActiveFiles[i]->master_header+sizeof("MPI_IO_Tag : ")-1, //-1 sizeof returns the size of the string+1 for "\0"
972 						sizeof(int) );                                                   // masterheader should look like "MPI_IO_Tag : 12180 " with 12180 in binary format
973 
974 				if ( magic_number != ENDIAN_TEST_NUMBER )
975 				{
976 					PhastaIOActiveFiles[i]->Wrong_Endian = true;
977 				}
978 
979 				memcpy( read_out_tag,
980 						PhastaIOActiveFiles[i]->master_header+MAX_FIELDS_NAME_LENGTH+1, // TODO: WHY +1???
981 						MAX_FIELDS_NAME_LENGTH );
982 
983 				// Read in # fields ...
984 				token = strtok ( read_out_tag, ":" );
985 				token = strtok( NULL," ,;<>" );
986 				PhastaIOActiveFiles[i]->nFields = atoi( token );
987 
988 				unsigned long **header_table;
989 				header_table = ( unsigned long ** )calloc(PhastaIOActiveFiles[i]->nFields, sizeof(unsigned long *));
990 
991 				for ( j = 0; j < PhastaIOActiveFiles[i]->nFields; j++ )
992 				{
993 					header_table[j]=( unsigned long * ) calloc( PhastaIOActiveFiles[i]->nPPF , sizeof( unsigned long));
994 				}
995 
996 				// Read in the offset table ...
997 				for ( j = 0; j < PhastaIOActiveFiles[i]->nFields; j++ )
998 				{
999                                         if ( PhastaIOActiveFiles[i]->local_myrank == 0 ) {
1000 						memcpy( header_table[j],
1001 							PhastaIOActiveFiles[i]->master_header +
1002 							VERSION_INFO_HEADER_SIZE +
1003 							j * PhastaIOActiveFiles[i]->nPPF * sizeof(unsigned long),
1004 							PhastaIOActiveFiles[i]->nPPF * sizeof(unsigned long) );
1005                                         }
1006 
1007 					MPI_Scatter( header_table[j],
1008 							PhastaIOActiveFiles[i]->nppp,
1009 							MPI_LONG_LONG_INT,
1010 							PhastaIOActiveFiles[i]->my_read_table[j],
1011 							PhastaIOActiveFiles[i]->nppp,
1012 							MPI_LONG_LONG_INT,
1013 							0,
1014 							PhastaIOActiveFiles[i]->local_comm );
1015 
1016 					// Swap byte order if endianess is different ...
1017 					if ( PhastaIOActiveFiles[i]->Wrong_Endian ) {
1018 						SwapArrayByteOrder( PhastaIOActiveFiles[i]->my_read_table[j],
1019 								sizeof(unsigned long),
1020 								PhastaIOActiveFiles[i]->nppp );
1021 					}
1022 				}
1023 
1024                                 for ( j = 0; j < PhastaIOActiveFiles[i]->nFields; j++ ) {
1025                                 	free ( header_table[j] );
1026                                 }
1027                                 free (header_table);
1028 
1029 			} // end of if MPI_IO_TAG
1030 			else //else not valid MPI file
1031 			{
1032 				*fileDescriptor = NOT_A_MPI_FILE;
1033 				printf("Error openfile: The file %s you opened is not in syncIO new format, please check again! File descriptor = %d, MasterHeaderSize = %d, read_out_tag = %s\n",fname,*fileDescriptor,MasterHeaderSize,read_out_tag);
1034 				//Printing MasterHeaderSize is useful to test a compiler bug on Intrepid BGP
1035                                 endTimer(&timer_end);
1036                                 printPerf("openfile", timer_start, timer_end, 0, 0, "");
1037 				return;
1038 			}
1039 		} // end of if "read"
1040 		else if( cscompare( "write", imode ) )
1041 		{
1042 			rc = MPI_File_open( PhastaIOActiveFiles[i]->local_comm,
1043 					fname,
1044 					MPI_MODE_WRONLY | MPI_MODE_CREATE,
1045 					MPI_INFO_NULL,
1046 					&(PhastaIOActiveFiles[i]->file_handle) );
1047 			if(rc)
1048 			{
1049 				*fileDescriptor = UNABLE_TO_OPEN_FILE;
1050 				return;
1051 			}
1052 		} // end of if "write"
1053 		free (fname);
1054 		free (imode);
1055 	} // end of if FileIndex != 0
1056 
1057 	endTimer(&timer_end);
1058 	printPerf("openfile", timer_start, timer_end, 0, 0, "");
1059 }
1060 
1061 /** close file for both POSIX and MPI-IO syncIO format.
1062  *
1063  * If it's old POSIX format, simply call posix fclose().
1064  *
1065  * If it's MPI-IO foramt:
1066  * in "read" mode, it simply close file with MPI-IO close routine.
1067  * in "write" mode, rank 0 in each group will re-assemble the master header and
1068  * offset table and write to the beginning of file, then close the file.
1069  */
1070 void closefile( int* fileDescriptor,
1071                  const char mode[] )
1072 {
1073 	double timer_start, timer_end;
1074 	startTimer(&timer_start);
1075 
1076 	int i = *fileDescriptor;
1077 	checkFileDescriptor("closefile",&i);
1078 
1079 	if ( PhastaIONextActiveIndex == 0 ) {
1080 		char* imode = StringStripper( mode );
1081 
1082 		if( cscompare( "write", imode )
1083 				|| cscompare( "append", imode ) ) {
1084 			fflush( fileArray[ *fileDescriptor - 1 ] );
1085 		}
1086 
1087 		fclose( fileArray[ *fileDescriptor - 1 ] );
1088 		free (imode);
1089      	}
1090 	else {
1091 		char* imode = StringStripper( mode );
1092 
1093 		//write master header here:
1094 		if ( cscompare( "write", imode ) ) {
1095 			//	      if ( PhastaIOActiveFiles[i]->nPPF * PhastaIOActiveFiles[i]->nFields < 2*ONE_MEGABYTE/8 )  //SHOULD BE CHECKED
1096 			//		MasterHeaderSize = 4*ONE_MEGABYTE;
1097 			//	      else
1098 			//		MasterHeaderSize = 4*ONE_MEGABYTE + PhastaIOActiveFiles[i]->nPPF * PhastaIOActiveFiles[i]->nFields * 8 - 2*ONE_MEGABYTE;
1099 
1100 			MasterHeaderSize = computeMHSize( PhastaIOActiveFiles[i]->nFields, PhastaIOActiveFiles[i]->nPPF, LATEST_WRITE_VERSION);
1101 			phprintf_0("Info closefile: myrank = %d, MasterHeaderSize = %d\n", PhastaIOActiveFiles[i]->myrank, MasterHeaderSize);
1102 
1103 			MPI_Status write_header_status;
1104 			char mpi_tag[MAX_FIELDS_NAME_LENGTH];
1105 			char version[MAX_FIELDS_NAME_LENGTH/4];
1106 			char mhsize[MAX_FIELDS_NAME_LENGTH/4];
1107 			int magic_number = ENDIAN_TEST_NUMBER;
1108 
1109 			if ( PhastaIOActiveFiles[i]->local_myrank == 0 )
1110 			{
1111 				bzero((void*)mpi_tag,MAX_FIELDS_NAME_LENGTH);
1112 				sprintf(mpi_tag, "MPI_IO_Tag : ");
1113 				memcpy(PhastaIOActiveFiles[i]->master_header,
1114 						mpi_tag,
1115 						MAX_FIELDS_NAME_LENGTH);
1116 
1117 				bzero((void*)version,MAX_FIELDS_NAME_LENGTH/4);
1118 				// this version is "1", print version in ASCII
1119 				sprintf(version, "version : %d",1);
1120 				memcpy(PhastaIOActiveFiles[i]->master_header + MAX_FIELDS_NAME_LENGTH/2,
1121 						version,
1122 						MAX_FIELDS_NAME_LENGTH/4);
1123 
1124 				// master header size is computed using the formula above
1125 				bzero((void*)mhsize,MAX_FIELDS_NAME_LENGTH/4);
1126 				sprintf(mhsize, "mhsize : ");
1127 				memcpy(PhastaIOActiveFiles[i]->master_header + MAX_FIELDS_NAME_LENGTH/4*3,
1128 						mhsize,
1129 						MAX_FIELDS_NAME_LENGTH/4);
1130 
1131 				bzero((void*)mpi_tag,MAX_FIELDS_NAME_LENGTH);
1132 				sprintf(mpi_tag,
1133 						"\nnFields : %d\n",
1134 						PhastaIOActiveFiles[i]->nFields);
1135 				memcpy(PhastaIOActiveFiles[i]->master_header+MAX_FIELDS_NAME_LENGTH,
1136 						mpi_tag,
1137 						MAX_FIELDS_NAME_LENGTH);
1138 
1139 				bzero((void*)mpi_tag,MAX_FIELDS_NAME_LENGTH);
1140 				sprintf(mpi_tag, "\nnPPF : %d\n", PhastaIOActiveFiles[i]->nPPF);
1141 				memcpy( PhastaIOActiveFiles[i]->master_header+
1142 						PhastaIOActiveFiles[i]->nFields *
1143 						MAX_FIELDS_NAME_LENGTH +
1144 						MAX_FIELDS_NAME_LENGTH * 2,
1145 						mpi_tag,
1146 						MAX_FIELDS_NAME_LENGTH);
1147 
1148 				memcpy( PhastaIOActiveFiles[i]->master_header+sizeof("MPI_IO_Tag : ")-1, //-1 sizeof returns the size of the string+1 for "\0"
1149 						&magic_number,
1150 						sizeof(int));
1151 
1152 				memcpy( PhastaIOActiveFiles[i]->master_header+sizeof("mhsize : ") -1 + MAX_FIELDS_NAME_LENGTH/4*3,
1153 						&MasterHeaderSize,
1154 						sizeof(int));
1155 			}
1156 
1157 			int j = 0;
1158 			unsigned long **header_table;
1159 			header_table = ( unsigned long ** )calloc(PhastaIOActiveFiles[i]->nFields, sizeof(unsigned long *));
1160 
1161 			for ( j = 0; j < PhastaIOActiveFiles[i]->nFields; j++ ) {
1162 				header_table[j]=( unsigned long * ) calloc( PhastaIOActiveFiles[i]->nPPF , sizeof( unsigned long));
1163 			}
1164 
1165 			//if( irank == 0 ) printf("gonna mpi_gather, myrank = %d\n", irank);
1166 			for ( j = 0; j < PhastaIOActiveFiles[i]->nFields; j++ ) {
1167 				MPI_Gather( PhastaIOActiveFiles[i]->my_offset_table[j],
1168 						PhastaIOActiveFiles[i]->nppp,
1169 						MPI_LONG_LONG_INT,
1170 						header_table[j],
1171 						PhastaIOActiveFiles[i]->nppp,
1172 						MPI_LONG_LONG_INT,
1173 						0,
1174 						PhastaIOActiveFiles[i]->local_comm );
1175 			}
1176 
1177 			if ( PhastaIOActiveFiles[i]->local_myrank == 0 ) {
1178 
1179 			        //if( irank == 0 ) printf("gonna memcpy for every procs, myrank = %d\n", irank);
1180 				for ( j = 0; j < PhastaIOActiveFiles[i]->nFields; j++ ) {
1181 					memcpy ( PhastaIOActiveFiles[i]->master_header +
1182 						VERSION_INFO_HEADER_SIZE +
1183 						j * PhastaIOActiveFiles[i]->nPPF * sizeof(unsigned long),
1184 						header_table[j],
1185 						PhastaIOActiveFiles[i]->nPPF * sizeof(unsigned long) );
1186 				}
1187 
1188 				//if( irank == 0 ) printf("gonna file_write_at(), myrank = %d\n", irank);
1189 				MPI_File_write_at( PhastaIOActiveFiles[i]->file_handle,
1190 						0,
1191 						PhastaIOActiveFiles[i]->master_header,
1192 						MasterHeaderSize,
1193 						MPI_CHAR,
1194 						&write_header_status );
1195 			}
1196 
1197 			////free(PhastaIOActiveFiles[i]->master_header);
1198 
1199 			for ( j = 0; j < PhastaIOActiveFiles[i]->nFields; j++ ) {
1200 				free ( header_table[j] );
1201 			}
1202 			free (header_table);
1203 		}
1204 
1205 		//if( irank == 0 ) printf("gonna file_close(), myrank = %d\n", irank);
1206 		MPI_File_close( &( PhastaIOActiveFiles[i]->file_handle ) );
1207 		free ( imode );
1208 	}
1209 
1210 	endTimer(&timer_end);
1211 	printPerf("closefile_", timer_start, timer_end, 0, 0, "");
1212 }
1213 
1214 int readHeader( FILE* f, const char phrase[],
1215     int* params, int numParams, const char* iotype) {
1216   isBinary(iotype);
1217   return readHeader(f,phrase,params,numParams);
1218 }
1219 
1220 void readheader( int* fileDescriptor,
1221                   const  char keyphrase[],
1222                   void* valueArray,
1223                   int*  nItems,
1224                   const char  datatype[],
1225                   const char  iotype[] )
1226 {
1227      	double timer_start, timer_end;
1228 
1229 	startTimer(&timer_start);
1230 
1231 	int i = *fileDescriptor;
1232 	checkFileDescriptor("readheader",&i);
1233 
1234 	if ( PhastaIONextActiveIndex == 0 ) {
1235 		int filePtr = *fileDescriptor - 1;
1236 		FILE* fileObject;
1237 		int* valueListInt;
1238 
1239 		if ( *fileDescriptor < 1 || *fileDescriptor > (int)fileArray.size() ) {
1240 			fprintf(stderr,"No file associated with Descriptor %d\n",*fileDescriptor);
1241 			fprintf(stderr,"openfile_ function has to be called before \n") ;
1242 			fprintf(stderr,"acessing the file\n ") ;
1243 			fprintf(stderr,"fatal error: cannot continue, returning out of call\n");
1244                         endTimer(&timer_end);
1245                         printPerf("readheader", timer_start, timer_end, 0, 0, "");
1246 			return;
1247 		}
1248 
1249                 LastHeaderKey[filePtr] = std::string(keyphrase);
1250 		LastHeaderNotFound = false;
1251 
1252 		fileObject = fileArray[ filePtr ] ;
1253 		Wrong_Endian = byte_order[ filePtr ];
1254 
1255 		isBinary( iotype );
1256 		typeSize( datatype );   //redundant call, just avoid a compiler warning.
1257 
1258 		// right now we are making the assumption that we will only write integers
1259 		// on the header line.
1260 
1261 		valueListInt = static_cast< int* >( valueArray );
1262 		int ierr = readHeader( fileObject ,
1263 				keyphrase,
1264 				valueListInt,
1265 				*nItems ) ;
1266 
1267 		byte_order[ filePtr ] = Wrong_Endian ;
1268 
1269 		if ( ierr ) LastHeaderNotFound = true;
1270 
1271 		//return ; // don't return, go to the end to print perf
1272 	}
1273 	else {
1274 		unsigned int skip_size;
1275 		int* valueListInt;
1276 		valueListInt = static_cast <int*>(valueArray);
1277 		char* token = NULL;
1278 		bool FOUND = false ;
1279 		isBinary( iotype );
1280 
1281 		MPI_Status read_offset_status;
1282 		char read_out_tag[MAX_FIELDS_NAME_LENGTH];
1283                 memset(read_out_tag, '\0', MAX_FIELDS_NAME_LENGTH);
1284 		char readouttag[MAX_FIELDS_NUMBER][MAX_FIELDS_NAME_LENGTH];
1285 		int j;
1286 
1287 		int string_length = strlen( keyphrase );
1288 		char* buffer = (char*) malloc ( string_length+1 );
1289 
1290 		strcpy ( buffer, keyphrase );
1291 		buffer[ string_length ] = '\0';
1292 
1293 		char* st2 = strtok ( buffer, "@" );
1294 		st2 = strtok (NULL, "@");
1295 		PhastaIOActiveFiles[i]->GPid = atoi(st2);
1296 		if ( char* p = strpbrk(buffer, "@") )
1297 			*p = '\0';
1298 
1299 		// Check if the user has input the right GPid
1300 		if ( ( PhastaIOActiveFiles[i]->GPid <=
1301 					PhastaIOActiveFiles[i]->myrank *
1302 					PhastaIOActiveFiles[i]->nppp )||
1303 				( PhastaIOActiveFiles[i]->GPid >
1304 					( PhastaIOActiveFiles[i]->myrank + 1 ) *
1305 					PhastaIOActiveFiles[i]->nppp ) )
1306 		{
1307 			*fileDescriptor = NOT_A_MPI_FILE;
1308 			printf("Error readheader: The file is not in syncIO new format, please check! myrank = %d, GPid = %d, nppp = %d, keyphrase = %s\n", PhastaIOActiveFiles[i]->myrank, PhastaIOActiveFiles[i]->GPid, PhastaIOActiveFiles[i]->nppp, keyphrase);
1309                         // It is possible atoi could not generate a clear integer from st2 because of additional garbage character in keyphrase
1310                         endTimer(&timer_end);
1311                         printPerf("readheader", timer_start, timer_end, 0, 0, "");
1312 			return;
1313 		}
1314 
1315 		// Find the field we want ...
1316 		//for ( j = 0; j<MAX_FIELDS_NUMBER; j++ )
1317 		for ( j = 0; j<PhastaIOActiveFiles[i]->nFields; j++ )
1318 		{
1319 			memcpy( readouttag[j],
1320 					PhastaIOActiveFiles[i]->master_header + j*MAX_FIELDS_NAME_LENGTH+MAX_FIELDS_NAME_LENGTH*2+1,
1321 					MAX_FIELDS_NAME_LENGTH-1 );
1322 		}
1323 
1324 		for ( j = 0; j<PhastaIOActiveFiles[i]->nFields; j++ )
1325 		{
1326 			token = strtok ( readouttag[j], ":" );
1327 
1328 			//if ( cscompare( buffer, token ) )
1329 			if ( cscompare( token , buffer ) && cscompare( buffer, token ) )
1330                         // This double comparison is required for the field "number of nodes" and all the other fields that start with "number of nodes" (i.g. number of nodes in the mesh").
1331                         // Would be safer to rename "number of nodes" by "number of nodes in the part" so that the name are completely unique. But much more work to do that (Nspre, phParAdapt, etc).
1332                         // Since the field name are unique in SyncIO (as it includes part ID), this should be safe and there should be no issue with the "?" trailing character.
1333 			{
1334 				PhastaIOActiveFiles[i]->read_field_count = j;
1335 				FOUND = true;
1336                                 //printf("buffer: %s | token: %s | j: %d\n",buffer,token,j);
1337 				break;
1338 			}
1339 		}
1340                 free(buffer);
1341 
1342 		if (!FOUND)
1343 		{
1344 			//if(irank==0) printf("Warning readheader: Not found %s \n",keyphrase); //PhastaIOActiveFiles[i]->myrank is certainly initialized here.
1345 			if(PhastaIOActiveFiles[i]->myrank == 0) printf("WARNING readheader: Not found %s\n",keyphrase);
1346                         endTimer(&timer_end);
1347                         printPerf("readheader", timer_start, timer_end, 0, 0, "");
1348 			return;
1349 		}
1350 
1351 		// Find the part we want ...
1352 		PhastaIOActiveFiles[i]->read_part_count = PhastaIOActiveFiles[i]->GPid -
1353 			PhastaIOActiveFiles[i]->myrank * PhastaIOActiveFiles[i]->nppp - 1;
1354 
1355 		PhastaIOActiveFiles[i]->my_offset =
1356 			PhastaIOActiveFiles[i]->my_read_table[PhastaIOActiveFiles[i]->read_field_count][PhastaIOActiveFiles[i]->read_part_count];
1357 
1358 		// 	  printf("****Rank %d offset is %d\n",PhastaIOActiveFiles[i]->myrank,PhastaIOActiveFiles[i]->my_offset);
1359 
1360 		// Read each datablock header here ...
1361 
1362 		MPI_File_read_at_all( PhastaIOActiveFiles[i]->file_handle,
1363 				PhastaIOActiveFiles[i]->my_offset+1,
1364 				read_out_tag,
1365 				MAX_FIELDS_NAME_LENGTH-1,
1366 				MPI_CHAR,
1367 				&read_offset_status );
1368 		token = strtok ( read_out_tag, ":" );
1369 
1370 		// 	  printf("&&&&Rank %d read_out_tag is %s\n",PhastaIOActiveFiles[i]->myrank,read_out_tag);
1371 
1372 		if( cscompare( keyphrase , token ) ) //No need to compare also token with keyphrase like above. We should already have the right one. Otherwise there is a problem.
1373 		{
1374 			FOUND = true ;
1375 			token = strtok( NULL, " ,;<>" );
1376 			skip_size = atoi( token );
1377 			for( j=0; j < *nItems && ( token = strtok( NULL," ,;<>") ); j++ )
1378 				valueListInt[j] = atoi( token );
1379 
1380 			if ( j < *nItems )
1381 			{
1382 				fprintf( stderr, "Expected # of ints not found for: %s\n", keyphrase );
1383 			}
1384 		}
1385                 else {
1386                   //if(irank==0)
1387 		  if(PhastaIOActiveFiles[i]->myrank == 0)
1388                   // If we enter this if, there is a problem with the name of some fields
1389                   {
1390                     printf("Error readheader: Unexpected mismatch between keyphrase = %s and token = %s\n",keyphrase,token);
1391                   }
1392                 }
1393 	}
1394 
1395 	endTimer(&timer_end);
1396 	printPerf("readheader", timer_start, timer_end, 0, 0, "");
1397 
1398 }
1399 
1400 void readDataBlock(
1401     FILE* fileObject,
1402     void* valueArray,
1403     int nItems,
1404     const char  datatype[],
1405     const char  iotype[] )
1406 {
1407   isBinary(iotype);
1408   size_t type_size = typeSize( datatype );
1409   if ( binary_format ) {
1410     char junk = '\0';
1411     fread( valueArray, type_size, nItems, fileObject );
1412     fread( &junk, sizeof(char), 1 , fileObject );
1413     if ( Wrong_Endian ) SwapArrayByteOrder( valueArray, type_size, nItems );
1414   } else {
1415     char* ts1 = StringStripper( datatype );
1416     if ( cscompare( "integer", ts1 ) ) {
1417       for( int n=0; n < nItems ; n++ )
1418         fscanf(fileObject, "%d\n",(int*)((int*)valueArray+n) );
1419     } else if ( cscompare( "double", ts1 ) ) {
1420       for( int n=0; n < nItems ; n++ )
1421         fscanf(fileObject, "%lf\n",(double*)((double*)valueArray+n) );
1422     }
1423     free (ts1);
1424   }
1425 }
1426 
1427 void readdatablock( int*  fileDescriptor,
1428                      const char keyphrase[],
1429                      void* valueArray,
1430                      int*  nItems,
1431                      const char  datatype[],
1432                      const char  iotype[] )
1433 {
1434 	//if(irank == 0) printf("entering readdatablock()\n");
1435 	unsigned long data_size = 0;
1436 	double timer_start, timer_end;
1437 	startTimer(&timer_start);
1438 
1439 	int i = *fileDescriptor;
1440 	checkFileDescriptor("readdatablock",&i);
1441 
1442 	if ( PhastaIONextActiveIndex == 0 ) {
1443 		int filePtr = *fileDescriptor - 1;
1444 		FILE* fileObject;
1445 		char junk;
1446 
1447 		if ( *fileDescriptor < 1 || *fileDescriptor > (int)fileArray.size() ) {
1448 			fprintf(stderr,"No file associated with Descriptor %d\n",*fileDescriptor);
1449 			fprintf(stderr,"openfile_ function has to be called before\n") ;
1450 			fprintf(stderr,"acessing the file\n ") ;
1451 			fprintf(stderr,"fatal error: cannot continue, returning out of call\n");
1452                         endTimer(&timer_end);
1453                         printPerf("readdatablock", timer_start, timer_end, 0, 0, "");
1454 			return;
1455 		}
1456 
1457 		// error check..
1458 		// since we require that a consistant header always preceed the data block
1459 		// let us check to see that it is actually the case.
1460 
1461 		if ( ! cscompare( LastHeaderKey[ filePtr ].c_str(), keyphrase ) ) {
1462 			fprintf(stderr, "Header not consistant with data block\n");
1463 			fprintf(stderr, "Header: %s\n", LastHeaderKey[ filePtr ].c_str() );
1464 			fprintf(stderr, "DataBlock: %s\n ", keyphrase );
1465 			fprintf(stderr, "Please recheck read sequence \n");
1466 			if( Strict_Error ) {
1467 				fprintf(stderr, "fatal error: cannot continue, returning out of call\n");
1468                                 endTimer(&timer_end);
1469                                 printPerf("readdatablock", timer_start, timer_end, 0, 0, "");
1470 				return;
1471 			}
1472 		}
1473 
1474 		if ( LastHeaderNotFound ) {
1475                         endTimer(&timer_end);
1476                         printPerf("readdatablock", timer_start, timer_end, 0, 0, "");
1477                         return;
1478                 }
1479 		fileObject = fileArray[ filePtr ];
1480 		Wrong_Endian = byte_order[ filePtr ];
1481                 LastHeaderKey.erase(filePtr);
1482                 readDataBlock(fileObject,valueArray,*nItems,datatype,iotype);
1483 
1484 		//return;
1485 	}
1486 	else {
1487 		// 	  printf("read data block\n");
1488 		MPI_Status read_data_status;
1489 		size_t type_size = typeSize( datatype );
1490 		int nUnits = *nItems;
1491 		isBinary( iotype );
1492 
1493 		// read datablock then
1494 		//MR CHANGE
1495 		//             if ( cscompare ( datatype, "double"))
1496 		char* ts2 = StringStripper( datatype );
1497 		if ( cscompare ( "double" , ts2))
1498 			//MR CHANGE END
1499 		{
1500 
1501 			MPI_File_read_at_all_begin( PhastaIOActiveFiles[i]->file_handle,
1502 					PhastaIOActiveFiles[i]->my_offset + DB_HEADER_SIZE,
1503 					valueArray,
1504 					nUnits,
1505 					MPI_DOUBLE );
1506 			MPI_File_read_at_all_end( PhastaIOActiveFiles[i]->file_handle,
1507 					valueArray,
1508 					&read_data_status );
1509 			data_size=8*nUnits;
1510 
1511 		}
1512 		//MR CHANGE
1513 		//             else if ( cscompare ( datatype, "integer"))
1514 		else if ( cscompare ( "integer" , ts2))
1515 			//MR CHANGE END
1516 		{
1517 			MPI_File_read_at_all_begin(PhastaIOActiveFiles[i]->file_handle,
1518 					PhastaIOActiveFiles[i]->my_offset + DB_HEADER_SIZE,
1519 					valueArray,
1520 					nUnits,
1521 					MPI_INT );
1522 			MPI_File_read_at_all_end( PhastaIOActiveFiles[i]->file_handle,
1523 					valueArray,
1524 					&read_data_status );
1525 			data_size=4*nUnits;
1526 		}
1527 		else
1528 		{
1529 			*fileDescriptor = DATA_TYPE_ILLEGAL;
1530 			printf("readdatablock - DATA_TYPE_ILLEGAL - %s\n",datatype);
1531                         endTimer(&timer_end);
1532                         printPerf("readdatablock", timer_start, timer_end, 0, 0, "");
1533 			return;
1534 		}
1535                 free(ts2);
1536 
1537 
1538 		// 	  printf("%d Read finishe\n",PhastaIOActiveFiles[i]->myrank);
1539 
1540 		// Swap data byte order if endianess is different ...
1541 		if ( PhastaIOActiveFiles[i]->Wrong_Endian )
1542 		{
1543 			SwapArrayByteOrder( valueArray, type_size, nUnits );
1544 		}
1545 	}
1546 
1547 	endTimer(&timer_end);
1548 	char extra_msg[1024];
1549 	memset(extra_msg, '\0', 1024);
1550 	char* key = StringStripper(keyphrase);
1551 	sprintf(extra_msg, " field is %s ", key);
1552 	printPerf("readdatablock", timer_start, timer_end, data_size, 1, extra_msg);
1553         free(key);
1554 
1555 }
1556 
1557 void writeHeader(FILE* f,
1558                  const char keyphrase[],
1559                  const void* valueArray,
1560                  const int nItems,
1561                  const int ndataItems,
1562                  const char datatype[],
1563                  const char iotype[])
1564 {
1565   isBinary( iotype );
1566 
1567   const int _newline =
1568     ( ndataItems > 0 ) ? sizeof( char ) : 0;
1569   int size_of_nextblock =
1570     ( binary_format ) ? typeSize(datatype) * ndataItems + _newline : ndataItems;
1571 
1572   fprintf( f, "%s : < %d > ", keyphrase, size_of_nextblock );
1573   for( int i = 0; i < nItems; i++ )
1574     fprintf( f, "%d ", *((int*)((int*)valueArray+i)));
1575   fprintf( f, "\n");
1576 }
1577 
1578 void writeheader(  const int* fileDescriptor,
1579                     const char keyphrase[],
1580                     const void* valueArray,
1581                     const int* nItems,
1582                     const int* ndataItems,
1583                     const char datatype[],
1584                     const char iotype[])
1585 {
1586 
1587 	//if(irank == 0) printf("entering writeheader()\n");
1588 
1589 	double timer_start, timer_end;
1590 	startTimer(&timer_start);
1591 
1592 	int i = *fileDescriptor;
1593 	checkFileDescriptor("writeheader",&i);
1594 
1595 	if ( PhastaIONextActiveIndex == 0 ) {
1596 		int filePtr = *fileDescriptor - 1;
1597 		if ( *fileDescriptor < 1 || *fileDescriptor > (int)fileArray.size() ) {
1598 			fprintf(stderr,"No file associated with Descriptor %d\n",*fileDescriptor);
1599 			fprintf(stderr,"openfile_ function has to be called before \n") ;
1600 			fprintf(stderr,"acessing the file\n ") ;
1601 			fprintf(stderr,"fatal error: cannot continue, returning out of call\n");
1602                         endTimer(&timer_end);
1603                         printPerf("writeheader", timer_start, timer_end, 0, 0, "");
1604 			return;
1605 		}
1606 
1607                 LastHeaderKey[filePtr] = std::string(keyphrase);
1608 		DataSize = *ndataItems;
1609                 FILE* fileObject = fileArray[ filePtr ] ;
1610                 header_type[ filePtr ] = typeSize( datatype );
1611                 writeHeader(fileObject,keyphrase,valueArray,*nItems,
1612                     *ndataItems,datatype,iotype);
1613 	}
1614 	else { // else it's parallel I/O
1615 		DataSize = *ndataItems;
1616 		size_t type_size = typeSize( datatype );
1617 		isBinary( iotype );
1618 		char mpi_tag[MAX_FIELDS_NAME_LENGTH];
1619 
1620 		int string_length = strlen( keyphrase );
1621 		char* buffer = (char*) malloc ( string_length+1 );
1622 
1623 		strcpy ( buffer, keyphrase);
1624 		buffer[ string_length ] = '\0';
1625 
1626 		char* st2 = strtok ( buffer, "@" );
1627 		st2 = strtok (NULL, "@");
1628 		PhastaIOActiveFiles[i]->GPid = atoi(st2);
1629 
1630 		if ( char* p = strpbrk(buffer, "@") )
1631 			*p = '\0';
1632 
1633 	        bzero((void*)mpi_tag,MAX_FIELDS_NAME_LENGTH);
1634 		sprintf(mpi_tag, "\n%s : %d\n", buffer, PhastaIOActiveFiles[i]->field_count);
1635 		unsigned long offset_value;
1636 
1637 		int temp = *ndataItems;
1638 		unsigned long number_of_items = (unsigned long)temp;
1639 		MPI_Barrier(PhastaIOActiveFiles[i]->local_comm);
1640 
1641 		MPI_Scan( &number_of_items,
1642 				&offset_value,
1643 				1,
1644 				MPI_LONG_LONG_INT,
1645 				MPI_SUM,
1646 				PhastaIOActiveFiles[i]->local_comm );
1647 
1648 		offset_value = (offset_value - number_of_items) * type_size;
1649 
1650 		offset_value += PhastaIOActiveFiles[i]->local_myrank *
1651 			DB_HEADER_SIZE +
1652 			PhastaIOActiveFiles[i]->next_start_address;
1653 		// This offset is the starting address of each datablock header...
1654 		PhastaIOActiveFiles[i]->my_offset = offset_value;
1655 
1656 		// Write in my offset table ...
1657 		PhastaIOActiveFiles[i]->my_offset_table[PhastaIOActiveFiles[i]->field_count][PhastaIOActiveFiles[i]->part_count] =
1658 			PhastaIOActiveFiles[i]->my_offset;
1659 
1660 		// Update the next-start-address ...
1661 		PhastaIOActiveFiles[i]->next_start_address = offset_value +
1662 			number_of_items * type_size +
1663 			DB_HEADER_SIZE;
1664 		MPI_Bcast( &(PhastaIOActiveFiles[i]->next_start_address),
1665 				1,
1666 				MPI_LONG_LONG_INT,
1667 				PhastaIOActiveFiles[i]->local_numprocs-1,
1668 				PhastaIOActiveFiles[i]->local_comm );
1669 
1670 		// Prepare datablock header ...
1671 		int _newline = (*ndataItems>0)?sizeof(char):0;
1672 		unsigned int size_of_nextblock = type_size * (*ndataItems) + _newline;
1673 
1674 		//char datablock_header[255];
1675 		//bzero((void*)datablock_header,255);
1676 		char datablock_header[DB_HEADER_SIZE];
1677 		bzero((void*)datablock_header,DB_HEADER_SIZE);
1678 
1679 		PhastaIOActiveFiles[i]->GPid = PhastaIOActiveFiles[i]->nppp*PhastaIOActiveFiles[i]->myrank+PhastaIOActiveFiles[i]->part_count;
1680 		sprintf( datablock_header,
1681 				"\n%s : < %u >",
1682 				keyphrase,
1683 				size_of_nextblock );
1684 
1685 		for ( int j = 0; j < *nItems; j++ )
1686 		{
1687 			sprintf( datablock_header,
1688 					"%s %d ",
1689 					datablock_header,
1690 					*((int*)((int*)valueArray+j)));
1691 		}
1692 		sprintf( datablock_header,
1693 				"%s\n ",
1694 				datablock_header );
1695 
1696 		// Write datablock header ...
1697 		//MR CHANGE
1698 		// 	if ( cscompare(datatype,"double") )
1699 		char* ts1 = StringStripper( datatype );
1700 		if ( cscompare("double",ts1) )
1701 			//MR CHANGE END
1702 		{
1703 			free ( PhastaIOActiveFiles[i]->double_chunk );
1704 			PhastaIOActiveFiles[i]->double_chunk = ( double * )malloc( (sizeof( double )*number_of_items+ DB_HEADER_SIZE));
1705 
1706 			double * aa = ( double * )datablock_header;
1707 			memcpy(PhastaIOActiveFiles[i]->double_chunk, aa, DB_HEADER_SIZE);
1708 		}
1709 		//MR CHANGE
1710 		// 	if  ( cscompare(datatype,"integer") )
1711 		else if ( cscompare("integer",ts1) )
1712 			//MR CHANGE END
1713 		{
1714 			free ( PhastaIOActiveFiles[i]->int_chunk );
1715 			PhastaIOActiveFiles[i]->int_chunk = ( int * )malloc( (sizeof( int )*number_of_items+ DB_HEADER_SIZE));
1716 
1717 			int * aa = ( int * )datablock_header;
1718 			memcpy(PhastaIOActiveFiles[i]->int_chunk, aa, DB_HEADER_SIZE);
1719 		}
1720 		else {
1721 			//             *fileDescriptor = DATA_TYPE_ILLEGAL;
1722 			printf("writeheader - DATA_TYPE_ILLEGAL - %s\n",datatype);
1723                         endTimer(&timer_end);
1724                         printPerf("writeheader", timer_start, timer_end, 0, 0, "");
1725 			return;
1726 		}
1727 		free(ts1);
1728 
1729 		PhastaIOActiveFiles[i]->part_count++;
1730 		if ( PhastaIOActiveFiles[i]->part_count == PhastaIOActiveFiles[i]->nppp ) {
1731 			//A new field will be written
1732 			if ( PhastaIOActiveFiles[i]->local_myrank == 0 ) {
1733 				memcpy( PhastaIOActiveFiles[i]->master_header +
1734 						PhastaIOActiveFiles[i]->field_count *
1735 						MAX_FIELDS_NAME_LENGTH +
1736 						MAX_FIELDS_NAME_LENGTH * 2,
1737 						mpi_tag,
1738 						MAX_FIELDS_NAME_LENGTH-1);
1739 			}
1740 			PhastaIOActiveFiles[i]->field_count++;
1741 			PhastaIOActiveFiles[i]->part_count=0;
1742 		}
1743                 free(buffer);
1744 	}
1745 
1746 	endTimer(&timer_end);
1747 	printPerf("writeheader", timer_start, timer_end, 0, 0, "");
1748 }
1749 
1750 void writeDataBlock( FILE* f,
1751                      const void* valueArray,
1752                      const int   nItems,
1753                      const char  datatype[],
1754                      const char  iotype[]  )
1755 {
1756   isBinary( iotype );
1757   size_t type_size = typeSize( datatype );
1758   if ( binary_format ) {
1759     fwrite( valueArray, type_size, nItems, f );
1760     fprintf( f,"\n");
1761   } else {
1762     char* ts1 = StringStripper( datatype );
1763     if ( cscompare( "integer", ts1 ) ) {
1764       for( int n=0; n < nItems ; n++ )
1765         fprintf(f,"%d\n",*((int*)((int*)valueArray+n)));
1766     } else if ( cscompare( "double", ts1 ) ) {
1767       for( int n=0; n < nItems ; n++ )
1768         fprintf(f,"%lf\n",*((double*)((double*)valueArray+n)));
1769     }
1770     free (ts1);
1771   }
1772 }
1773 
1774 void writedatablock( const int* fileDescriptor,
1775                       const char keyphrase[],
1776                       const void* valueArray,
1777                       const int* nItems,
1778                       const char datatype[],
1779                       const char iotype[] )
1780 {
1781 	//if(irank == 0) printf("entering writedatablock()\n");
1782 
1783 	unsigned long data_size = 0;
1784 	double timer_start, timer_end;
1785 	startTimer(&timer_start);
1786 
1787 	int i = *fileDescriptor;
1788 	checkFileDescriptor("writedatablock",&i);
1789 
1790 	if ( PhastaIONextActiveIndex == 0 ) {
1791 		int filePtr = *fileDescriptor - 1;
1792 
1793 		if ( *fileDescriptor < 1 || *fileDescriptor > (int)fileArray.size() ) {
1794 			fprintf(stderr,"No file associated with Descriptor %d\n",*fileDescriptor);
1795 			fprintf(stderr,"openfile_ function has to be called before \n") ;
1796 			fprintf(stderr,"acessing the file\n ") ;
1797 			fprintf(stderr,"fatal error: cannot continue, returning out of call\n");
1798                         endTimer(&timer_end);
1799                         printPerf("writedatablock", timer_start, timer_end, 0, 0, "");
1800 			return;
1801 		}
1802 		// since we require that a consistant header always preceed the data block
1803 		// let us check to see that it is actually the case.
1804 
1805 		if ( ! cscompare( LastHeaderKey[ filePtr ].c_str(), keyphrase ) ) {
1806 			fprintf(stderr, "Header not consistant with data block\n");
1807 			fprintf(stderr, "Header: %s\n", LastHeaderKey[ filePtr ].c_str() );
1808 			fprintf(stderr, "DataBlock: %s\n ", keyphrase );
1809 			fprintf(stderr, "Please recheck write sequence \n");
1810 			if( Strict_Error ) {
1811 				fprintf(stderr, "fatal error: cannot continue, returning out of call\n");
1812                                 endTimer(&timer_end);
1813                                 printPerf("writedatablock", timer_start, timer_end, 0, 0, "");
1814 				return;
1815 			}
1816 		}
1817 
1818 		FILE* fileObject =  fileArray[ filePtr ] ;
1819 		size_t type_size=typeSize( datatype );
1820 		isBinary( iotype );
1821 
1822                 LastHeaderKey.erase(filePtr);
1823 
1824 		if ( header_type[filePtr] != (int)type_size ) {
1825 			fprintf(stderr,"header and datablock differ on typeof data in the block for\n");
1826 			fprintf(stderr,"keyphrase : %s\n", keyphrase);
1827 			if( Strict_Error ) {
1828 				fprintf(stderr,"fatal error: cannot continue, returning out of call\n" );
1829                                 endTimer(&timer_end);
1830                                 printPerf("writedatablock", timer_start, timer_end, 0, 0, "");
1831 				return;
1832 			}
1833 		}
1834 
1835 		int nUnits = *nItems;
1836 
1837 		if ( nUnits != DataSize ) {
1838 			fprintf(stderr,"header and datablock differ on number of data items for\n");
1839 			fprintf(stderr,"keyphrase : %s\n", keyphrase);
1840 			if( Strict_Error ) {
1841 				fprintf(stderr,"fatal error: cannot continue, returning out of call\n" );
1842                                 endTimer(&timer_end);
1843                                 printPerf("writedatablock", timer_start, timer_end, 0, 0, "");
1844 				return;
1845 			}
1846 		}
1847                 writeDataBlock(fileObject,valueArray,*nItems,datatype,iotype);
1848 	}
1849 	else {  // syncIO case
1850 		MPI_Status write_data_status;
1851 		isBinary( iotype );
1852 		int nUnits = *nItems;
1853 
1854 		//MR CHANGE
1855 		// 	if ( cscompare(datatype,"double") )
1856 		char* ts1 = StringStripper( datatype );
1857 		if ( cscompare("double",ts1) )
1858 			//MR CHANGE END
1859 		{
1860 			memcpy((PhastaIOActiveFiles[i]->double_chunk+DB_HEADER_SIZE/sizeof(double)), valueArray, nUnits*sizeof(double));
1861 			MPI_File_write_at_all_begin( PhastaIOActiveFiles[i]->file_handle,
1862 					PhastaIOActiveFiles[i]->my_offset,
1863 					PhastaIOActiveFiles[i]->double_chunk,
1864 					//BLOCK_SIZE/sizeof(double),
1865 					nUnits+DB_HEADER_SIZE/sizeof(double),
1866 					MPI_DOUBLE );
1867 			MPI_File_write_at_all_end( PhastaIOActiveFiles[i]->file_handle,
1868 					PhastaIOActiveFiles[i]->double_chunk,
1869 					&write_data_status );
1870 			data_size=8*nUnits;
1871 		}
1872 		//MR CHANGE
1873 		// 	else if ( cscompare ( datatype, "integer"))
1874 		else if ( cscompare("integer",ts1) )
1875 			//MR CHANGE END
1876 		{
1877 			memcpy((PhastaIOActiveFiles[i]->int_chunk+DB_HEADER_SIZE/sizeof(int)), valueArray, nUnits*sizeof(int));
1878 			MPI_File_write_at_all_begin( PhastaIOActiveFiles[i]->file_handle,
1879 					PhastaIOActiveFiles[i]->my_offset,
1880 					PhastaIOActiveFiles[i]->int_chunk,
1881 					nUnits+DB_HEADER_SIZE/sizeof(int),
1882 					MPI_INT );
1883 			MPI_File_write_at_all_end( PhastaIOActiveFiles[i]->file_handle,
1884 					PhastaIOActiveFiles[i]->int_chunk,
1885 					&write_data_status );
1886 			data_size=4*nUnits;
1887 		}
1888 		else {
1889 			printf("Error: writedatablock - DATA_TYPE_ILLEGAL - %s\n",datatype);
1890                         endTimer(&timer_end);
1891                         printPerf("writedatablock", timer_start, timer_end, 0, 0, "");
1892 			return;
1893 		}
1894                 free(ts1);
1895 	}
1896 
1897 	endTimer(&timer_end);
1898 	char extra_msg[1024];
1899 	memset(extra_msg, '\0', 1024);
1900 	char* key = StringStripper(keyphrase);
1901 	sprintf(extra_msg, " field is %s ", key);
1902 	printPerf("writedatablock", timer_start, timer_end, data_size, 1, extra_msg);
1903         free(key);
1904 
1905 }
1906 
1907 void
1908 SwapArrayByteOrder( void* array,
1909                      int   nbytes,
1910                      int   nItems )
1911 {
1912 	/* This swaps the byte order for the array of nItems each
1913 		 of size nbytes , This will be called only locally  */
1914 	int i,j;
1915 	unsigned char* ucDst = (unsigned char*)array;
1916 
1917 	for(i=0; i < nItems; i++) {
1918 		for(j=0; j < (nbytes/2); j++)
1919 			std::swap( ucDst[j] , ucDst[(nbytes - 1) - j] );
1920 		ucDst += nbytes;
1921 	}
1922 }
1923 
1924 void
1925 writestring( int* fileDescriptor,
1926               const char inString[] )
1927 {
1928 
1929 	int filePtr = *fileDescriptor - 1;
1930 	FILE* fileObject = fileArray[filePtr] ;
1931 	fprintf(fileObject,"%s",inString );
1932 	return;
1933 }
1934 
1935 void
1936 Gather_Headers( int* fileDescriptor,
1937                 vector< string >& headers )
1938 {
1939 
1940 	FILE* fileObject;
1941 	char Line[1024];
1942 
1943 	fileObject = fileArray[ (*fileDescriptor)-1 ];
1944 
1945 	while( !feof(fileObject) ) {
1946 		fgets( Line, 1024, fileObject);
1947 		if ( Line[0] == '#' ) {
1948 			headers.push_back( Line );
1949 		} else {
1950 			break;
1951 		}
1952 	}
1953 	rewind( fileObject );
1954 	clearerr( fileObject );
1955 }
1956 void
1957 isWrong( void ) { (Wrong_Endian) ? fprintf(stdout,"YES\n"): fprintf(stdout,"NO\n") ; }
1958 
1959 void
1960 togglestrictmode( void ) { Strict_Error = !Strict_Error; }
1961 
1962 int
1963 isLittleEndian( void )
1964 {
1965 	// this function returns a 1 if the current running architecture is
1966 	// LittleEndian Byte Ordered, else it returns a zero
1967 
1968 	union  {
1969 		long a;
1970 		char c[sizeof( long )];
1971 	} endianUnion;
1972 
1973 	endianUnion.a = 1 ;
1974 
1975 	if ( endianUnion.c[sizeof(long)-1] != 1 ) return 1 ;
1976 	else return 0;
1977 }
1978 
1979 namespace PHASTA {
1980 	const char* const PhastaIO_traits<int>::type_string = "integer";
1981 	const char* const PhastaIO_traits<double>::type_string = "double";
1982 }
1983