1 subroutine genbkb (ibksz) 2c 3c---------------------------------------------------------------------- 4c 5c This routine reads the boundary elements, reorders them and 6c generates traces for the gather/scatter operations. 7c 8c Zdenek Johan, Fall 1991. 9c---------------------------------------------------------------------- 10c 11 use dtnmod 12 use pointer_data 13c 14 include "common.h" 15!MR CHANGE 16 include "mpif.h" !Required to determine the max for itpblk 17!MR CHANGE END 18c 19 20 integer, allocatable :: ientp(:,:),iBCBtp(:,:) 21 real*8, allocatable :: BCBtp(:,:) 22 integer materb(ibksz) 23 integer intfromfile(50) ! integers read from headers 24 character*255 fname1 25 26cccccccccccccc New Phasta IO starts here cccccccccccccccccccccccccccccc 27 28 integer :: descriptor, descriptorG, GPID, color, nfiles, nfields 29 integer :: numparts, nppf, nppp, nprocs, writeLock 30 integer :: ierr_io, numprocs, itmp, itmp2 31! integer :: num_local_loop, num_global_loop 32!MR CHANGE 33 integer :: itpblktot,ierr 34!MR CHANGE END 35 36 character*255 fnamer, fname2, temp2 37 character*64 temp1, temp3 38 39 nfiles = nsynciofiles 40 nfields = nsynciofieldsreadgeombc 41 numparts = numpe !This is the common settings. Beware if you try to compute several parts per process 42 43 nppp = numparts/numpe 44 nppf = numparts/nfiles 45 46 color = int(myrank/(numparts/nfiles)) 47 itmp2 = int(log10(float(color+1)))+1 48 write (temp2,"('(''geombc-dat.'',i',i1,')')") itmp2 49 temp2=trim(temp2) 50 write (fnamer,temp2) (color+1) 51 fnamer=trim(fnamer) 52 53 ione=1 54 itwo=2 55 ieight=8 56 ieleven=11 57 itmp = int(log10(float(myrank+1)))+1 58 59! num_global_loop = 4 60! num_local_loop = nelblb 61 62cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc 63 64 iel=1 65 itpblk=nelblb 66 67 68!MR CHANGE 69 ! Get the total number of different interior topologies in the whole domain. 70 ! Try to read from a field. If the field does not exist, scan the geombc file. 71 itpblktot=-1 72 write(temp1, 73 & "('(''total number of boundary tpblocks@'',i',i1,',A1)')") itmp 74 write (fname2,temp1) (myrank+1),'?' 75 call readheader(igeom,fname2 // char(0) ,itpblktot,ione, 76 & 'integer' // char(0),iotype) 77 78! write (*,*) 'Rank: ',myrank,' boundary itpblktot intermediate:', 79! & itpblktot 80 81 if (itpblktot == -1) then 82 ! The field 'total number of different boundary tpblocks' was not found in the geombc file. 83 ! Scan all the geombc file for the 'connectivity interior' fields to get this information. 84 iblk=0 85 neltp=0 86 do while(neltp .ne. -1) 87 88 ! intfromfile is reinitialized to -1 every time. 89 ! If connectivity boundary@xxx is not found, then 90 ! readheader will return intfromfile unchanged 91 92 intfromfile(:)=-1 93 iblk = iblk+1 94 write (temp1,"('connectivity boundary',i1)") iblk 95 temp1 = trim(temp1) 96 write (temp3,"('(''@'',i',i1,',A1)')") itmp 97 write (fname2, temp3) (myrank+1), '?' 98 fname2 = trim(temp1)//trim(fname2) 99 !write(*,*) 'rank, fname2',myrank, trim(adjustl(fname2)) 100 call readheader(igeom,fname2 // char(0),intfromfile, 101 & ieight,'integer' // char(0),iotype) 102 neltp = intfromfile(1) ! -1 if fname2 was not found, >=0 otherwise 103 end do 104 itpblktot = iblk-1 105 end if 106 107 if (myrank == 0) then 108 write(*,*) 'Number of boundary topologies: ',itpblktot 109 endif 110! write (*,*) 'Rank: ',myrank,' boundary itpblktot final:', 111! & itpblktot 112 113!MR CHANGE END 114 nelblb=0 115 mattyp=0 116 ndofl = ndof 117!MR CHANGE 118! call initphmpiio( nfields, nppf, nfiles, igeom ) 119! call openfile( fnamer, 'read', igeom ) 120 121 do iblk = 1, itpblktot 122 writeLock=0; 123!MR CHANGE END 124 125 fname1='connectivity boundary?' 126 127! print *, "Loop ",iblk, myrank, itpblk, trim(fnamer) 128 129ccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc 130 131 write (temp1,"('connectivity boundary',i1)") iblk 132 temp1 = trim(temp1) 133 write (temp3,"('(''@'',i',i1,',A1)')") itmp 134 write (fname2, temp3) (myrank+1), '?' 135 fname2 = trim(temp1)//trim(fname2) 136 fname2 = trim(fname2) 137! write(*,*) 'rank, fname2',myrank, trim(adjustl(fname2)) 138 139ccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc 140 141 ! Synchronization for performance monitoring, as some parts do not include some topologies 142 call MPI_Barrier(MPI_COMM_WORLD,ierr) 143 call readheader(igeom,fname2 // char(0),intfromfile,ieight, 144 & 'integer' // char(0),iotype) 145 neltp =intfromfile(1) 146 nenl =intfromfile(2) 147 ipordl=intfromfile(3) 148 nshl =intfromfile(4) 149 nshlb =intfromfile(5) 150 nenbl =intfromfile(6) 151 lcsyst=intfromfile(7) 152 numnbc=intfromfile(8) 153 154!MR CHANGE 155! write (temp1,"('connectivityBoundaryHeader_',i1,'_',i1)") 156! & iblk,myrank 157! temp1=trim(temp1) 158! open(unit=14,file=temp1,status='unknown') 159! write(14,*) intfromfile(:) 160! close(14) 161!MR CHANGE END 162 163 allocate (ientp(neltp,nshl)) 164 allocate (iBCBtp(neltp,ndiBCB)) 165 allocate (BCBtp(neltp,ndBCB)) 166 iientpsiz=neltp*nshl 167 168 if (neltp==0) then 169 writeLock=1; 170 endif 171 172! print *, "neltp is ", neltp 173 174 call readdatablock(igeom,fname2 // char(0),ientp,iientpsiz, 175 & 'integer' // char(0),iotype) 176 177!MR CHANGE 178! write (temp1,"('connectivityBoundaryDatablock_',i1,'_',i1)") 179! & iblk,myrank 180! temp1=trim(temp1) 181 182! open(unit=14,file=temp1,status='unknown') 183! do i=1,neltp 184! do j=1,nshl 185! write(14,*) ientp(i,j) 186! enddo 187! enddo 188! write(14,*) iientpsiz 189! close(14) 190!MR CHANGE END 191 192 193c 194c.... Read the boundary flux codes 195c 196 197! print *,"connectivity [] is ", trim(fname2),ientp(0,0) 198 199 200 fname1='nbc codes?' 201 202ccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc 203 204 write (temp1,"('nbc codes',i1)") iblk 205 temp1=trim(temp1) 206 write (temp3,"('(''@'',i',i1,',A1)')") itmp 207 write (fname2, temp3) (myrank+1), '?' 208 fname2 = trim(temp1)//trim(fname2) 209! write(*,*) 'rank, fname2',myrank, trim(adjustl(fname2)) 210 call MPI_BARRIER(MPI_COMM_WORLD, ierr) 211 212ccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc 213 214 call readheader(igeom,fname2 // char(0) ,intfromfile, 215 & ieight,'integer' // char(0),iotype) 216 iiBCBtpsiz=neltp*ndiBCB 217 call readdatablock(igeom,fname2 // char(0) ,iBCBtp, 218 & iiBCBtpsiz,'integer' // char(0),iotype) 219 220!MR CHANGE 221! print *, "ndiBCB is ",ndiBCB 222! write (temp1,"('nbcCodesDatablock_',i1,'_',i1)") 223! & iblk,myrank 224! temp1=trim(temp1) 225! 226! open(unit=13,file=temp1,status='unknown') 227! do i=1,neltp 228! do j=1,ndiBCB 229! write(13,*) iBCBtp(i,j) 230! enddo 231! enddo 232! write(13,*) iiBCBtpsiz 233! close(13) 234!! iBCBtp(:,:) = 0 ! JUST FOR TEST 235!MR CHANGE END 236 237c 238c.... read the boundary condition data 239c 240 fname1='nbc values?' 241 242ccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc 243 244 write (temp1,"('nbc values',i1)") iblk 245 temp1=trim(temp1) 246 write (temp3,"('(''@'',i',i1,',A1)')") itmp 247 write (fname2, temp3) (myrank+1), '?' 248 fname2 = trim(temp1)//trim(fname2) 249 call MPI_BARRIER(MPI_COMM_WORLD, ierr) 250ccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc 251 252 call readheader(igeom,fname2 // char(0) ,intfromfile, 253 & ieight,'integer' // char(0) ,iotype) 254 BCBtp = zero 255 iBCBtpsiz=neltp*ndBCB 256 call readdatablock(igeom,fname2 // char(0),BCBtp,iBCBtpsiz, 257 & 'double' // char(0) ,iotype) 258 259!MR CHANGE 260! write (temp1,"('nbcValuesDatablock_',i1,'_',i1)") 261! & iblk,myrank 262! temp1=trim(temp1) 263! open(unit=12,file=temp1,status='unknown') 264! do i=1,neltp 265! do j=1,ndBCB 266! write(12,*) BCBtp(i,j) 267! enddo 268! enddo 269! write(12,*) iBCBtpsiz 270! close(12) 271!MR CHANGE END 272 273c 274c This is a temporary fix until NSpre properly zeros this array where it 275c is not set. DEC has indigestion with these arrays though the 276c result is never used (never effects solution). 277c 278 279 280!MR CHANGE 281 if(writeLock==0) then 282!MR CHANGE 283! print *,"In ASSIGN ASSIGN",myrank 284 where(.not.btest(iBCBtp(:,1),0)) BCBtp(:,1)=zero 285 where(.not.btest(iBCBtp(:,1),1)) BCBtp(:,2)=zero 286 where(.not.btest(iBCBtp(:,1),3)) BCBtp(:,6)=zero 287 if(ndBCB.gt.6) then 288 do i=6,ndof 289 where(.not.btest(iBCBtp(:,1),i-1)) BCBtp(:,i+1)=zero 290 enddo 291 endif 292 where(.not.btest(iBCBtp(:,1),2)) 293 BCBtp(:,3)=zero 294 BCBtp(:,4)=zero 295 BCBtp(:,5)=zero 296 endwhere 297 298 do n=1,neltp,ibksz 299 nelblb=nelblb+1 300 npro= min(IBKSZ, neltp - n + 1) 301c 302 lcblkb(1,nelblb) = iel 303c lcblkb(2,nelblb) = iopen ! available for later use 304 lcblkb(3,nelblb) = lcsyst 305 lcblkb(4,nelblb) = ipordl 306 lcblkb(5,nelblb) = nenl 307 lcblkb(6,nelblb) = nenbl 308 lcblkb(7,nelblb) = mattyp 309 lcblkb(8,nelblb) = ndofl 310 lcblkb(9,nelblb) = nshl 311 lcblkb(10,nelblb) = nshlb ! # of shape functions per elt 312c 313c.... save the element block 314c 315 n1=n 316 n2=n+npro-1 317 materb=1 ! all one material for now 318c 319c.... allocate memory for stack arrays 320c 321 322 allocate (mienb(nelblb)%p(npro,nshl)) 323c 324 allocate (miBCB(nelblb)%p(npro,ndiBCB)) 325c 326 allocate (mBCB(nelblb)%p(npro,nshlb,ndBCB)) 327c 328 allocate (mmatb(nelblb)%p(npro)) 329c 330c.... save the boundary element block 331c 332 call gensvb (ientp(n1:n2,1:nshl), 333 & iBCBtp(n1:n2,:), BCBtp(n1:n2,:), 334 & materb, mienb(nelblb)%p, 335 & miBCB(nelblb)%p, mBCB(nelblb)%p, 336 & mmatb(nelblb)%p) 337c 338 iel=iel+npro 339 enddo 340 341!MR CHANGE 342 endif 343!MR CHANGE 344 deallocate(ientp) 345 deallocate(iBCBtp) 346 deallocate(BCBtp) 347 348 enddo 349 lcblkb(1,nelblb+1) = iel 350 351! call closefile( igeom, "read" ) 352! call finalizephmpiio( igeom ) 353 354c 355c.... return 356c 357 return 358c 359c.... end of file error handling 360c 361 911 call error ('genbcb ','end file',igeomBAK) 362c 3631000 format(a80,//, 364 & ' B o u n d a r y E l e m e n t C o n n e c t i v i t y',//, 365 & ' Elem BC codes',/, 366 & ' Number C P V H ',5x,27('Node',i1,:,2x)) 3671100 format(2x,i5,2x,4i2,3x,27i7) 368c$$$2000 format(a80,//, 369c$$$ & ' B o u n d a r y E l e m e n t B C D a t a ',//, 370c$$$ & ' Node ',3x,'mass',/, 371c$$$ & ' Number ',3x,'flux',6x,'Pressure',6x,'Heat',6x, 372c$$$ & 3('Viscous',i1,:,4x)) 3732100 format(2x,i5,1p,1x,6e12.4) 374c 375 end 376 377