xref: /phasta/phSolver/common/genblk.f (revision fe88b52d35537d3736cc3d33d9db7e51353aba7d)
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