metis_mod.F90 coverage: 0.00 %func 0.00 %block
1) #ifdef HAVE_CONFIG_H
2) #include "config.h"
3) #endif
4)
5) module metis_mod
6) use kinds, only : iulog
7) use parallel_mod, only : abortmp
8) implicit none
9) private
10) integer, parameter :: VertexWeight = 1000
11) integer, parameter :: EdgeWeight = 1000
12)
13) public :: genmetispart
14) private :: PartitionGraph
15) private :: CreateMeshGraph
16) private :: PrintMeshGraph
17) private :: genLocal2Global
18) private :: sort
19) contains
20)
21) subroutine genmetispart(GridEdge, GridVertex)
22) use gridgraph_mod, only : GridVertex_t, GridEdge_t, freegraph, createsubgridgraph, printgridvertex
23) use kinds, only : int_kind
24) use parallel_mod , only : iam, FrameWeight, PartitionForNodes
25) use dimensions_mod , only : nmpi_per_node, npart, nnodes, nelem
26) use control_mod, only: partmethod
27) use params_mod, only : wrecursive
28)
29) implicit none
30)
31) type (GridVertex_t), intent(inout) :: GridVertex(:)
32) type (GridEdge_t), intent(inout) :: GridEdge(:)
33)
34) integer , target, allocatable :: xadj(:),adjncy(:) ! Adjacency structure for METIS
35) integer , target, allocatable :: vwgt(:),adjwgt(:) ! Weights for the adj struct for METIS
36)
37) integer , target, allocatable :: xadj_nl(:),adjncy_nl(:)
38) integer , target, allocatable :: vwgt_nl(:),adjwgt_nl(:)
39)
40) type (GridVertex_t), allocatable :: SubVertex(:)
41)
42) real(kind=4), allocatable :: tpwgts(:)
43) real(kind=4) :: tmp
44)
45) integer(kind=int_kind), allocatable :: part(:)
46) integer, allocatable :: part_nl(:),local2global_nl(:)
47) integer, allocatable :: part_fl(:),local2global_fl(:)
48)
49)
50) integer, allocatable :: cnt(:),newnum(:),oldnum(:)
51)
52) integer :: nelem_edge,numflag,edgecut,wgtflag
53) integer :: head_part,tail_part
54) integer :: options(5)
55) integer :: i,j,ii,ig,in,ip,if
56) integer :: nelem_nl,nelem_fl,newPartition
57) integer :: partitionmethod,numpartitions
58) integer :: nodes_per_frame
59) logical , parameter :: Debug = .true.
60) real (kind=4) :: dummy(1)
61) nelem_edge = SIZE(GridEdge)
62)
63) print *, "nelem = ", nelem
64) print *, "nelem_edge = ", nelem_edge
65) allocate(tpwgts(npart))
66) allocate(part(nelem))
67) allocate(xadj(nelem+1))
68) allocate(vwgt(nelem))
69) allocate(adjncy(nelem_edge))
70) allocate(adjwgt(nelem_edge))
71)
72) if(Debug) write(iulog,*)'genmetispart: point #1'
73) ! =============================================
74) ! Generate Graph for METIS
75) ! =============================================
76) !DBG call PrintGridVertex(GridVertex)
77) call CreateMeshGraph(GridVertex,xadj,adjncy,adjwgt)
78) ! Add weights to all the vertices
79) vwgt(:)=VertexWeight
80) if(Debug) write(iulog,*)'genmetispart: point #2'
81)
82) ! ====================================================
83) ! Some cruff from the Weighted partitioning
84) ! experimentation... Remove soon
85) ! ====================================================
86) do i=1,npart
87) tpwgts(i) = i
88) enddo
89) tmp = SUM(tpwgts)
90) tpwgts = tpwgts/tmp
91) if(Debug) write(iulog,*)'genmetispart: point #3'
92)
93) !================================================
94) ! Setup some options for the METIS partitioning
95) !================================================
96) options(1) = 0 ! Use Default METIS options
97) ! wgtflag = 1 ! Weights on the Edges only
98) wgtflag = 3 ! Weights on the edges and vertices
99) numflag = 1 ! Use Fortran based numbering
100)
101) if(npart > 1) then
102)
103) if (PartitionForNodes) then
104) numPartitions = nnodes
105) PartitionMethod = partmethod
106) else
107) numPartitions = npart
108) PartitionMethod = partmethod
109) endif
110) !===========================================
111) ! Generate the METIS partitioning
112) !===========================================
113) if(iam .eq. 1) write(6,100) nelem,PartitionMethod,numPartitions
114) if(Debug) write(iulog,*)'genmetispart: point #4'
115) if (PartitionMethod==WRECURSIVE) then
116) call PartitionGraph(PartitionMethod,nelem,xadj,adjncy, &
117) vwgt,adjwgt,wgtflag,numflag, &
118) numPartitions,FrameWeight,options,edgecut,part)
119) else
120) ! FrameWeight has not been allocated in this case, so replace
121) ! with dummy argument:
122) call PartitionGraph(PartitionMethod,nelem,xadj,adjncy, &
123) vwgt,adjwgt,wgtflag,numflag, &
124) numPartitions,dummy,options,edgecut,part)
125) endif
126) if(Debug) write(iulog,*)'genmetispart: point #5'
127) ! ===============================
128) ! Do not partitiion for nodes
129) ! ===============================
130) ! PartitionForNodes=.FALSE.
131) if(PartitionForNodes) then
132) ! =========================================
133) ! Partition the graph on each compute node
134) ! =========================================
135) if(nmpi_per_node .ne. 1 ) then
136)
137) do ip = 1,nelem
138) part(ip) = nmpi_per_node*(part(ip)-1) + 1
139) enddo
140)
141) if(Debug) write(iulog,*)'genmetispart: point #14'
142)
143) ! ====================
144) ! Loop over each node
145) ! ====================
146) do in = 1,nnodes
147)
148) newPartition = nmpi_per_node*(in -1) + 1
149) nelem_nl = COUNT(part .eq. newPartition)
150)
151) ! ===================================================
152) ! Allocate all the memory for node level partitioning
153) ! ===================================================
154)
155) allocate(local2global_nl(nelem_nl))
156) allocate(part_nl(nelem_nl))
157)
158) if(Debug) write(iulog,*)'genmetispart: point #15'
159) allocate(xadj_nl(nelem_nl+1))
160) allocate(vwgt_nl(nelem_nl))
161) allocate(adjncy_nl(8*nelem_nl))
162) allocate(adjwgt_nl(8*nelem_nl))
163) adjncy_nl(:) = 0
164) adjwgt_nl(:) = 0
165) part_nl(:) = 0
166)
167) if(Debug) write(iulog,*)'genmetispart: point #16'
168) !Add vertex weights to the subgraphs
169) vwgt_nl(:) = VertexWeight
170) allocate(SubVertex(nelem_nl))
171)
172) ! =====================================
173) ! Setup the index translation arrays
174) ! =====================================
175) call genLocal2Global(local2global_nl,part,newPartition)
176)
177) if(Debug) write(iulog,*)'genmetispart: point #17'
178) ! ======================================================
179) ! Create a set of the Vertices that represent a subgraph
180) ! ======================================================
181) if(Debug) write(iulog,*)'genmetispart: local2global_nl',local2global_nl
182) if(Debug) call PrintGridVertex(GridVertex)
183) call CreateSubGridGraph(GridVertex,SubVertex,local2global_nl)
184) !JMD call CheckGridNeighbors(SubVertex)
185) if(Debug) call PrintGridVertex(SubVertex)
186) if(Debug) write(iulog,*)'genmetispart: point #18'
187)
188) ! ==============================
189) ! Convert grep to Metis format
190) ! ==============================
191) call CreateMeshGraph(SubVertex,xadj_nl,adjncy_nl,adjwgt_nl)
192)
193) !debug call PrintMetisgraph(xadj_nl,adjncy_nl,adjwgt_nl)
194) ! =======================
195) ! Partition the subgraph
196) ! =======================
197) if(iam .eq. 1) write(6,100) nelem_nl,partmethod,nmpi_per_node
198) call PartitionGraph(partmethod,nelem_nl,xadj_nl,adjncy_nl, &
199) vwgt_nl,adjwgt_nl,wgtflag,numflag, &
200) nmpi_per_node,tpwgts,options,edgecut,part_nl)
201) if(Debug) write(iulog,*)'genmetispart: point #19'
202)
203) ! =========================================================
204) ! Apply the node partitioning the the overall partitioning
205) ! =========================================================
206) do i=1,nelem_nl
207) ig = local2global_nl(i)
208) part(ig) = part(ig) + (part_nl(i)-1)
209) enddo
210)
211) ! ======================================
212) ! Free up the temporaries that were used
213) ! ======================================
214) deallocate(local2global_nl)
215) deallocate(part_nl)
216) deallocate(xadj_nl)
217) deallocate(vwgt_nl)
218) deallocate(adjncy_nl)
219) deallocate(adjwgt_nl)
220)
221)
222) if(Debug) write(iulog,*)'genmetispart: point #20'
223) call FreeGraph(SubVertex)
224) deallocate(SubVertex)
225) if(Debug) write(iulog,*)'genmetispart: point #21'
226)
227) enddo
228)
229)
230) endif ! if node based partitioning
231)
232) endif ! if multilevel
233)
234) else ! else no partitioning needed
235) !======================================================================
236) ! It appears that Metis will set part(:) = 2 if nnodes == 1,
237) ! which messes stuff up, so set it myself
238) !======================================================================
239) part(:)=1
240) endif
241)
242) if(Debug) write(iulog,*)'genmetispart: point #22'
243) ! ===========================================
244) ! Add the partitioning information into the
245) ! Grid Vertex and Grid Edge structures
246) ! ===========================================
247)
248) GridVertex(:)%processor_number = part
249) if(Debug) write(iulog,*)'genmetispart: point #23'
250) !=================================================
251) ! Output the partitioning information
252) !=================================================
253) #if 0
254)
255) write(iulog,*)'Metis Parititioning: '
256) write(iulog,*)part
257)
258) write(iulog,*)'Metis Parititioning: '
259) do k=1,nelem
260) write(iulog,*) k,part(k)
261) enddo
262) call abortmp(' at the end of genmetispart')
263)
264) if(OutputFiles) then
265) do k=1,nelem
266) write(10,*) GridVertex(k)%processor_number
267) end do
268) close(10)
269) end if
270) #endif
271)
272) 100 format("genmetispart: Partitioning ",i4," elements, using method: ",i2," into ",i4," pieces")
273)
274) end subroutine genmetispart
275)
276) subroutine NewPartitionNumber(newnum,oldnum,cnt)
277)
278) implicit none
279)
280) integer, intent(in) :: oldnum(:)
281) integer, intent(in) :: cnt(:)
282) integer, intent(out) :: newnum(:)
283)
284) integer :: i,n
285)
286) n = SIZE(cnt)
287)
288) newnum(1) = oldnum(1)
289) do i=2,n
290) newnum(i) = newnum(i-1) + cnt(i-1)
291) enddo
292)
293) end subroutine NewPartitionNumber
294)
295) subroutine genLocal2Global(local2Global,part,newP)
296)
297) implicit none
298)
299) integer,intent(inout) :: local2Global(:)
300) integer,intent(in) :: part(:)
301) integer,intent(in) :: newP
302)
303) integer :: i,ii,nelem
304)
305) nelem = SIZE(part)
306) ii = 1
307) do i=1,nelem
308) if( part(i) .eq. newP) then
309) local2global(ii) = i
310) ii = ii + 1
311) endif
312) enddo
313)
314) end subroutine genLocal2Global
315)
316) subroutine PartitionGraph(partmethod,nelem,xadj,adjncy,vwgt,adjwgt, &
317) wgtflag,numflag,npart,tpwgts,options,edgecut,part)
318)
319) implicit none
320)
321) integer :: partmethod,nelem
322) integer :: xadj(:),adjncy(:)
323) integer :: vwgt(:),adjwgt(:)
324) integer :: wgtflag,numflag,edgecut
325) real(kind=4) :: tpwgts(:)
326) integer :: options(5)
327) integer :: npart
328) integer :: part(:)
329)
330)
331) #ifdef _USEMETIS
332)
333) if (partmethod .eq. KWAY) then
334) call metis_partgraphkway(nelem,xadj,adjncy, &
335) vwgt,adjwgt,wgtflag,numflag,npart,options,edgecut,part)
336) else if (partmethod .eq. RECURSIVE) then
337) call METIS_PartGraphRecursive(nelem,xadj,adjncy, &
338) vwgt,adjwgt,wgtflag,numflag,npart,options,edgecut,part)
339) else if (partmethod .eq. WRECURSIVE) then
340) call METIS_WPartGraphRecursive(nelem,xadj,adjncy, &
341) vwgt,adjwgt,wgtflag,numflag,npart,tpwgts,options,edgecut,part)
342) else if (partmethod .eq. VOLUME) then
343) call METIS_PartGraphVKway(nelem,xadj,adjncy, &
344) vwgt,adjwgt,wgtflag,numflag,npart,options,edgecut,part)
345) endif
346) #endif
347)
348) end subroutine PartitionGraph
349)
350) subroutine CreateMeshGraph(GridVertex,xadj,adjncy,adjwgt)
351) use gridgraph_mod, only : GridVertex_t, num_neighbors
352) use kinds, only : real_kind, int_kind
353) type (GridVertex_t), intent(in) :: GridVertex(:)
354) integer,intent(out) :: xadj(:),adjncy(:),adjwgt(:)
355)
356) integer :: i,j,k,ii,jj
357) integer :: degree,nelem
358) !integer(kind=int_kind),allocatable :: neigh_list(:),sort_indices(:)
359) !real(kind=REAL_KIND),allocatable :: neigh_wgt(:)
360) integer(kind=int_kind) :: neigh_list(num_neighbors), &
361) sort_indices(num_neighbors)
362) real(kind=REAL_KIND) :: neigh_wgt(num_neighbors)
363) integer :: max_neigh
364)
365) integer :: start, cnt
366)
367) nelem = SIZE(GridVertex)
368)
369) degree = 0
370) ii = 1
371)
372) max_neigh = num_neighbors
373)
374) do i=1,nelem
375) print *, "i = ", i
376) xadj(i) = ii
377) degree = 0
378) neigh_list=0
379)
380)
381) do j=1,num_neighbors
382) cnt = GridVertex(i)%nbrs_ptr(j+1) - GridVertex(i)%nbrs_ptr(j)
383) start = GridVertex(i)%nbrs_ptr(j)
384) print *, "j,cnt,start = ", j, cnt, start
385) do k=0, cnt-1
386) print *, "k,wgt = ", k, GridVertex(i)%nbrs_wgt(start+k)
387) if(GridVertex(i)%nbrs_wgt(start+k) .gt. 0) then
388) degree = degree + 1
389) !adjncy(ii+degree-1) = GridVertex(i)%nbrs(start+k)
390) neigh_list(degree) = GridVertex(i)%nbrs(start+k)
391) !adjwgt(ii+degree-1) = GridVertex(i)%nbrs_wgt(start+k)*EdgeWeight
392) neigh_wgt(degree) = GridVertex(i)%nbrs_wgt(start+k)*EdgeWeight
393) endif
394) enddo
395) enddo
396) if (degree > max_neigh) call abortmp( "number of neighbors found exceeds expected max")
397)
398) call sort(degree,neigh_list,sort_indices)
399) !degree = COUNT(GridVertex(i)%nbrs_wgt(:) .gt. 0)
400) ! Copy the sorted adjncy list in
401) do j=1,degree
402) adjncy(ii+j-1) = neigh_list(sort_indices(j))
403) adjwgt(ii+j-1) = neigh_wgt(sort_indices(j))
404) enddo
405) print *, "degree = ", degree
406) !print *, "wgts = ", GridVertex(i)%nbrs_wgt(:)
407) ii = ii + degree
408) enddo
409) xadj(nelem+1) = ii
410) !open(unit=11, file="csr.txt")
411) !write(11,*) xadj
412) !write(11,*) adjncy
413) !write(11,*) adjwgt
414) !close(11)
415)
416) end subroutine CreateMeshGraph
417)
418) subroutine sort(n,list,index)
419) use kinds, only : int_kind
420) implicit none
421) integer, intent(in) :: n
422) integer(kind=int_kind), intent(in) :: list(n)
423) integer(kind=int_kind), intent(inout) :: index(n)
424)
425) ! Local variables
426) integer :: i,iloc,nct,ii,lmax,lmin
427) logical :: msk(n)
428) logical,parameter :: Debug =.FALSE.
429)
430) msk=.TRUE.
431) if(Debug) write(iulog,*)'sort: point #1'
432) nct=0
433) do i=1,n
434) if(list(i) .eq. 0) then
435) msk(i)=.FALSE.
436) else
437) nct = nct + 1
438) endif
439) enddo
440) if(Debug) write(iulog,*)'sort: point #2 list: ',list
441) if(Debug) write(iulog,*)'sort: point #2.1 msk: ',msk
442)
443) lmax=maxval(list)
444) do i=1,nct
445) if(Debug) write(iulog,*)'sort: point #3'
446)
447) ! pgf90 Rel 3.1-4i: minloc() with mask is buggy
448) ! iloc = minloc(list,dim=1,mask=msk)
449) lmin=lmax
450) iloc=-1
451) do ii=1,n
452) if (msk(ii) .and. list(ii)<= lmin) iloc=ii
453) enddo
454) if (iloc==-1) call abortmp( "sort() error")
455)
456) index(i)=iloc
457) if(Debug) write(iulog,*)'sort: point #4'
458) msk(iloc)=.FALSE.
459) if(Debug) write(iulog,*)'sort: point #5'
460) if(Debug) write(iulog,*)'sort: i, msk',i,msk
461) enddo
462) if(Debug) write(iulog,*)'sort: point #6'
463) !DBG write(iulog,*)'sort: list is:',list
464) !DBG write(iulog,*)'sort: index is:',index
465) !DBG stop
466)
467) end subroutine sort
468) !------------------------------------------------------------------------
469) subroutine PrintMeshGraph(xadj,adjncy,adjwgt)
470) integer, intent(in) :: xadj(:),adjncy(:),adjwgt(:)
471) integer :: nelem,nadj,i,j
472)
473) nelem = SIZE(xadj)
474)
475) do i=1,nelem
476) write(iulog,*)'xadj(i):= ',xadj(i)
477) enddo
478) nadj = xadj(nelem)
479) do j=1,nadj-1
480) write(iulog,*)'{adjncy(i),adjwgt(i)}:=',adjncy(j),adjwgt(j)
481) enddo
482)
483) end subroutine PrintMeshGraph
484) !------------------------------------------------------------------------
485) end module metis_mod
