1c readnblk.f (pronounce "Reed and Block Dot Eff") contains: 2c 3c module readarrays ("Red Arrays") -- contains the arrays that 4c are read in from binary files but not immediately blocked 5c through pointers. 6c 7c subroutine readnblk ("Reed and Block") -- allocates space for 8c and reads data to be contained in module readarrays. Reads 9c all remaining data and blocks them with pointers. 10c 11 12 13 module readarrays 14 15 real*8, allocatable :: point2x(:,:) 16 real*8, allocatable :: qold(:,:) 17 real*8, allocatable :: uold(:,:) 18 real*8, allocatable :: acold(:,:) 19 integer, allocatable :: iBCtmp(:) 20 real*8, allocatable :: BCinp(:,:) 21 22 integer, allocatable :: point2ilwork(:) 23 integer, allocatable :: nBC(:) 24 integer, allocatable :: point2iper(:) 25 integer, allocatable :: point2ifath(:) 26 integer, allocatable :: point2nsons(:) 27 28 end module 29 30 31 32 subroutine readnblk 33c 34 use readarrays 35 include "common.h" 36c 37 real*8, allocatable :: xread(:,:), qread(:,:), acread(:,:) 38 real*8, allocatable :: uread(:,:) 39 real*8, allocatable :: BCinpread(:,:) 40 integer, allocatable :: iperread(:), iBCtmpread(:) 41 integer, allocatable :: ilworkread(:), nBCread(:) 42 character*10 cname2 43 character*8 mach2 44!MR CHANGE 45! character*20 fmt1 46 character*30 fmt1 47!MR CHANGE END 48 character*255 fname1,fnamer,fnamelr 49 character*255 warning 50 integer igeomBAK, ibndc, irstin, ierr 51 integer intfromfile(50) ! integers read from headers 52 53cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc 54ccccccccccccccccccccccccccccccccccccccccc New PhastaIO Definition Part ccccccccccccccccccccccccccccccccccccccccc 55 56 integer :: descriptor, descriptorG, GPID, color, nfiles, nfields 57 integer :: numparts, nppf 58 integer :: ierr_io, numprocs, itmp, itmp2 59 character*255 fname2, temp2 60 character*64 temp1 61 62cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc 63cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc 64 65 66c 67c.... determine the step number to start with 68c 69 open(unit=72,file='numstart.dat',status='old') 70 read(72,*) irstart 71 close(72) 72 lstep=irstart ! in case restart files have no fields 73 74c 75 fname1='geombc.dat' 76 fname1= trim(fname1) // cname2(myrank+1) 77c fnamelr='restart.latest' 78 79 itmp=1 80 if (irstart .gt. 0) itmp = int(log10(float(irstart)))+1 81 write (fmt1,"('(''restart.'',i',i1,',1x)')") itmp 82 write (fnamer,fmt1) irstart 83 fnamer = trim(fnamer) // cname2(myrank+1) 84c fnamelr = trim(fnamelr) // cname2(myrank+1) 85 86c 87c.... open input files 88c 89 90 91c call openfile( fname1, 'read?', igeomBAK ); 92 93 94c 95c.... try opening restart.latest.proc before trying restart.stepno.proc 96c 97c call openfile( fnamelr, 'read?', irstin ); 98c if ( irstin .eq. 0 ) 99 100!MR CHANGE 101! call openfile( fnamer, 'read?', irstin ); 102!MR CHANGE END 103 104! either one will work 105c 106c.... input the geometry parameters 107c 108 109cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc 110!MR CHANGE 111 112! nfiles = 2 113! nfields = 31 114! numparts = 8 115! nppp = numparts/numpe 116! nppf = numparts/nfiles 117 118 nfiles = nsynciofiles 119! nfields = nsynciofieldsreadgeombc 120 numparts = numpe !This is the common settings. Beware if you try to compute several parts per process 121 122! nppp = numparts/numpe 123! nppf = numparts/nfiles 124!MR CHANGE END 125 126c fnamer="/home/nliu/develop/PIG4/512-procs_case/geombc-dat" 127 128 color = int(myrank/(numparts/nfiles)) !Should call the color routine in SyncIO here 129 itmp2 = int(log10(float(color+1)))+1 130 write (temp2,"('(''geombc-dat.'',i',i1,')')") itmp2 131 write (fnamer,temp2) (color+1) 132 fnamer = trim(fnamer) // char(0) 133 134c fnamer="/home/nliu/develop/test-case/512-procs_case/geombc-dat" 135 136 itwo=2 137 ione=1 138 ieleven=11 139 itmp = int(log10(float(myrank+1)))+1 140 141 call queryphmpiio(fnamer, nfields, nppf); 142 if (myrank == 0) then 143 write(*,*) 'Number of fields in geombc-dat: ',nfields 144 write(*,*) 'Number of parts per file geombc-dat: ',nppf 145 endif 146 call initphmpiio( nfields, nppf, nfiles, igeom, 147 & 'read' // char(0)) 148 call openfile( fnamer, 'read' // char(0), igeom ) 149 150 call readheader(igeom,'number of nodes' // char(0),numnp,ione, 151 & 'integer' // char(0), iotype) 152 153 call readheader(igeom,'number of modes' // char(0),nshg,ione, 154 & 'integer' // char(0), iotype) 155 156 call readheader(igeom,'number of interior elements' // char(0), 157 & numel,ione,'integer' // char(0), iotype) 158 159 call readheader(igeom,'number of boundary elements' // char(0), 160 & numelb,ione,'integer' // char(0),iotype) 161 162 call readheader(igeom,'maximum number of element nodes' // char(0), 163 & nen,ione,'integer' // char(0),iotype) 164 165 call readheader(igeom,'number of interior tpblocks' // char(0), 166 & nelblk,ione,'integer' // char(0) ,iotype) 167 168 call readheader(igeom,'number of boundary tpblocks' // char(0), 169 & nelblb,ione,'integer' // char(0), iotype) 170 171 call readheader(igeom,'number of nodes with Dirichlet BCs' // char(0), 172 & numpbc,ione,'integer' // char(0),iotype) 173 174 call readheader(igeom,'number of shape functions' // char(0), 175 & ntopsh,ione,'integer' // char(0),iotype) 176 177c call closefile( igeom, "read" ) 178c call finalizephmpiio( igeom ) 179 180! if(myrank==0) then 181! print *, "Reading Finished, ********* : ", trim(fname2) 182! endif 183 184 185ccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc 186 187c ieleven=11 188c ione=1 189c fname1='number of nodes?' 190c call readheader(igeomBAK,fname1,numnp,ione,'integer', iotype) 191c fname1='number of modes?' 192c call readheader(igeomBAK,fname1,nshg,ione,'integer', iotype) 193cc fname1='number of global modes?' 194cc call readheader(igeomBAK,fname1,nshgt,ione,'integer', iotype) 195c fname1='number of interior elements?' 196c call readheader(igeomBAK,fname1,numel,ione,'integer', iotype) 197c fname1='number of boundary elements?' 198c call readheader(igeomBAK,fname1,numelb,ione,'integer', iotype) 199c fname1='maximum number of element nodes?' 200c call readheader(igeomBAK,fname1,nen,ione,'integer', iotype) 201c fname1='number of interior tpblocks?' 202c call readheader(igeomBAK,fname1,nelblk,ione,'integer', iotype) 203c fname1='number of boundary tpblocks?' 204c call readheader(igeomBAK,fname1,nelblb,ione,'integer', iotype) 205c fname1='number of nodes with Dirichlet BCs?' 206c call readheader(igeomBAK,fname1,numpbc,ione,'integer', iotype) 207c fname1='number of shape functions?' 208c call readheader(igeomBAK,fname1,ntopsh,ione,'integer', iotype) 209 210c 211c.... calculate the maximum number of boundary element nodes 212c 213 nenb = 0 214 do i = 1, melCat 215 if (nen .eq. nenCat(i,nsd)) nenb = max(nenCat(i,nsd-1), nenb) 216 enddo 217c 218 if (myrank == master) then 219 if (nenb .eq. 0) call error ('input ','nen ',nen) 220 endif 221c 222c.... setup some useful constants 223c 224 I3nsd = nsd / 3 ! nsd=3 integer flag 225 E3nsd = float(I3nsd) ! nsd=3 real flag 226c 227 if(matflg(1,1).lt.0) then 228 nflow = nsd + 1 229 else 230 nflow = nsd + 2 231 endif 232 ndof = nsd + 2 233 nsclr=impl(1)/100 234 ndof=ndof+nsclr ! number of sclr transport equations to solve 235 236 ndofBC = ndof + I3nsd ! dimension of BC array 237 ndiBCB = 2 ! dimension of iBCB array 238 ndBCB = ndof + 1 ! dimension of BCB array 239c 240 nsymdf = (ndof*(ndof + 1)) / 2 ! symm. d.o.f.'s 241c 242c.... ----------------------> Communication tasks <-------------------- 243c 244 245cc still read in new 246 247 if(numpe > 1) then 248 249cc MR CHANGE 250 call readheader(igeom,'size of ilwork array' // char(0), 251 & nlwork,ione,'integer' // char(0) ,iotype) 252 253 call readheader(igeom,'ilwork' //char(0) ,nlwork,ione, 254 & 'integer' // char(0) , iotype) 255 256 allocate( point2ilwork(nlwork) ) 257 allocate( ilworkread(nlwork) ) 258 call readdatablock(igeom,'ilwork' // char(0),ilworkread, 259 & nlwork,'integer' // char(0) , iotype) 260 261c call closefile( igeom, "read" ) 262c call finalizephmpiio( igeom ) 263 264c fname1='size of ilwork array?' 265c call readheader(igeomBAK,fname1,nlwork,ione,'integer', iotype) 266 267c ione=1 268c fname1='ilwork?' 269c call readheader(igeomBAK,fname1,nlwork,ione,'integer', iotype) 270 271c allocate( point2ilwork(nlwork) ) 272c allocate( ilworkread(nlwork) ) 273c call readdatablock(igeomBAK,fname1,ilworkread, 274c & nlwork,'integer', iotype) 275cc MR CHANGE 276 277 278 point2ilwork = ilworkread 279 call ctypes (point2ilwork) 280 else 281 nlwork=1 282 allocate( point2ilwork(1)) 283 nshg0 = nshg 284 endif 285 286cccccccccccccccccccccccccccccccccccccccccc 287 288 itwo=2 289 290c print *, "fname2 is", fname2 291 292cc MR CHANGE 293c call initphmpiio(nfields,nppf,nfiles,igeom,'read') 294c call openfile( fnamer, 'read', igeom ) 295CC MR CHANGE 296 297 call readheader(igeom,'co-ordinates' // char(0),intfromfile,itwo, 298 & 'double' // char(0), iotype) 299 numnp=intfromfile(1) 300c print *, "read out @@@@@@ is ", numnp 301 allocate( point2x(numnp,nsd) ) 302 allocate( xread(numnp,nsd) ) 303 ixsiz=numnp*nsd 304 call readdatablock(igeom,'co-ordinates' // char(0),xread,ixsiz, 305 & 'double' // char(0), iotype) 306 point2x = xread 307 308! call closefile( igeom, "read" ) 309! call finalizephmpiio( igeom ) 310 311! print *, "Finished finalize" 312 313c deallocate (point2x) 314c deallocate (xread) 315 316cccccccccccccccccccccccccccccccccccccccccc 317 318c 319c.... read the node coordinates 320c 321 322c itwo=2 323c fname1='co-ordinates?' 324c call readheader(igeomBAK,fname1,intfromfile,itwo, 'double', iotype) 325c numnp=intfromfile(1) 326cc nsd=intfromfile(2) 327c allocate( point2x(numnp,nsd) ) 328c allocate( xread(numnp,nsd) ) 329c ixsiz=numnp*nsd 330c call readdatablock(igeomBAK,fname1,xread,ixsiz, 'double',iotype) 331c point2x = xread 332 333c 334c.... read in and block out the connectivity 335c 336 337! !MR CHANGE 338! This is not necessary but this avoids to have the geombc files opend two times. 339! A better way consists in pasisng the file handler to genblk or make it global or use igeomBAK instead of igeom 340! call closefile( igeom, "read" ) 341! call finalizephmpiio( igeom ) 342! !MR CHANGE END 343 344 call genblk (IBKSIZ) 345 346! !MR CHANGE 347! call initphmpiio( nfields, nppf, nfiles, igeom, 'read') 348! call openfile( fnamer, 'read', igeom ) 349! !MR CHANGE END 350 351c 352c.... read the boundary condition mapping array 353c 354 355cc MR CHANGE 356! call initphmpiio(nfields,nppf,nfiles,igeom, 'read') 357! call openfile( fnamer, 'read', igeom ) 358cc MR CHANGE 359 360 ione=1 361 call readheader(igeom,'bc mapping array' // char(0),nshg,ione, 362 & 'integer' // char(0),iotype) 363 364c fname1='bc mapping array?' 365c call readheader(igeomBAK,fname1,nshg, 366c & ione,'integer', iotype) 367 368 allocate( nBC(nshg) ) 369 370 allocate( nBCread(nshg) ) 371 372c call readdatablock(igeomBAK,fname1,nBCread,nshg,'integer',iotype) 373 call readdatablock(igeom,'bc mapping array' // char(0), 374 & nBCread,nshg,'integer' // char(0),iotype) 375 376 nBC=nBCread 377 378c 379c.... read the temporary iBC array 380c 381 ione=1 382 call readheader(igeom,'bc codes array' // char(0) ,numpbc,ione, 383 & 'integer' // char(0),iotype) 384 385c ione = 1 386c fname1='bc codes array?' 387c call readheader(igeomBAK,fname1,numpbc, 388c & ione, 'integer', iotype) 389 390!MR CHANGE 391! if ( numpbc > 0 ) then 392! allocate( iBCtmp(numpbc) ) 393! allocate( iBCtmpread(numpbc) ) 394! c call readdatablock(igeomBAK,fname1,iBCtmpread,numpbc, 395! c & 'integer',iotype) 396! call readdatablock(igeom,fname2,iBCtmpread,numpbc, 397! & 'integer',iotype) 398! iBCtmp=iBCtmpread 399! else ! sometimes a partition has no BC's 400! allocate( iBCtmp(1) ) 401! iBCtmp=0 402! endif 403 404 if ( numpbc > 0 ) then 405 allocate( iBCtmp(numpbc) ) 406 allocate( iBCtmpread(numpbc) ) 407 else 408 allocate( iBCtmp(1) ) 409 allocate( iBCtmpread(1) ) 410 endif 411c call readdatablock(igeomBAK,fname1,iBCtmpread,numpbc, 412c & 'integer',iotype) 413 call readdatablock(igeom,'bc codes array' // char(0), 414 & iBCtmpread,numpbc,'integer' // char(0),iotype) 415 416 if ( numpbc > 0 ) then 417 iBCtmp=iBCtmpread 418 else ! sometimes a partition has no BC's 419 deallocate( iBCtmpread) 420 iBCtmp=0 421 endif 422!MR CHANGE END 423 424c 425c.... read boundary condition data 426c 427 428 ione=1 429 430c ione=1 431c fname1='boundary condition array?' 432c call readheader(igeomBAK,fname1,intfromfile, 433c & ione, 'double', iotype) 434 call readheader(igeom,'boundary condition array' // char(0), 435 & intfromfile,ione, 'double' // char(0), iotype) 436c here intfromfile(1) contains (ndof+7)*numpbc 437!MR CHANGE 438! if ( numpbc > 0 ) then 439! if(intfromfile(1).ne.(ndof+7)*numpbc) then 440! warning='WARNING more data in BCinp than needed: keeping 1st' 441! write(*,*) warning, ndof+7 442! endif 443! allocate( BCinp(numpbc,ndof+7) ) 444! nsecondrank=intfromfile(1)/numpbc 445! allocate( BCinpread(numpbc,nsecondrank) ) 446! iBCinpsiz=intfromfile(1) 447! c call readdatablock(igeomBAK,fname1,BCinpread,iBCinpsiz, 448! c & 'double',iotype) 449! call readdatablock(igeom,fname2,BCinpread,iBCinpsiz, 450! & 'double',iotype) 451! BCinp(:,1:(ndof+7))=BCinpread(:,1:(ndof+7)) 452! else ! sometimes a partition has no BC's 453! allocate( BCinp(1,ndof+7) ) 454! BCinp=0 455! endif 456 457 if ( numpbc > 0 ) then 458! if(intfromfile(1).ne.(ndof+7)*numpbc) then 459! warning='WARNING more data in BCinp than needed: keeping 1st' 460! write(*,*) warning, ndof+7 461! endif 462 allocate( BCinp(numpbc,ndof+7) ) 463 nsecondrank=intfromfile(1)/numpbc 464 allocate( BCinpread(numpbc,nsecondrank) ) 465 iBCinpsiz=intfromfile(1) 466 else 467 allocate( BCinp(1,ndof+7) ) 468 allocate( BCinpread(0,0) ) !dummy 469 iBCinpsiz=intfromfile(1) 470 endif 471c call readdatablock(igeomBAK,fname1,BCinpread,iBCinpsiz, 472c & 'double',iotype) 473 474 call readdatablock(igeom,'boundary condition array' // char(0), 475 & BCinpread,iBCinpsiz,'double' // char(0) ,iotype) 476 477 if ( numpbc > 0 ) then 478 BCinp(:,1:(ndof+7))=BCinpread(:,1:(ndof+7)) 479 else ! sometimes a partition has no BC's 480 deallocate(BCinpread) 481 BCinp=0 482 endif 483!MR CHANGE END 484 485c 486c.... read periodic boundary conditions 487c 488 489 ione=1 490c fname1='periodic masters array?' 491c call readheader(igeomBAK,fname1,nshg, 492c & ione, 'integer', iotype) 493 call readheader(igeom,'periodic masters array' // char(0) ,nshg, 494 & ione,'integer' // char(0), iotype) 495 allocate( point2iper(nshg) ) 496 allocate( iperread(nshg) ) 497c call readdatablock(igeomBAK,fname1,iperread,nshg, 498c & 'integer',iotype) 499 call readdatablock(igeom,'periodic masters array' // char(0), 500 & iperread,nshg,'integer' // char(0),iotype) 501 point2iper=iperread 502 503 504! !MR CHANGE 505! call closefile( igeom, "read" ) 506! call finalizephmpiio( igeom ) 507! !MR CHANGE END 508 509c 510c.... generate the boundary element blocks 511c 512 call genbkb (ibksiz) 513 514 515! !MR CHANGE 516! write(*,*) 'HELLO 12 from ', myrank 517! !MR CHANGE END 518 519c 520c Read in the nsons and ifath arrays if needed 521c 522c There is a fundamental shift in the meaning of ifath based on whether 523c there exist homogenous directions in the flow. 524c 525c HOMOGENOUS DIRECTIONS EXIST: Here nfath is the number of inhomogenous 526c points in the TOTAL mesh. That is to say that each partition keeps a 527c link to ALL inhomogenous points. This link is furthermore not to the 528c sms numbering but to the original structured grid numbering. These 529c inhomogenous points are thought of as fathers, with their sons being all 530c the points in the homogenous directions that have this father's 531c inhomogeneity. The array ifath takes as an arguement the sms numbering 532c and returns as a result the father. 533c 534c In this case nsons is the number of sons that each father has and ifath 535c is an array which tells the 536c 537c NO HOMOGENOUS DIRECTIONS. In this case the mesh would grow to rapidly 538c if we followed the above strategy since every partition would index its 539c points to the ENTIRE mesh. Furthermore, there would never be a need 540c to average to a node off processor since there is no spatial averaging. 541c Therefore, to properly account for this case we must recognize it and 542c inerrupt certain actions (i.e. assembly of the average across partitions). 543c This case is easily identified by noting that maxval(nsons) =1 (i.e. no 544c father has any sons). Reiterating to be clear, in this case ifath does 545c not point to a global numbering but instead just points to itself. 546c 547 nfath=1 ! some architectures choke on a zero or undeclared 548 ! dimension variable. This sets it to a safe, small value. 549 if(((iLES .lt. 20) .and. (iLES.gt.0)) 550 & .or. (itwmod.gt.0) ) then ! don't forget same 551 ! conditional in proces.f 552 553c read (igeomBAK) nfath ! nfath already read in input.f, 554 ! needed for alloc 555 ione=1 556c call creadlist(igeomBAK,ione,nfath) 557c fname1='keyword sonfath?' 558 if(nohomog.gt.0) then 559 560 561! fname1='number of father-nodes?' 562! call readheader(igeomBAK,fname1,nfath,ione,'integer', iotype) 563 564 call readheader(igeom,'number of father-nodes' // char(0), 565 & nfath,ione,'integer' // char(0), iotype) 566 567c 568c fname1='keyword nsons?' 569! fname1='number of son-nodes for each father?' 570! call readheader(igeomBAK,fname1,nfath,ione,'integer', iotype) 571 572 call readheader(igeom, 573 & 'number of son-nodes for each father' // char(0), 574 & nfath,ione,'integer' // char(0), iotype) 575 576 allocate (point2nsons(nfath)) 577 578! call readdatablock(igeomBAK,fname1,point2nsons,nfath, 579! & 'integer',iotype) 580 call readdatablock(igeom, 581 & 'number of son-nodes for each father' // char(0), 582 & point2nsons,nfath,'integer' // char(0), iotype) 583 584c 585! fname1='keyword ifath?' 586! call readheader(igeomBAK,fname1,nshg,ione,'integer', iotype) 587 588 call readheader(igeom,'keyword ifath' // char(0),nshg,ione, 589 & 'integer' // char(0), iotype) 590 591 allocate (point2ifath(nshg)) 592 593! call readdatablock(igeomBAK,fname1,point2ifath,nshg, 594! & 'integer',iotype) 595 call readdatablock(igeom,'keyword ifath' // char(0),point2ifath, 596 & nshg,'integer' // char(0) , iotype) 597 598c 599 nsonmax=maxval(point2nsons) 600c 601 else ! this is the case where there is no homogeneity 602 ! therefore ever node is a father (too itself). sonfath 603 ! (a routine in NSpre) will set this up but this gives 604 ! you an option to avoid that. 605 nfath=nshg 606 allocate (point2nsons(nfath)) 607 point2nsons=1 608 allocate (point2ifath(nshg)) 609 do i=1,nshg 610 point2ifath(i)=i 611 enddo 612 nsonmax=1 613c 614 endif 615 else 616 allocate (point2nsons(1)) 617 allocate (point2ifath(1)) 618 endif 619 620! !MR CHANGE 621 call closefile( igeom, "read" // char(0) ) 622 call finalizephmpiio( igeom ) 623! !MR CHANGE END 624 625! !MR CHANGE 626! write(*,*) 'HELLO 13 from ', myrank 627! !MR CHANGE END 628 629c 630c renumber the master partition for SPEBC 631c 632c if((myrank.eq.master).and.(irscale.ge.0)) then 633c call setSPEBC(numnp, nfath, nsonmax) 634c call renum(point2x,point2ifath,point2nsons) 635c endif 636c 637c.... Read restart files 638 639c$$$ccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc 640c 641c nfields = 11 642c numparts = 512 643c nppp = numparts/numpe 644c startpart = myrank * nppp +1 645c int endpart = startpart + nppp - 1 646c nppf = numparts/nfiles 647cc fnamer = "/users/nliu/PIG4/4-procs_case/restart-file" 648cc fnamer="./4-procs_case/restart-file" 649 650 itmp=1 651 if (irstart .gt. 0) itmp = int(log10(float(irstart+1)))+1 652 653 write (fmt1,"('(''restart-dat.'',i',i1,',1x)')") itmp 654 655 write (fnamer,fmt1) irstart 656 fnamer = trim(fnamer) // cname2(color+1) 657 658 call queryphmpiio(fnamer // char(0), nfields, nppf); 659 if (myrank == 0) then 660 write(*,*) 'Number of fields in restart-dat: ',nfields 661 write(*,*) 'Number of parts per file restart-dat: ',nppf 662 endif 663 call initphmpiio(nfields,nppf,nfiles,descriptor, 664 & 'read' // char(0)) 665 call openfile( fnamer // char(0) , 666 & 'read' // char(0), descriptor ) 667 668 ithree=3 669c call creadlist(irstin,ithree,nshg2,ndof2,lstep) 670 671 itmp = int(log10(float(myrank+1)))+1 672 673c print *, "Solution is : ", fname1 674 675 intfromfile=0 676 call readheader(descriptor,'solution' // char(0) ,intfromfile, 677 & ithree,'integer' // char(0), iotype) 678c 679c.... read the values of primitive variables into q 680c 681 682c print *, "intfromfile(1) is ", intfromfile(1) 683c print *, "intfromfile(2) is ", intfromfile(2) 684c print *, "intfromfile(3) is ", intfromfile(3) 685 686 allocate( qold(nshg,ndof) ) 687 if(intfromfile(1).ne.0) then 688 nshg2=intfromfile(1) 689 ndof2=intfromfile(2) 690 lstep=intfromfile(3) 691 if(ndof2.ne.ndof) then 692 693 endif 694c 695 if (nshg2 .ne. nshg) 696 & call error ('restar ', 'nshg ', nshg) 697 allocate( qread(nshg,ndof2) ) 698 iqsiz=nshg*ndof2 699 call readdatablock(descriptor,fname1 // char(0),qread,iqsiz, 700 & 'double' // char(0),iotype) 701 qold(:,1:ndof)=qread(:,1:ndof) 702 deallocate(qread) 703 else 704 if (myrank.eq.master) then 705 if (matflg(1,1).eq.0) then ! compressible 706 warning='Solution is set to zero (with p and T to one)' 707 else 708 warning='Solution is set to zero' 709 endif 710 write(*,*) warning 711 endif 712 qold=zero 713 if (matflg(1,1).eq.0) then ! compressible 714 qold(:,1)=one ! avoid zero pressure 715 qold(:,nflow)=one ! avoid zero temperature 716 endif 717 endif 718 719 720! !MR CHANGE 721! write(*,*) 'HELLO 16-8 from ', myrank 722! !MR CHANGE END 723 724c itmp=1 725c if (myrank .gt. 0) itmp = int(log10(float(myrank)))+1 726 intfromfile=0 727 call readheader(descriptor,'time derivative of solution' // char(0), 728 & intfromfile,ithree,'integer' // char(0),iotype) 729 allocate( acold(nshg,ndof) ) 730 if(intfromfile(1).ne.0) then 731 nshg2=intfromfile(1) 732 ndof2=intfromfile(2) 733 lstep=intfromfile(3) 734 735c print *, "intfromfile(1) is ", intfromfile(1) 736c print *, "intfromfile(2) is ", intfromfile(2) 737c print *, "intfromfile(3) is ", intfromfile(3) 738 739 if (nshg2 .ne. nshg) 740 & call error ('restar ', 'nshg ', nshg) 741c 742 allocate( acread(nshg,ndof2) ) 743 acread=zero 744 iacsiz=nshg*ndof2 745 call readdatablock(descriptor, 746 & 'time derivative of solution' // char(0),acread, 747 & iacsiz, 'double' // char(0),iotype) 748 acold(:,1:ndof)=acread(:,1:ndof) 749 deallocate(acread) 750 else 751 if (myrank.eq.master) then 752 warning='Time derivative of solution is set to zero (SAFE)' 753 write(*,*) warning 754 endif 755 acold=zero 756 endif 757 758cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc 759cc 760c 761cc 762cc.... read the header and check it against the run data 763cc 764 765 766c ithree=3 767ccc call creadlist(irstin,ithree,nshg2,ndof2,lstep) 768c fname1='solution?' 769c intfromfile=0 770c call readheader(irstin,fname1,intfromfile, 771c & ithree,'integer', iotype) 772cc 773cc.... read the values of primitive variables into q 774cc 775c allocate( qold(nshg,ndof) ) 776c if(intfromfile(1).ne.0) then 777c nshg2=intfromfile(1) 778c ndof2=intfromfile(2) 779c lstep=intfromfile(3) 780c if(ndof2.ne.ndof) then 781c warning='WARNING more data in restart than needed: keeping 1st ' 782c write(*,*) warning , ndof 783c endif 784cc 785c if (nshg2 .ne. nshg) 786c & call error ('restar ', 'nshg ', nshg) 787c allocate( qread(nshg,ndof2) ) 788 789c iqsiz=nshg*ndof2 790c call readdatablock(irstin,fname1,qread,iqsiz, 791c & 'double',iotype) 792c qold(:,1:ndof)=qread(:,1:ndof) 793c deallocate(qread) 794c else 795c if (myrank.eq.master) then 796c if (matflg(1,1).eq.0) then ! compressible 797c warning='Solution is set to zero (with p and T to one)' 798c else 799c warning='Solution is set to zero' 800c endif 801c write(*,*) warning 802c endif 803c qold=zero 804c if (matflg(1,1).eq.0) then ! compressible 805c qold(:,1)=one ! avoid zero pressure 806c qold(:,nflow)=one ! avoid zero temperature 807c endif 808c endif 809cc 810c fname1='time derivative of solution?' 811c intfromfile=0 812c call readheader(irstin,fname1,intfromfile, 813c & ithree,'integer', iotype) 814c allocate( acold(nshg,ndof) ) 815c if(intfromfile(1).ne.0) then 816c nshg2=intfromfile(1) 817c ndof2=intfromfile(2) 818c lstep=intfromfile(3) 819c 820c if (nshg2 .ne. nshg) 821c & call error ('restar ', 'nshg ', nshg) 822cc 823c allocate( acread(nshg,ndof2) ) 824c acread=zero 825c 826c iacsiz=nshg*ndof2 827c call readdatablock(irstin,fname1,acread,iacsiz, 828c & 'double',iotype) 829c acold(:,1:ndof)=acread(:,1:ndof) 830c deallocate(acread) 831c else 832c if (myrank.eq.master) then 833c warning='Time derivative of solution is set to zero (SAFE)' 834c write(*,*) warning 835c endif 836c acold=zero 837c endif 838 839c call creadlist(irstin,ithree,nshg2,ndisp,lstep) 840 841 if (ideformwall.eq.1) then 842! fname1='displacement?' 843! call readheader(irstin,fname1,intfromfile, 844! & ithree,'integer', iotype) 845 846 intfromfile=0 847 call readheader(descriptor,'displacement' // char(0), 848 & intfromfile,ithree,'integer' // char(0),iotype) 849 850 nshg2=intfromfile(1) 851 ndisp=intfromfile(2) 852 lstep=intfromfile(3) 853 if(ndisp.ne.nsd) then 854 warning='WARNING ndisp not equal nsd' 855 write(*,*) warning , ndisp 856 endif 857c 858 if (nshg2 .ne. nshg) 859 & call error ('restar ', 'nshg ', nshg) 860c 861c.... read the values of primitive variables into uold 862c 863 864 allocate( uold(nshg,nsd) ) 865 allocate( uread(nshg,nsd) ) 866 867 iusiz=nshg*nsd 868 869! call readdatablock(irstin,fname1,uread,iusiz, 870! & 'double',iotype) 871 call readdatablock(descriptor,fname1 // char(0) ,uread,iusiz, 872 & 'double' // char(0),iotype) 873 874 uold(:,1:nsd)=uread(:,1:nsd) 875 else 876 allocate( uold(nshg,nsd) ) 877 uold(:,1:nsd) = zero 878 endif 879 880c 881c.... close c-binary files 882c 883!MR CHANGE 884! call closefile( irstin, "read" ) 885 886 call closefile( descriptor, "read" // char(0) ) 887 call finalizephmpiio( descriptor ) 888 889!MR CHANGE 890! call closefile( igeomBAK, "read" ) 891c 892 deallocate(xread) 893 if ( numpbc > 0 ) then 894 deallocate(bcinpread) 895 deallocate(ibctmpread) 896 endif 897 deallocate(iperread) 898 if(numpe.gt.1) 899 & deallocate(ilworkread) 900 deallocate(nbcread) 901 902 return 903c 904 994 call error ('input ','opening ', igeomBAK) 905 995 call error ('input ','opening ', igeomBAK) 906 997 call error ('input ','end file', igeomBAK) 907 998 call error ('input ','end file', igeomBAK) 908c 909 end 910 911c 912c No longer called but kept around in case.... 913c 914 subroutine genpzero(iBC) 915 916 use pointer_data 917c 918 include "common.h" 919 integer iBC(nshg) 920c 921c.... check to see if any of the nodes have a dirichlet pressure 922c 923 pzero=1 924 if (any(btest(iBC,2))) pzero=0 925c 926 do iblk = 1, nelblb 927 npro = lcblkb(1,iblk+1)-lcblkb(1,iblk) 928 do i=1, npro 929 iBCB1=miBCB(iblk)%p(i,1) 930c 931c.... check to see if any of the nodes have a Neumann pressure 932c but not periodic (note that 933c 934 if(btest(iBCB1,1)) pzero=0 935 enddo 936c 937c.... share results with other processors 938c 939 pzl=pzero 940 if (numpe .gt. 1) 941 & call MPI_ALLREDUCE (pzl, pzero, 1, 942 & MPI_DOUBLE_PRECISION,MPI_MIN, MPI_COMM_WORLD,ierr) 943 944 enddo 945c 946c.... return 947c 948 return 949c 950 end 951 952