Line data Source code
1 : module metagraph_mod
2 : use cam_logfile, only: iulog
3 : use gridgraph_mod, only : gridvertex_t, gridedge_t, &
4 : allocate_gridvertex_nbrs, assignment ( = )
5 :
6 : implicit none
7 : private
8 :
9 : type, public :: MetaEdge_t
10 : type (GridEdge_t),pointer :: members(:)
11 : integer ,pointer :: edgeptrP(:)
12 : integer ,pointer :: edgeptrP_ghost(:)
13 : integer ,pointer :: edgeptrS(:)
14 : integer :: number
15 : integer :: type
16 : integer :: wgtP ! sum of lengths of all messages to pack for edges
17 : integer :: wgtP_ghost ! sum of lengths of all messages to pack for ghost cells
18 : integer :: wgtS
19 : integer :: HeadVertex ! processor number to send to
20 : integer :: TailVertex ! processor number to send from
21 : integer :: nmembers ! number of messages to (un)pack (out)into this buffer
22 : integer :: padding ! just to quite compiler
23 : end type MetaEdge_t
24 :
25 : type, public :: MetaVertex_t ! one for each processor
26 : integer :: number ! USELESS just the local processor number
27 : integer :: nmembers ! number of elements on this processor
28 : type (GridVertex_t),pointer :: members(:) ! array of elements on this processor
29 : type (MetaEdge_t),pointer :: edges(:) ! description of messages to send/receive
30 : integer :: nedges ! number of processors to communicate with (length of edges)
31 : integer :: padding ! just to quite compiler
32 : end type MetaVertex_t
33 :
34 :
35 : public :: edge_uses_vertex
36 : public :: PrintMetaEdge, PrintMetaVertex
37 : public :: LocalElemCount
38 : public :: initMetaGraph
39 :
40 : interface assignment ( = )
41 : module procedure copy_metaedge
42 : end interface
43 :
44 : CONTAINS
45 :
46 : ! =====================================
47 : ! copy vertex:
48 : ! copy device for overloading = sign.
49 : ! =====================================
50 :
51 : recursive subroutine copy_metaedge(edge2,edge1)
52 :
53 : type (MetaEdge_t), intent(out) :: edge2
54 : type (MetaEdge_t), intent(in) :: edge1
55 :
56 : integer i
57 :
58 : edge2%number = edge1%number
59 : edge2%type = edge1%type
60 : edge2%wgtP = edge1%wgtP
61 : edge2%wgtP_ghost = edge1%wgtP_ghost
62 : edge2%nmembers = edge1%nmembers
63 :
64 : if (associated(edge1%members)) then
65 : allocate(edge2%members(edge2%nmembers))
66 : do i=1,edge2%nmembers
67 : edge2%members(i)=edge1%members(i)
68 : end do
69 : end if
70 :
71 : if (associated(edge1%edgeptrP)) then
72 : allocate(edge2%edgeptrP(edge2%nmembers))
73 : allocate(edge2%edgeptrS(edge2%nmembers))
74 : allocate(edge2%edgeptrP_ghost(edge2%nmembers))
75 : do i=1,edge2%nmembers
76 : edge2%edgeptrP(i)=edge1%edgeptrP(i)
77 : edge2%edgeptrS(i)=edge1%edgeptrS(i)
78 : edge2%edgeptrP_ghost(i)=edge1%edgeptrP_ghost(i)
79 : end do
80 : end if
81 :
82 : edge2%HeadVertex = edge1%HeadVertex
83 : edge2%TailVertex = edge1%TailVertex
84 :
85 : end subroutine copy_metaedge
86 :
87 1536 : function LocalElemCount(Vertex) result(nelemd)
88 :
89 : type (MetaVertex_t),intent(in) :: Vertex
90 : integer :: nelemd
91 :
92 1536 : nelemd = Vertex%nmembers
93 :
94 1536 : end function LocalElemCount
95 :
96 0 : function edge_uses_vertex(Vertex,Edge) result(log)
97 :
98 : type(MetaVertex_t), intent(in) :: Vertex
99 : type(MetaEdge_t), intent(in) :: Edge
100 : logical :: log
101 : integer :: number
102 :
103 0 : number = Vertex%number
104 0 : if(number == Edge%HeadVertex .or. number == Edge%TailVertex) then
105 : log = .TRUE.
106 : else
107 0 : log = .FALSE.
108 : endif
109 :
110 0 : end function edge_uses_vertex
111 :
112 0 : subroutine PrintMetaEdge(Edge)
113 : use gridgraph_mod, only : PrintGridEdge
114 : implicit none
115 : type (MetaEdge_t), intent(in) :: Edge(:)
116 : integer :: i,nedge
117 :
118 0 : nedge = SIZE(Edge)
119 0 : do i=1,nedge
120 0 : print *
121 0 : write(iulog,90) Edge(i)%number,Edge(i)%type,Edge(i)%wgtP,Edge(i)%nmembers, &
122 0 : Edge(i)%TailVertex, Edge(i)%HeadVertex
123 0 : if(associated(Edge(i)%members)) then
124 : call PrintGridEdge(Edge(i)%members)
125 : endif
126 : enddo
127 : 90 format('METAEDGE #',I4,2x,'TYPE ',I1,2x,'WGT ',I4,2x,'NUM ',I6,2x,'Processors ',I4,' ---> ',I4)
128 :
129 0 : end subroutine PrintMetaEdge
130 :
131 0 : subroutine PrintMetaVertex(Vertex)
132 : use gridgraph_mod, only : PrintGridVertex
133 : implicit none
134 : type (MetaVertex_t), intent(in),target :: Vertex
135 :
136 : integer :: j
137 :
138 :
139 0 : write(iulog,*)
140 0 : write(iulog,95) Vertex%nmembers
141 0 : call PrintGridVertex(Vertex%members)
142 0 : write(iulog,96) Vertex%nedges
143 0 : if(associated(Vertex%edges)) then
144 0 : do j=1,Vertex%nedges
145 0 : write(iulog,97) Vertex%edges(j)%number, Vertex%edges(j)%type, &
146 0 : Vertex%edges(j)%wgtP, Vertex%edges(j)%HeadVertex, &
147 0 : Vertex%edges(j)%TailVertex
148 : enddo
149 : endif
150 :
151 : 95 format(5x,I2,' Member Grid Vertices')
152 : 96 format(5x,I2,' Incident Meta Edges ')
153 : 97 format(10x,'METAEDGE #',I2,2x,'TYPE ',I1,2x,'WGT ',I4,2x,'Processors ',I2,' ---> ',I2)
154 :
155 0 : end subroutine PrintMetaVertex
156 :
157 1536 : subroutine initMetaGraph(ThisProcessorNumber,MetaVertex,GridVertex,GridEdge)
158 : use ll_mod, only : root_t, LLSetEdgeCount, LLFree, LLInsertEdge, LLGetEdgeCount, LLFindEdge
159 : use gridgraph_mod, only : GridEdge_type, printGridVertex
160 : !------------------
161 : !------------------
162 : implicit none
163 :
164 : integer, intent(in) :: ThisProcessorNumber
165 : type (MetaVertex_t), intent(out) :: MetaVertex
166 : type (GridVertex_t), intent(in),target :: GridVertex(:)
167 : type (GridEdge_t), intent(in),target :: GridEdge(:)
168 :
169 : !type (MetaEdge_t), allocatable :: MetaEdge(:)
170 : integer :: nelem,nelem_edge, nedges
171 1536 : integer,allocatable :: icount(:)
172 : integer :: ic,i,j,ii
173 : integer :: npart
174 : integer :: head_processor_number
175 : integer :: tail_processor_number
176 : integer :: nedge_active,enum
177 : logical :: found
178 : integer iTail, iHead, wgtP,wgtS
179 :
180 : type (root_t) :: mEdgeList ! root_t = C++ std::set<std::pair<int,int> >
181 :
182 : logical :: Verbose = .FALSE.
183 : logical :: Debug = .FALSE.
184 :
185 :
186 0 : if(Debug) write(iulog,*)'initMetagraph: point #1'
187 : ! Number of grid vertices
188 1536 : nelem = SIZE(GridVertex)
189 : ! Number of grid edges
190 1536 : nelem_edge = SIZE(GridEdge)
191 :
192 1536 : mEdgeList%number = ThisProcessorNumber
193 1536 : NULLIFY(mEdgeList%first)
194 1536 : call LLSetEdgeCount(0)
195 :
196 66319872 : do i=1,nelem_edge
197 66318336 : tail_processor_number = GridEdge(i)%tail%processor_number
198 66318336 : head_processor_number = GridEdge(i)%head%processor_number
199 66318336 : if(tail_processor_number .eq. ThisProcessorNumber .or. &
200 1536 : head_processor_number .eq. ThisProcessorNumber ) then
201 134168 : call LLInsertEdge(mEdgeList,tail_processor_number,head_processor_number,eNum)
202 : endif
203 : enddo
204 :
205 1536 : call LLGetEdgeCount(nedges)
206 :
207 1536 : NULLIFY(MetaVertex%edges)
208 :
209 4608 : allocate(MetaVertex%edges(nedges))
210 :
211 : ! Initalize the Meta Vertices to zero... probably should be done
212 : ! in a separate routine
213 1536 : MetaVertex%nmembers=0
214 1536 : MetaVertex%number=0
215 1536 : MetaVertex%nedges=0
216 1536 : if(Debug) write(iulog,*)'initMetagraph: point #2'
217 :
218 :
219 : ! Give some identity to the Meta_vertex
220 1536 : MetaVertex%number = ThisProcessorNumber
221 1536 : if(Debug) write(iulog,*)'initMetagraph: point #3'
222 :
223 : ! Look through all the small_vertices and determine the number of
224 : ! member vertices
225 1536 : if(Debug) call PrintGridVertex(GridVertex)
226 1536 : if(Debug) write(iulog,*)'initMetagraph: After call to PrintGridVertex point #3.1'
227 1536 : if(Debug) write(iulog,*)'initMetaGraph: ThisProcessorNumber is ',ThisProcessorNumber
228 :
229 8295936 : do j=1,nelem ! count number of elements on this processor
230 8295936 : if(GridVertex(j)%processor_number .eq. ThisProcessorNumber) then
231 10800 : MetaVertex%nmembers = MetaVertex%nmembers + 1
232 : endif
233 : enddo
234 :
235 1536 : if(Debug) write(iulog,*)'initMetagraph: point #4 '
236 : ! Allocate space for the members of the MetaVertices
237 1536 : if(Debug) write(iulog,*)'initMetagraph: point #4.1 i,MetaVertex%nmembers',i,MetaVertex%nmembers
238 15408 : allocate(MetaVertex%members(MetaVertex%nmembers))
239 :
240 12336 : do j=1, MetaVertex%nmembers
241 12336 : call allocate_gridvertex_nbrs(MetaVertex%members(j))
242 : end do
243 :
244 1536 : if(Debug) write(iulog,*)'initMetagraph: point #5'
245 :
246 : ! Set the identity of the members of the MetaVertices
247 1536 : ic=1
248 8295936 : do j=1,nelem
249 8295936 : if( GridVertex(j)%processor_number .eq. ThisProcessorNumber) then
250 10800 : MetaVertex%members(ic) = GridVertex(j)
251 10800 : ic=ic+1
252 : endif
253 : enddo
254 :
255 1536 : nedges = SIZE(MetaVertex%edges)
256 1536 : if(Debug) write(iulog,*)'initMetagraph: point #6 nedges',nedges
257 : ! Zero out all the edge numbers ... this should probably be
258 : ! move to some initalization routine
259 22856 : MetaVertex%edges%number = 0
260 22856 : MetaVertex%edges%nmembers = 0
261 22856 : MetaVertex%edges%wgtP = 0
262 22856 : MetaVertex%edges%wgtS = 0
263 22856 : MetaVertex%edges%wgtP_ghost = 0
264 22856 : do i=1,nedges
265 22856 : NULLIFY(MetaVertex%edges(i)%members)
266 : enddo
267 :
268 1536 : if(Debug) write(iulog,*)'initMetagraph: point #7'
269 :
270 : ! Insert all the grid edges into the Meta Edges
271 66319872 : do i=1, nelem_edge
272 : ! Which Meta Edge does this grid edge belong
273 66318336 : head_processor_number = GridEdge(i)%head%processor_number
274 66318336 : tail_processor_number = GridEdge(i)%tail%processor_number
275 66318336 : call LLFindEdge(mEdgeList,tail_processor_number,head_processor_number,j,found)
276 66319872 : if(found) then
277 :
278 : ! Increment the number of grid edges contained in the grid edge
279 : ! and setup the pointers
280 134168 : if(Debug) write(iulog,*)'initMetagraph: point #8'
281 134168 : ii=GridEdge(i)%tail_face
282 :
283 134168 : wgtP=Gridedge(i)%tail%nbrs_wgt(ii)
284 134168 : wgtS=1
285 :
286 134168 : MetaVertex%edges(j)%nmembers = MetaVertex%edges(j)%nmembers+1
287 134168 : MetaVertex%edges(j)%wgtP = MetaVertex%edges(j)%wgtP + wgtP
288 134168 : MetaVertex%edges(j)%wgtS = MetaVertex%edges(j)%wgtS + wgtS
289 :
290 134168 : MetaVertex%edges(j)%wgtP_ghost = MetaVertex%edges(j)%wgtP_ghost + Gridedge(i)%tail%nbrs_wgt_ghost(ii)
291 :
292 134168 : if(Debug) write(iulog,*)'initMetagraph: point #9'
293 :
294 : ! If this the first grid edge to be inserted into the Meta Edge
295 : ! do some more stuff
296 :
297 134168 : if(MetaVertex%edges(j)%nmembers .eq. 1) then
298 :
299 21320 : if(Debug) write(iulog,*)'initMetagraph: point #10'
300 21320 : MetaVertex%edges(j)%number = j ! its identity
301 21320 : MetaVertex%edges(j)%type = gridedge_type(GridEdge(i)) ! Type of grid edge
302 :
303 21320 : if(Debug) write(iulog,*)'initMetagraph: point #11'
304 :
305 : ! Setup the pointer to the head and tail of the Vertex
306 21320 : MetaVertex%edges(j)%HeadVertex = head_processor_number
307 21320 : MetaVertex%edges(j)%TailVertex = tail_processor_number
308 21320 : if(Debug) write(iulog,*)'initMetagraph: point #12'
309 :
310 : ! Determine the number of edges for the Meta_Vertex
311 : ! This is the number of processors to communicate with
312 21320 : MetaVertex%nedges = MetaVertex%nedges + 1
313 21320 : if(Debug) write(iulog,*)'initMetagraph: point #13'
314 : endif
315 : endif
316 : enddo
317 :
318 22856 : do i=1,nedges
319 : ! Allocate space for the member edges and edge index
320 198128 : allocate(MetaVertex%edges(i)%members (MetaVertex%edges(i)%nmembers))
321 63960 : allocate(MetaVertex%edges(i)%edgeptrP(MetaVertex%edges(i)%nmembers))
322 63960 : allocate(MetaVertex%edges(i)%edgeptrS(MetaVertex%edges(i)%nmembers))
323 63960 : allocate(MetaVertex%edges(i)%edgeptrP_ghost(MetaVertex%edges(i)%nmembers))
324 155488 : MetaVertex%edges(i)%edgeptrP(:)=0
325 155488 : MetaVertex%edges(i)%edgeptrS(:)=0
326 157024 : MetaVertex%edges(i)%edgeptrP_ghost(:)=0
327 : enddo
328 1536 : if(Debug) write(iulog,*)'initMetagraph: point #14'
329 :
330 : ! Insert the edges into the proper meta edges
331 4608 : allocate(icount(nelem_edge))
332 66319872 : icount=1
333 66319872 : do i=1,nelem_edge
334 66318336 : head_processor_number = GridEdge(i)%head%processor_number
335 66318336 : tail_processor_number = GridEdge(i)%tail%processor_number
336 66318336 : call LLFindEdge(mEdgeList,tail_processor_number,head_processor_number,j,found)
337 66319872 : if(found) then
338 134168 : MetaVertex%edges(j)%members(icount(j)) = GridEdge(i)
339 134168 : if(icount(j)+1 .le. MetaVertex%edges(j)%nmembers) then
340 :
341 112848 : ii=GridEdge(i)%tail_face
342 :
343 112848 : wgtP=Gridedge(i)%tail%nbrs_wgt(ii)
344 112848 : MetaVertex%edges(j)%edgeptrP(icount(j)+1) = MetaVertex%edges(j)%edgeptrP(icount(j)) + wgtP
345 :
346 112848 : wgtS = 1
347 112848 : MetaVertex%edges(j)%edgeptrS(icount(j)+1) = MetaVertex%edges(j)%edgeptrS(icount(j)) + wgtS
348 :
349 112848 : wgtP=Gridedge(i)%tail%nbrs_wgt_ghost(ii)
350 112848 : MetaVertex%edges(j)%edgeptrP_ghost(icount(j)+1) = MetaVertex%edges(j)%edgeptrP_ghost(icount(j)) + wgtP
351 : endif
352 134168 : if(Debug) write(iulog,*)'initMetagraph: point #15'
353 134168 : icount(j)=icount(j)+1
354 : endif
355 : enddo
356 1536 : deallocate(icount)
357 1536 : if(Debug) write(iulog,*)'initMetagraph: point #16'
358 :
359 1536 : if(Verbose) then
360 0 : print *
361 0 : write(iulog,*)"edge bundle list:(INITMETAGRAPH)"
362 0 : call PrintMetaEdge( MetaVertex%edges)
363 0 : write(iulog,*)'initmetagrap: Before last call to PrintMetaVertex'
364 0 : call PrintMetaVertex(MetaVertex)
365 : endif
366 :
367 1536 : call LLFree(mEdgeList)
368 :
369 : 90 format('EDGE #',I2,2x,'TYPE ',I1,2x,'Processor Numbers ',I2,' ---> ',I2)
370 : 100 format(10x,I2,1x,'(',I1,') ---> ',I2,1x,'(',I1,')')
371 :
372 3072 : end subroutine initMetaGraph
373 :
374 :
375 0 : end module metagraph_mod
|