|
DOUG 0.2
|
00001 ! DOUG - Domain decomposition On Unstructured Grids 00002 ! Copyright (C) 1998-2006 Faculty of Computer Science, University of Tartu and 00003 ! Department of Mathematics, University of Bath 00004 ! 00005 ! This library is free software; you can redistribute it and/or 00006 ! modify it under the terms of the GNU Lesser General Public 00007 ! License as published by the Free Software Foundation; either 00008 ! version 2.1 of the License, or (at your option) any later version. 00009 ! 00010 ! This library is distributed in the hope that it will be useful, 00011 ! but WITHOUT ANY WARRANTY; without even the implied warranty of 00012 ! MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU 00013 ! Lesser General Public License for more details. 00014 ! 00015 ! You should have received a copy of the GNU Lesser General Public 00016 ! License along with this library; if not, write to the Free Software 00017 ! Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA 00018 ! or contact the authors (University of Tartu, Faculty of Computer Science, Chair 00019 ! of Distributed Systems, Liivi 2, 50409 Tartu, Estonia, http://dougdevel.org, 00020 ! mailto:info(at)dougdevel.org) 00021 00022 module Graph_class 00023 00024 use DOUG_utils 00025 00026 implicit none 00027 00028 integer, parameter :: D_GRAPH_DUAL = 1 00029 integer, parameter :: D_GRAPH_NODAL = 2 00030 00031 !-------------------------------------------------------------- 00032 ! Graph type (with parameters for partitioning and its result 00033 ! (although one may wish to separate them I tend 00034 ! to keep them together)). 00035 ! Adjacency structure of the sparse graph is stored using 00036 ! the compressed storage format (CSR). Structure of the 00037 ! graph with 'nvtx' vertices and 'nedges' edges is represented using 00038 ! two arrays 'xadj' and 'adjncy'. The 'xadj' array is of size 00039 ! n+1 whereas 'adjncy' array is of size 2*m. 00040 ! 00041 ! data fields: 00042 ! wgtflag -- type of graph 00043 ! type vertices weighted? edges weighted? 00044 ! 0 no no 00045 ! 1 yes no 00046 ! 2 no yes 00047 ! 3 yes yes 00048 ! 00049 !-------------------------------------------------------------- 00050 type Graph 00051 ! Number of vertices in the graph 00052 integer :: nvtx = -1 00053 ! Number of edges in the graph 00054 integer :: nedges = -1 00055 integer, dimension(:), pointer :: xadj ! size: nvtx+1 00056 integer, dimension(:), pointer :: adjncy ! size: 2*nedges 00057 ! Type of the graph : 00058 ! dual: D_GRAPH_DUAL = 1, nodal: D_GRAPH_NODAL = 2 00059 integer :: type = -1 00060 00061 ! Partitioning: 00062 logical :: parted = .false. 00063 real(kind=rk), dimension(:), pointer :: vwgt => NULL() 00064 real(kind=rk), dimension(:), pointer :: adjwgt => NULL() 00065 integer :: wgtflag = 0 00066 integer :: numflag = 1 ! one based Fortran 00067 ! Number of parts to partition the graph 00068 integer :: nparts = 1 00069 !integer, dimension(5) :: options = (/0, 3, 1, 1, 0/) 00070 integer :: edgecut 00071 integer, dimension(:), pointer :: part 00072 end type Graph 00073 00074 ! METIS_PartGraphRecursive() 00075 ! Objective: minimize the edgecut. 00076 ! (preferable to partition into smaller than 8 number of partitions) 00077 integer, parameter :: D_PART_PMETIS = 1 00078 ! METIS_PartGraphKway() 00079 ! Objective: minimize the edgecut. 00080 ! (use to partition in 8 or more parts) 00081 integer, parameter :: D_PART_KMETIS = 2 00082 ! METIS_PartGraphVKway() 00083 ! Objective: minimize the total communication volume 00084 integer, parameter :: D_PART_VKMETIS = 3 00085 00086 ! Mesh partitioning: 00087 integer, parameter :: D_PART_OMETIS = 4 ! 00088 00089 00090 ! Private methods 00091 private :: & 00092 Graph_partng 00093 00094 contains 00095 00096 !------------------------------------------ 00097 ! Graph_new() 00098 ! If arguments are given the space will be allocated: 00099 ! n - Number of vertices of the graph 00100 ! m - Number of edges of the graph 00101 ! Result: Graph 00102 !------------------------------------------ 00103 function Graph_New(n, m, type) result(G) 00104 implicit none 00105 00106 integer, intent(in), optional :: n, m, type 00107 type(Graph) :: G ! Graph (not filled, just allocated 00108 ! if n and m are present) 00109 00110 if (present(n).and.present(m)) then 00111 allocate(G%xadj(n+1), G%adjncy(2*m)) 00112 G%nvtx = n 00113 G%nedges = m 00114 else 00115 G%nvtx = 0 00116 G%nedges = 0 00117 G%xadj => NULL() 00118 G%adjncy => NULL() 00119 end if 00120 00121 G%wgtflag = 0 00122 G%numflag = 1 00123 G%vwgt => NULL() 00124 G%adjwgt => NULL() 00125 G%part => NULL() 00126 00127 if (present(type)) then 00128 if ((type == D_GRAPH_DUAL).or.(type == D_GRAPH_NODAL)) then 00129 G%type = type 00130 else 00131 call DOUG_abort('[Graph_New] : Wrong graph type given.',-1) 00132 end if 00133 else 00134 ! Default value 00135 write(stream,*) 'Type of Graph set to default D_GRAPH_DUAL' 00136 G%type = D_GRAPH_DUAL 00137 end if 00138 end function Graph_New 00139 00140 00141 !----------------------------------------------- 00142 ! Graph_newInit() 00143 !----------------------------------------------- 00144 function Graph_newInit(n, m, xadj, adjncy, type) result(G) 00145 00146 implicit none 00147 00148 integer, intent(in) :: n, m 00149 integer, dimension(:), intent(in) :: xadj 00150 integer, dimension(:), intent(in) :: adjncy 00151 integer, intent(in), optional :: type 00152 type(Graph) :: G ! Filled Graph 00153 00154 ! Check for consistancy 00155 if (size(xadj,1) /= n+1) & 00156 call DOUG_abort('[Graph_Init] : size(xadj) /= nell+1') 00157 if (size(adjncy,1) /= 2*m) & 00158 call DOUG_abort('[Graph_Init] : size(adjncy,1) /= 2*(graph edges)') 00159 00160 G = Graph_New(n, m, type) 00161 G%xadj = xadj 00162 G%adjncy = adjncy 00163 00164 end function Graph_newInit 00165 00166 00167 !---------------------------------------------- 00168 ! Graph_Init() 00169 !---------------------------------------------- 00170 subroutine Graph_Init(G, xadj, adjncy, type) 00171 type(Graph), intent(in out) :: G 00172 integer, dimension(:), intent(in) :: xadj 00173 integer, dimension(:), intent(in) :: adjncy 00174 integer, intent(in), optional :: type 00175 00176 ! Check for consistancy 00177 if (size(xadj,1) /= G%nvtx+1) & 00178 call DOUG_abort('[Graph_Init] : size(xadj) /= nverts+1',-1) 00179 if (size(adjncy,1) /= 2*G%nedges) & 00180 call DOUG_abort('[Graph_Init] : size(adjncy,1) /= 2*(graph edges)',-1) 00181 00182 if (.not.associated(G%xadj)) allocate(G%xadj(G%nvtx+1)) 00183 if (.not.associated(G%adjncy)) allocate(G%adjncy(2*G%nedges)) 00184 G%xadj = xadj 00185 G%adjncy = adjncy 00186 00187 if (present(type)) then 00188 if ((type == D_GRAPH_DUAL).or.(type == D_GRAPH_NODAL)) then 00189 G%type = type 00190 else 00191 call DOUG_abort('[Graph_New] : Wrong graph type given.',-1) 00192 end if 00193 else 00194 ! Default value 00195 write(stream,*) 'Type of Graph set to default D_GRAPH_DUAL' 00196 G%type = D_GRAPH_DUAL 00197 end if 00198 end subroutine Graph_Init 00199 00200 00201 !---------------------------------------------- 00202 ! Graph_Destroy() - Graph destrucor 00203 ! Arguments: 00204 ! G - Graph 00205 ! Result: dealocated arrays and zeroed fields 00206 !---------------------------------------------- 00207 subroutine Graph_Destroy(G) 00208 implicit none 00209 00210 type(Graph), intent(in out) :: G ! Graph 00211 00212 if (associated(G%xadj)) deallocate(G%xadj) 00213 if (associated(G%adjncy)) deallocate(G%adjncy) 00214 00215 G%nvtx = 0 00216 G%nedges = 0 00217 00218 ! And much more here for partitioning... 00219 if (associated(G%part)) deallocate(G%part) 00220 00221 end subroutine Graph_Destroy 00222 00223 00224 !------------------------------------------------------------ 00225 ! Graph_Partition() - graph partitioning (interface to METIS) 00226 !------------------------------------------------------------ 00227 subroutine Graph_Partition(G, nparts, method, options) 00228 00229 use globals, only: stream 00230 00231 implicit none 00232 00233 !!$ include 'globals_partng.F90' 00234 00235 type(Graph), intent(in out) :: G 00236 integer, intent(in) :: nparts 00237 integer, optional :: method 00238 integer, dimension(5), optional :: options 00239 00240 ! Default: recursive graph partitionig 00241 ! pmetis() -> METIS_PartGraphRecursive() 00242 integer :: part_method = 0 00243 integer, dimension(5) :: part_options = (/0, 3, 2, 2, 0/) 00244 ! Defaults: 00245 ! ... 00246 00247 ! Check for arguments. 00248 if (present(method)) then 00249 part_method = method 00250 else 00251 ! Switch to Kway-partitioning if number of desired 00252 ! partitions is greater than 8. (Significantly speeds up.) 00253 if (G%nparts > 8) part_method = D_PART_KMETIS 00254 end if 00255 if (present(options)) part_options = options 00256 00257 ! Call partitioner wrapper 00258 call Graph_partng(G, nparts, part_method, part_options) 00259 00260 end subroutine Graph_Partition 00261 00262 00263 !-------------------------------------------------- 00264 ! Graph_partng() 00265 !-------------------------------------------------- 00266 subroutine Graph_partng(G, nparts, method, options) 00267 00268 !use globals_partng 00269 use globals, only : stream 00270 00271 implicit none 00272 00273 !!$ include 'globals_partng.F90' 00274 00275 type(Graph), intent(in out) :: G 00276 integer, intent(in) :: nparts 00277 integer, intent(in) :: method 00278 integer, dimension(5), intent(in) :: options 00279 00280 write(stream, '(/a)', advance='no') ' Graph partitioning: ' 00281 00282 G%nparts = nparts 00283 allocate(G%part(G%nvtx)) 00284 00285 if (nparts > 1) then ! METIS dislikes partitioning into one partition 00286 00287 ! Select appropriate method. 00288 select case(method) 00289 case (D_PART_PMETIS) 00290 00291 ! Multilevel recursive bisection 00292 write(stream, *) 'multilevel recursive bisection' 00293 00294 !metis_partgraphrecursive 00295 call METIS_PartGraphRecursive(G%nvtx, G%xadj, G%adjncy, & 00296 G%vwgt, G%adjwgt, G%wgtflag, G%numflag, & 00297 nparts, options, G%edgecut, G%part) 00298 00299 case (D_PART_KMETIS) 00300 00301 ! Multilevel K-way partitioning 00302 write(stream, *) 'multilevel K-way partitioning' 00303 00304 call METIS_PartGraphKway(G%nvtx, G%xadj, G%adjncy, & 00305 G%vwgt, G%adjwgt, G%wgtflag, G%numflag, & 00306 nparts, options, G%edgecut, G%part) 00307 00308 case (D_PART_VKMETIS) 00309 00310 ! Multilevel K-way (min. communication) 00311 write(stream, *) 'multilevel K-way (min. communication)' 00312 call METIS_PartGraphVKway(G%nvtx, G%xadj, G%adjncy, & 00313 !G%vwgt, G%adjwgt, G%wgtflag, G%numflag, & ! currently gfortran fails on G%vwgt=>NULL() and G%adjwgt=>NULL(), compiler bug? 00314 NULL(), NULL(), G%wgtflag, G%numflag, & 00315 nparts, options, G%edgecut, G%part) 00316 00317 case default 00318 call DOUG_abort('[Graph_partng] : Wrong Graph partitioning'//& 00319 ' method specified') 00320 end select 00321 00322 else 00323 G%nparts = 1 00324 G%part = 1 00325 end if 00326 00327 G%parted = .true. 00328 00329 end subroutine Graph_partng 00330 00331 end module Graph_class
1.7.3-20110217