1 #include <stdio.h> 2 #include <iostream> 3 #include <string.h> 4 #include <stdlib.h> 5 //#define OMPI_SKIP_MPICXX 1 //Added in the CMakeList.txt file 6 #include <mpi.h> 7 #include <math.h> 8 #include <sys/stat.h> 9 #include <sys/types.h> 10 #include "phastaIO.h" 11 12 inline int 13 cscompare( const char teststring[], 14 const char targetstring[] ) 15 { 16 char* s1 = const_cast<char*>(teststring); 17 char* s2 = const_cast<char*>(targetstring); 18 19 while( *s1 == ' ') s1++; 20 while( *s2 == ' ') s2++; 21 while( ( *s1 ) 22 && ( *s2 ) 23 && ( *s2 != '?') 24 && ( tolower( *s1 )==tolower( *s2 ) ) ) { 25 s1++; 26 s2++; 27 while( *s1 == ' ') s1++; 28 while( *s2 == ' ') s2++; 29 } 30 if ( !( *s1 ) || ( *s1 == '?') ) return 1; 31 else return 0; 32 } 33 34 int main(int argc, char *argv[]) { 35 36 MPI_Init(&argc,&argv); 37 38 int myrank, N_procs; 39 MPI_Comm_rank(MPI_COMM_WORLD, &myrank); 40 MPI_Comm_size(MPI_COMM_WORLD, &N_procs); 41 42 FILE * pFile; 43 char target[1024], pool[256]; 44 char * temp, * token; 45 int i, j, k, N_restart_integer, N_restart_double; 46 int N_geombc_double, N_geombc_integer; 47 int N_steps, N_parts, N_files; 48 49 pFile = fopen("./IO.N2O.input","r"); 50 if (pFile == NULL) 51 printf("Error openning\n"); 52 53 fgets( target, 1024, pFile ); 54 token = strtok ( target, ";" );strcpy(pool,token); 55 temp = strtok ( pool, ":" );temp = strtok ( NULL, ":" ); 56 N_geombc_double = atoi(temp); 57 58 fgets( target, 1024, pFile ); 59 token = strtok ( target, ";" );strcpy(pool,token); 60 temp = strtok ( pool, ":" );temp = strtok ( NULL, ":" ); 61 N_geombc_integer = atoi(temp); 62 63 fgets( target, 1024, pFile ); 64 token = strtok ( target, ";" );strcpy(pool,token); 65 temp = strtok ( pool, ":" );temp = strtok ( NULL, ":" ); 66 N_restart_double = atoi(temp); 67 68 fgets( target, 1024, pFile ); 69 token = strtok ( target, ";" );strcpy(pool,token); 70 temp = strtok ( pool, ":" );temp = strtok ( NULL, ":" ); 71 N_restart_integer = atoi(temp); 72 73 fgets( target, 1024, pFile ); 74 token = strtok ( target, ";" );strcpy(pool,token); 75 temp = strtok ( pool, ":" );temp = strtok ( NULL, ":" ); 76 N_steps = atoi(temp); 77 78 fgets( target, 1024, pFile ); 79 token = strtok ( target, ";" );strcpy(pool,token); 80 temp = strtok ( pool, ":" );temp = strtok ( NULL, ":" ); 81 N_parts = atoi(temp); 82 83 if(myrank==0){ 84 printf("numpe is %d and start is %d\n",N_parts,N_steps); 85 } 86 87 fgets( target, 1024, pFile ); 88 token = strtok ( target, ";" );strcpy(pool,token); 89 temp = strtok ( pool, ":" );temp = strtok ( NULL, ":" ); 90 N_files = atoi(temp); 91 92 double ***Dfield; int ***Ifield; 93 int ***paraD, ***paraI, *expectD, *expectI; 94 char **fieldNameD, **fileTypeD, **dataTypeD, **headerTypeD; 95 char **fieldNameI, **fileTypeI, **dataTypeI, **headerTypeI; 96 97 int WriteLockD[N_restart_double]; 98 int WriteLockI[N_restart_integer]; 99 100 int nppp = N_parts/N_procs; 101 int startpart = myrank * nppp +1; 102 int endpart = startpart + nppp - 1; 103 char gfname[64], numTemp[128]; 104 int iarray[10], igeom; 105 long isize; 106 107 108 if (N_parts != N_procs) { 109 printf("Input error: number of parts should be equal to the number of procs!\n"); 110 printf("Please modify the IO.O2N.input file!\n"); 111 return 0; 112 } 113 114 115 116 ///////////////////// reading /////////////////////////////// 117 118 int nppf = N_parts/N_files; 119 int N_geombc = N_geombc_double + N_geombc_integer; 120 int readHandle, GPID; 121 char fname[255],fieldtag[255]; 122 123 int irestart; 124 125 Dfield = new double**[N_restart_double]; 126 Ifield = new int**[N_restart_integer]; 127 128 paraD = new int**[N_restart_double]; 129 paraI = new int**[N_restart_integer]; 130 131 expectD = new int[N_restart_double]; 132 expectI = new int[N_restart_integer]; 133 134 fieldNameD = new char*[N_restart_double]; 135 fileTypeD = new char*[N_restart_double]; 136 dataTypeD = new char*[N_restart_double]; 137 headerTypeD = new char*[N_restart_double]; 138 139 fieldNameI = new char*[N_restart_integer]; 140 fileTypeI = new char*[N_restart_integer]; 141 dataTypeI = new char*[N_restart_integer]; 142 headerTypeI = new char*[N_restart_integer]; 143 144 if (N_restart_double>0) 145 for ( i = 0; i < N_restart_double; i++ ) 146 { 147 WriteLockD[i]=0; 148 Dfield[i] = new double*[nppp]; 149 150 paraD[i] = new int*[nppp]; 151 152 fieldNameD[i] = new char[128]; 153 fileTypeD[i] = new char[128]; 154 dataTypeD[i] = new char[128]; 155 headerTypeD[i] = new char[128]; 156 } 157 158 if (N_restart_integer>0) 159 for ( i = 0; i < N_restart_integer; i++ ) 160 { 161 WriteLockI[i]=0; 162 Ifield[i] = new int*[nppp]; 163 164 paraI[i] = new int*[nppp]; 165 166 fieldNameI[i] = new char[128]; 167 fileTypeI[i] = new char[128]; 168 dataTypeI[i] = new char[128]; 169 headerTypeI[i] = new char[128]; 170 } 171 172 173 ///////////////////// reading /////////////////////////////// 174 175 int N_restart = N_restart_double + N_restart_integer; 176 int readHandle1; 177 178 bzero((void*)fname,255); 179 sprintf(fname,"./%d-procs_case/restart-dat.%d.%d",N_parts,N_steps,((int)(myrank/(N_procs/N_files))+1)); 180 181 //if(myrank==0){ 182 // printf("Myrank is %d - Filename is %s \n",myrank,fname); 183 //} 184 185 int nfields; 186 queryphmpiio(fname, &nfields, &nppf); 187 //initphmpiio(&N_restart, &nppf, &N_files,&readHandle1, "write") ;//WRONG 188 initphmpiio(&nfields, &nppf, &N_files, &readHandle1, "read"); 189 openfile(fname, "read", &readHandle1); 190 191 for ( i = 0; i < N_restart_double; i++ ) 192 { 193 fgets( target, 1024, pFile ); 194 temp = strtok( target, ";" ); 195 token = strtok( temp, "," ); 196 strcpy( fileTypeD[i], token ); 197 token = strtok ( NULL, "," ); 198 strcpy( fieldNameD[i], token ); 199 token = strtok ( NULL, "," ); 200 strcpy( dataTypeD[i], token ); 201 token = strtok ( NULL, "," ); 202 strcpy( headerTypeD[i], token ); 203 token = strtok ( NULL, "," ); 204 strcpy( numTemp, token ); 205 expectD[i] = atoi (numTemp); 206 207 for ( j = 0; j < nppp; j++ ) 208 { 209 paraD[i][j] = new int[expectD[i]]; 210 211 for ( k = 0; k < 10; k++ ) 212 iarray[k]=0; 213 214 GPID = startpart + j; 215 bzero((void*)fieldtag,255); 216 sprintf(fieldtag,"%s@%d",fieldNameD[i],GPID); 217 218 //printf("myrank %d - filedtag %s\n",myrank,fieldtag); 219 220 iarray[0]=-1; 221 readheader( &readHandle1, 222 fieldtag, 223 (void*)iarray, 224 &expectD[i], 225 "double", 226 "binary" ); 227 228 if ( iarray[0]==-1 ) 229 WriteLockD[i]=1; 230 if ( WriteLockD[i]==0 ) 231 { 232 for ( k = 0; k < expectD[i]; k++ ) 233 paraD[i][j][k] = iarray[k]; 234 235 if ( cscompare("block",headerTypeD[i]) ) 236 { 237 if ( expectD[i]==1) 238 isize = paraD[i][j][0]; 239 else 240 isize = paraD[i][j][0] * paraD[i][j][1]; 241 242 Dfield[i][j] = new double[isize]; 243 readdatablock( &readHandle1, 244 fieldtag, 245 (void*)Dfield[i][j], 246 &isize, 247 "double", 248 "binary" ); 249 } 250 251 } 252 } 253 } 254 255 for ( i = 0; i < N_restart_integer; i++ ) 256 { 257 fgets( target, 1024, pFile ); 258 temp = strtok( target, ";" ); 259 token = strtok( temp, "," ); 260 strcpy( fileTypeI[i], token ); 261 token = strtok ( NULL, "," ); 262 strcpy( fieldNameI[i], token ); 263 token = strtok ( NULL, "," ); 264 strcpy( dataTypeI[i], token ); 265 token = strtok ( NULL, "," ); 266 strcpy( headerTypeI[i], token ); 267 token = strtok ( NULL, "," ); 268 strcpy( numTemp, token ); 269 expectI[i] = atoi (numTemp); 270 271 for ( j = 0; j < nppp; j++ ) 272 { 273 paraI[i][j] = new int[expectI[i]]; 274 275 for ( k = 0; k < 10; k++ ) 276 iarray[k]=0; 277 278 GPID = startpart + j; 279 bzero((void*)fieldtag,255); 280 sprintf(fieldtag,"%s@%d",fieldNameI[i],GPID); 281 iarray[0]=-1; 282 283 //printf("Rank %d, fieldname is %s \n",myrank,fieldtag); 284 285 readheader( &readHandle1, 286 fieldtag, 287 (void*)iarray, 288 &expectI[i], 289 "integer", 290 "binary" ); 291 292 if ( iarray[0]==-1) 293 WriteLockI[i]=1; 294 if ( WriteLockI[i]==0 ) 295 { 296 for ( k = 0; k < expectI[i]; k++ ) 297 paraI[i][j][k] = iarray[k]; 298 299 if ( cscompare("block",headerTypeI[i]) ) 300 { 301 if ( expectI[i]==1) 302 isize = paraI[i][j][0]; 303 else 304 isize = paraI[i][j][0] * paraI[i][j][1]; 305 306 Ifield[i][j] = new int[isize]; 307 readdatablock( &readHandle1, 308 fieldtag, 309 (void*)Ifield[i][j], 310 &isize, 311 "integer", 312 "binary" ); 313 } 314 } 315 } 316 317 } 318 319 closefile(&readHandle1, "write"); 320 finalizephmpiio(&readHandle1); 321 322 //////////////////////////writing//////////////////////////// 323 324 int irstou; 325 int magic_number = 362436; 326 int* mptr = &magic_number; 327 int nitems = 1; 328 329 //MR CHANGE 330 bzero((void*)fname,255); 331 sprintf(fname,"./%d-procs_case-1PPP",N_parts); 332 if(0<mkdir(fname,0777)) { printf("ERROR - Could not create procs_case-1PPP directory\n"); return 1; } 333 //MR CHANGE END 334 335 bzero((void*)fname,255); 336 sprintf(fname,"./%d-procs_case-1PPP/restart.%d.%d",N_parts,N_steps,myrank+1); 337 openfile(fname,"write", &irstou); 338 339 /* writing the top ascii header for the restart file */ 340 341 writestring( &irstou,"# PHASTA Input File Version 2.0\n"); 342 writestring( &irstou, 343 "# format \"keyphrase : sizeofnextblock usual headers\"\n"); 344 345 bzero( (void*)fname, 255 ); 346 writestring( &irstou, fname ); 347 348 writestring( &irstou, fname ); 349 writestring( &irstou,"\n"); 350 351 352 isize = 1; 353 nitems = 1; 354 iarray[ 0 ] = 1; 355 writeheader( &irstou, "byteorder magic number ", 356 (void*)iarray, &nitems, &isize, "integer", "binary" ); 357 358 nitems = 1; 359 writedatablock( &irstou, "byteorder magic number ", 360 (void*)mptr, &nitems, "integer", "binary" ); 361 362 for ( i = 0; i < N_restart_double; i++ ) 363 { 364 for ( j = 0; j < nppp; j++ ) 365 { 366 if ( WriteLockD[i] == 0 ) 367 { 368 if ( cscompare("header",headerTypeD[i]) ) 369 { 370 bzero( (void*)fname, 255 ); 371 sprintf(fname,"%s : < 0 > %d\n", fieldNameD[i],paraD[i][j][0]); 372 writestring( &irstou, fname ); 373 } 374 375 if ( cscompare("block",headerTypeD[i]) ) 376 { 377 if ( expectD[i]==1 ) 378 isize = paraD[i][j][0]; 379 else 380 isize = paraD[i][j][0] * paraD[i][j][1]; 381 382 for ( k = 0; k < expectD[i]; k++ ) 383 iarray[k] = paraD[i][j][k]; 384 385 if ( cscompare("header",headerTypeD[i]) ) 386 isize = 0; 387 388 writeheader( &irstou, 389 fieldNameD[i], 390 (void*)iarray, 391 &expectD[i], 392 &isize, 393 "double", 394 "binary"); 395 writedatablock( &irstou, 396 fieldNameD[i], 397 (void*)Dfield[i][j], 398 &isize, 399 "double", 400 "binary"); 401 } 402 403 404 if ( cscompare("block",headerTypeD[i]) ) 405 delete [] Dfield[i][j]; 406 } 407 delete [] paraD[i][j]; 408 } 409 } 410 411 for ( i = 0; i < N_restart_integer; i++ ) 412 { 413 for ( j = 0; j < nppp; j++ ) 414 { 415 416 if ( WriteLockI[i] == 0 ) 417 { 418 419 if ( cscompare("header",headerTypeI[i]) ) 420 { 421 bzero( (void*)fname, 255 ); 422 sprintf(fname,"%s : < 0 > %d\n", fieldNameI[i],paraI[i][j][0]); 423 writestring( &irstou, fname ); 424 } 425 426 if ( cscompare("block",headerTypeI[i]) ) 427 { 428 if ( expectI[i]==1 ) 429 isize = paraI[i][j][0]; 430 else 431 isize = paraI[i][j][0] * paraI[i][j][1]; 432 433 for ( k = 0; k < expectI[i]; k++ ) 434 iarray[k] = paraI[i][j][k]; 435 436 writeheader( &irstou, 437 fieldNameI[i], 438 (void*)iarray, 439 &expectI[i], 440 &isize, 441 "integer", 442 "binary"); 443 writedatablock( &irstou, 444 fieldNameI[i], 445 (void*)Ifield[i][j], 446 &isize, 447 "integer", 448 "binary"); 449 } 450 451 if ( cscompare("block",headerTypeI[i]) ) 452 delete [] Ifield[i][j]; 453 } 454 delete [] paraI[i][j]; 455 } 456 } 457 458 459 closefile( &irstou, "write" ); 460 MPI_Barrier(MPI_COMM_WORLD); 461 462 463 if (N_restart_double>0) 464 for ( i = 0; i < N_restart_double; i++ ) 465 { 466 delete [] Dfield[i]; 467 delete [] paraD[i]; 468 469 delete [] fieldNameD[i]; 470 delete [] fileTypeD[i]; 471 delete [] dataTypeD[i]; 472 delete [] headerTypeD[i]; 473 } 474 475 if (N_restart_integer>0) 476 for ( i = 0; i < N_restart_integer; i++ ) 477 { 478 delete [] Ifield[i]; 479 delete [] paraI[i]; 480 481 delete [] fieldNameI[i]; 482 delete [] fileTypeI[i]; 483 delete [] dataTypeI[i]; 484 delete [] headerTypeI[i]; 485 } 486 487 delete [] Dfield; 488 delete [] Ifield; 489 490 delete [] paraD; 491 delete [] paraI; 492 493 delete [] expectD; 494 delete [] expectI; 495 496 delete [] fieldNameD; 497 delete [] fileTypeD; 498 delete [] dataTypeD; 499 delete [] headerTypeD; 500 501 delete [] fieldNameI; 502 delete [] fileTypeI; 503 delete [] dataTypeI; 504 delete [] headerTypeI; 505 506 fclose(pFile); 507 508 if (myrank==0) 509 { 510 printf("\nFinished transfer, please check data using:\n"); 511 printf(" grep -a ': <' filename \n\n"); 512 } 513 514 MPI_Finalize(); 515 516 } 517 518 519