xref: /phasta/phSolver/common/genblk.f (revision dc9538425e8e2132ad0d36a4b7a0472f4a4e7110)
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, nfiles
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        nfiles = nsynciofiles
29        numparts = numpe !This is the common settings. Beware if you try to compute several parts per process
30        ione=1
31        itwo=2
32        iseven=7
33        ieleven=11
34        iel=1
35        itpblk=nelblk
36
37        ! Get the total number of different interior topologies in the whole domain.
38        ! Try to read from a field. If the field does not exist, scan the geombc file.
39          itpblktot=1  ! hardwired to montopology for now
40        call phio_readheader(fhandle,
41     &   c_char_'total number of interior tpblocks' // char(0),
42     &   c_loc(itpblktot), ione, dataInt, iotype)
43
44        if (itpblktot == -1) then
45          ! The field 'total number of different interior tpblocks' was not found in the geombc file.
46          ! Scan all the geombc file for the 'connectivity interior' fields to get this information.
47          iblk=0
48          neltp=0
49          do while(neltp .ne. -1)
50
51            ! intfromfile is reinitialized to -1 every time.
52            ! If connectivity interior@xxx is not found, then
53            ! readheader will return intfromfile unchanged
54
55            intfromfile(:)=-1
56            iblk = iblk+1
57            if(nfiles.gt.0) then
58              write (fname2,"('connectivity interior',i1)") iblk
59            else
60              write (fname2,"('connectivity interior linear tetrahedron')")
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(nfiles.gt.0) then
83              write (fname2,"('connectivity interior',i1)") iblk
84            else
85              write (fname2,"('connectivity interior linear tetrahedron')")
86            endif
87
88           ! Synchronization for performance monitoring, as some parts do not include some topologies
89           call MPI_Barrier(MPI_COMM_WORLD,ierr)
90           call phio_readheader(fhandle, fname2 // char(0),
91     &      c_loc(intfromfile), iseven, dataInt, iotype)
92           neltp  =intfromfile(1)
93           nenl   =intfromfile(2)
94           ipordl =intfromfile(3)
95           nshl   =intfromfile(4)
96           ijunk  =intfromfile(5)
97           ijunk  =intfromfile(6)
98           lcsyst =intfromfile(7)
99           allocate (ientp(neltp,nshl))
100           iientpsiz=neltp*nshl
101
102           if (neltp==0) then
103              writeLock=1;
104           endif
105
106           call phio_readdatablock(fhandle,fname2 // char(0),
107     &      c_loc(ientp), iientpsiz, dataInt, iotype)
108
109           if(writeLock==0) then
110             do n=1,neltp,ibksz
111                nelblk=nelblk+1
112                npro= min(IBKSZ, neltp - n + 1)
113                lcblk(1,nelblk)  = iel
114                lcblk(3,nelblk)  = lcsyst
115                lcblk(4,nelblk)  = ipordl
116                lcblk(5,nelblk)  = nenl
117                lcblk(6,nelblk)  = nfacel
118                lcblk(7,nelblk)  = mattyp
119                lcblk(8,nelblk)  = ndofl
120                lcblk(9,nelblk)  = nsymdl
121                lcblk(10,nelblk) = nshl ! # of shape functions per elt
122c
123c.... allocate memory for stack arrays
124c
125                allocate (mmat(nelblk)%p(npro))
126c
127                allocate (mien(nelblk)%p(npro,nshl))
128                allocate (mxmudmi(nelblk)%p(npro,maxsh))
129c
130c.... save the element block
131c
132                n1=n
133                n2=n+npro-1
134                mater=1   ! all one material for now
135                call gensav (ientp(n1:n2,1:nshl),
136     &                       mater,           mien(nelblk)%p,
137     &                       mmat(nelblk)%p)
138                iel=iel+npro
139             enddo
140           endif
141           deallocate(ientp)
142        enddo
143
144        lcblk(1,nelblk+1) = iel
145        return
1461000    format(a80,//,
147     &  ' N o d a l   C o n n e c t i v i t y',//,
148     &  '   Elem  ',/,
149     &  '  Number  ',7x,27('Node',i2,:,2x))
1501100    format(2x,i5,6x,27i8)
151        end
152