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