1 subroutine genblk (IBKSZ) 2c 3c---------------------------------------------------------------------- 4c 5c This routine reads the interior elements and generates the 6c appropriate blocks. 7c 8c Zdenek Johan, Fall 1991. 9c---------------------------------------------------------------------- 10c 11 use pointer_data 12 use phio 13 use iso_c_binding 14c 15 include "common.h" 16!MR CHANGE 17 include "mpif.h" !Required to determine the max for itpblk 18!MR CHANGE END 19c 20 integer, target, allocatable :: ientp(:,:) 21 integer mater(ibksz) 22 integer, target :: intfromfile(50) ! integers read from headers 23 character*255 fname1 24 25cccccccccccccc New Phasta IO starts here ccccccccccccccccccccccccc 26 27 integer :: descriptor, descriptorG, GPID, color, nfiles 28 integer :: numparts, writeLock 29 integer :: ierr_io, numprocs 30!MR CHANGE 31 integer, target :: itpblktot,ierr,iseven 32!MR CHANGE END 33 character*255 fname2 34 35 character(len=30) :: dataInt 36 dataInt = c_char_'integer'//c_null_char 37 38!THIS NEEDS TO BE CLEANED - MR 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 44 ione=1 45 itwo=2 46 iseven=7 47 ieleven=11 48 49cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc 50 51c 52 iel=1 53 itpblk=nelblk 54!MR CHANGE 55 56 ! Get the total number of different interior topologies in the whole domain. 57 ! Try to read from a field. If the field does not exist, scan the geombc file. 58 itpblktot=1 ! hardwired to montopology for now 59 call phio_readheader(fhandle, 60 & c_char_'total number of interior tpblocks' // char(0), 61 & c_loc(itpblktot), ione, dataInt, iotype) 62 63! write (*,*) 'Rank: ',myrank,' interior itpblktot intermediate:', 64! & itpblktot 65 66 if (itpblktot == -1) then 67 ! The field 'total number of different interior tpblocks' was not found in the geombc file. 68 ! Scan all the geombc file for the 'connectivity interior' fields to get this information. 69 iblk=0 70 neltp=0 71 do while(neltp .ne. -1) 72 73 ! intfromfile is reinitialized to -1 every time. 74 ! If connectivity interior@xxx is not found, then 75 ! readheader will return intfromfile unchanged 76 77 intfromfile(:)=-1 78 iblk = iblk+1 79 if(nfiles.gt.0) then 80 write (fname2,"('connectivity interior',i1)") iblk 81 else 82 write (fname2,"('connectivity interior linear tetrahedron')") 83 endif 84 85 !write(*,*) 'rank, fname2',myrank, trim(adjustl(fname2)) 86 call phio_readheader(fhandle, fname2 // char(0), 87 & c_loc(intfromfile), iseven, dataInt, iotype) 88 neltp = intfromfile(1) ! -1 if fname2 was not found, >=0 otherwise 89 end do 90 itpblktot = iblk-1 91 end if 92 93 if (myrank == 0) then 94 write(*,*) 'Number of interior topologies: ',itpblktot 95 endif 96! write (*,*) 'Rank: ',myrank,' interior itpblktot final:', 97! & itpblktot 98 99!MR CHANGE END 100 101 nelblk=0 102 mattyp = 0 103 ndofl = ndof 104 nsymdl = nsymdf 105 106! call initphmpiio( nfields, nppf, nfiles, igeom ) 107! call openfile( fnamer, 'read', igeom ) 108 109! do iblk = 1, itpblk 110 do iblk = 1, itpblktot 111 writeLock=0; 112!MR CHANGE END 113c 114c read(igeomBAK) neltp,nenl,ipordl,nshl, ijunk, ijunk, lcsyst 115c call creadlist(igeomBAK,iseven, 116c & neltp,nenl,ipordl,nshl, ijunk, ijunk, lcsyst) 117 118ccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc 119 120 if(nfiles.gt.0) then 121 write (fname2,"('connectivity interior',i1)") iblk 122 else 123 write (fname2,"('connectivity interior linear tetrahedron')") 124 endif 125 126ccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc 127 128c fname1='connectivity interior?' 129 130 ! Synchronization for performance monitoring, as some parts do not include some topologies 131 call MPI_Barrier(MPI_COMM_WORLD,ierr) 132 call phio_readheader(fhandle, fname2 // char(0), 133 & c_loc(intfromfile), iseven, dataInt, iotype) 134 neltp =intfromfile(1) 135 nenl =intfromfile(2) 136 ipordl =intfromfile(3) 137 nshl =intfromfile(4) 138 ijunk =intfromfile(5) 139 ijunk =intfromfile(6) 140 lcsyst =intfromfile(7) 141 allocate (ientp(neltp,nshl)) 142c read(igeomBAK) ientp 143 iientpsiz=neltp*nshl 144 145 if (neltp==0) then 146 writeLock=1; 147 endif 148 149 call phio_readdatablock(fhandle,fname2 // char(0), 150 & c_loc(ientp), iientpsiz, dataInt, iotype) 151 152! call closefile( igeom, "read" // char(0) ) 153! call finalizephmpiio( igeom ) 154 155!MR CHANGE 156 if(writeLock==0) then 157!MR CHANGE 158 159 do n=1,neltp,ibksz 160 nelblk=nelblk+1 161 npro= min(IBKSZ, neltp - n + 1) 162c 163 lcblk(1,nelblk) = iel 164c lcblk(2,nelblk) = iopen ! available for later use 165 lcblk(3,nelblk) = lcsyst 166 lcblk(4,nelblk) = ipordl 167 lcblk(5,nelblk) = nenl 168 lcblk(6,nelblk) = nfacel 169 lcblk(7,nelblk) = mattyp 170 lcblk(8,nelblk) = ndofl 171 lcblk(9,nelblk) = nsymdl 172 lcblk(10,nelblk) = nshl ! # of shape functions per elt 173c 174c.... allocate memory for stack arrays 175c 176 allocate (mmat(nelblk)%p(npro)) 177c 178 allocate (mien(nelblk)%p(npro,nshl)) 179 allocate (mxmudmi(nelblk)%p(npro,maxsh)) 180c 181c.... save the element block 182c 183 n1=n 184 n2=n+npro-1 185 mater=1 ! all one material for now 186 call gensav (ientp(n1:n2,1:nshl), 187 & mater, mien(nelblk)%p, 188 & mmat(nelblk)%p) 189 iel=iel+npro 190c 191 enddo 192!MR CHANGE 193 endif 194!MR CHANGE 195 deallocate(ientp) 196 enddo 197 198! call closefile( igeom, "read" // char(0) ) 199! call finalizephmpiio( igeom ) 200 201 lcblk(1,nelblk+1) = iel 202c 203c.... return 204c 205CAD call timer ('Back ') 206c 207 return 208c 2091000 format(a80,//, 210 & ' N o d a l C o n n e c t i v i t y',//, 211 & ' Elem ',/, 212 & ' Number ',7x,27('Node',i2,:,2x)) 2131100 format(2x,i5,6x,27i8) 214 end 215