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 14 include "common.h" 15 include "mpif.h" !Required to determine the max for itpblk 16 17 integer, target, allocatable :: ientp(:,:) 18 integer mater(ibksz) 19 integer, target :: intfromfile(50) ! integers read from headers 20 character*255 fname1 21 integer :: descriptor, descriptorG, GPID, color 22 integer :: numparts, writeLock 23 integer :: ierr_io, numprocs 24 integer, target :: itpblktot,ierr,iseven 25 character*255 fname2 26 character(len=30) :: dataInt 27 dataInt = c_char_'integer'//c_null_char 28 numparts = numpe !This is the common settings. Beware if you try to compute several parts per process 29 ione=1 30 itwo=2 31 iseven=7 32 ieleven=11 33 iel=1 34 itpblk=nelblk 35 36 ! Get the total number of different interior topologies in the whole domain. 37 ! Try to read from a field. If the field does not exist, scan the geombc file. 38 itpblktot=1 ! hardwired to montopology for now 39 call phio_readheader(fhandle, 40 & c_char_'total number of interior tpblocks' // char(0), 41 & c_loc(itpblktot), ione, dataInt, iotype) 42 43 if (itpblktot == -1) then 44 ! The field 'total number of different interior tpblocks' was not found in the geombc file. 45 ! Scan all the geombc file for the 'connectivity interior' fields to get this information. 46 iblk=0 47 neltp=0 48 do while(neltp .ne. -1) 49 50 ! intfromfile is reinitialized to -1 every time. 51 ! If connectivity interior@xxx is not found, then 52 ! readheader will return intfromfile unchanged 53 54 intfromfile(:)=-1 55 iblk = iblk+1 56 if(input_mode.ge.1) then 57 write (fname2,"('connectivity interior',i1)") iblk 58 else 59 write (fname2,"('connectivity interior linear tetrahedron')") 60! write (fname2,"('connectivity interior?')") 61 endif 62 63 !write(*,*) 'rank, fname2',myrank, trim(adjustl(fname2)) 64 call phio_readheader(fhandle, fname2 // char(0), 65 & c_loc(intfromfile), iseven, dataInt, iotype) 66 neltp = intfromfile(1) ! -1 if fname2 was not found, >=0 otherwise 67 end do 68 itpblktot = iblk-1 69 end if 70 71 if (myrank == 0) then 72 write(*,*) 'Number of interior topologies: ',itpblktot 73 endif 74 75 nelblk=0 76 mattyp = 0 77 ndofl = ndof 78 nsymdl = nsymdf 79 80 do iblk = 1, itpblktot 81 writeLock=0; 82 if(input_mode.ge.1) then 83 write (fname2,"('connectivity interior',i1)") iblk 84 else 85 write (fname2,"('connectivity interior linear tetrahedron')") 86! write (fname2,"('connectivity interior?')") 87 endif 88 89 ! Synchronization for performance monitoring, as some parts do not include some topologies 90 call MPI_Barrier(MPI_COMM_WORLD,ierr) 91 call phio_readheader(fhandle, fname2 // char(0), 92 & c_loc(intfromfile), iseven, dataInt, iotype) 93 neltp =intfromfile(1) 94 nenl =intfromfile(2) 95 ipordl =intfromfile(3) 96 nshl =intfromfile(4) 97 ijunk =intfromfile(5) 98 ijunk =intfromfile(6) 99 lcsyst =intfromfile(7) 100 allocate (ientp(neltp,nshl)) 101 iientpsiz=neltp*nshl 102 103 if (neltp==0) then 104 writeLock=1; 105 endif 106 107 call phio_readdatablock(fhandle,fname2 // char(0), 108 & c_loc(ientp), iientpsiz, dataInt, iotype) 109 110 if(writeLock==0) then 111 do n=1,neltp,ibksz 112 nelblk=nelblk+1 113 npro= min(IBKSZ, neltp - n + 1) 114 lcblk(1,nelblk) = iel 115 lcblk(3,nelblk) = lcsyst 116 lcblk(4,nelblk) = ipordl 117 lcblk(5,nelblk) = nenl 118 lcblk(6,nelblk) = nfacel 119 lcblk(7,nelblk) = mattyp 120 lcblk(8,nelblk) = ndofl 121 lcblk(9,nelblk) = nsymdl 122 lcblk(10,nelblk) = nshl ! # of shape functions per elt 123c 124c.... allocate memory for stack arrays 125c 126 allocate (mmat(nelblk)%p(npro)) 127c 128 allocate (mien(nelblk)%p(npro,nshl)) 129 allocate (mxmudmi(nelblk)%p(npro,maxsh)) 130 if(usingpetsc.eq.0) then 131 allocate (mienG(nelblk)%p(1,1)) 132 else 133 allocate (mienG(nelblk)%p(npro,nshl)) 134 endif 135 ! note mienG will be passed to gensav but nothing filled if not 136 ! using PETSc so this is safe 137c 138c.... save the element block 139c 140 n1=n 141 n2=n+npro-1 142 mater=1 ! all one material for now 143 call gensav (ientp(n1:n2,1:nshl), 144 & mater, mien(nelblk)%p, 145 & mienG(nelblk)%p, 146 & mmat(nelblk)%p) 147 iel=iel+npro 148 enddo 149 endif 150 deallocate(ientp) 151 enddo 152 153 lcblk(1,nelblk+1) = iel 154 return 1551000 format(a80,//, 156 & ' N o d a l C o n n e c t i v i t y',//, 157 & ' Elem ',/, 158 & ' Number ',7x,27('Node',i2,:,2x)) 1591100 format(2x,i5,6x,27i8) 160 end 161