1! 2! Copyright (C) by Argonne National Laboratory 3! See COPYRIGHT in top-level directory 4! 5 6! This program is Fortran version of dgraph_unwgt.c 7! Specify a distributed graph of a bidirectional ring of the MPI_COMM_WORLD, 8! i.e. everyone only talks to left and right neighbors. 9 10 logical function validate_dgraph(dgraph_comm) 11 use mpi_f08 12 13 type(MPI_Comm) dgraph_comm 14 integer comm_topo 15 integer src_sz, dest_sz 16 integer ierr; 17 logical wgt_flag; 18 integer srcs(2), dests(2) 19 20 integer world_rank, world_size; 21 integer idx, nbr_sep 22 23 comm_topo = MPI_UNDEFINED 24 call MPI_Topo_test(dgraph_comm, comm_topo, ierr); 25 if (comm_topo .ne. MPI_DIST_GRAPH) then 26 validate_dgraph = .false. 27 write(6,*) "dgraph_comm is NOT of type MPI_DIST_GRAPH." 28 return 29 endif 30 31 call MPI_Dist_graph_neighbors_count(dgraph_comm, & 32 & src_sz, dest_sz, wgt_flag, & 33 & ierr) 34 if (ierr .ne. MPI_SUCCESS) then 35 validate_dgraph = .false. 36 write(6,*) "MPI_Dist_graph_neighbors_count() fails!" 37 return 38 endif 39 if (wgt_flag) then 40 validate_dgraph = .false. 41 write(6,*) "dgraph_comm is NOT created with MPI_UNWEIGHTED." 42 return 43 endif 44 45 if (src_sz .ne. 2 .or. dest_sz .ne. 2) then 46 validate_dgraph = .false. 47 write(6,*) "source or destination edge array is not size 2." 48 write(6,"('src_sz = ',I3,', dest_sz = ',I3)") src_sz, dest_sz 49 return 50 endif 51 52 call MPI_Dist_graph_neighbors(dgraph_comm, & 53 & src_sz, srcs, MPI_UNWEIGHTED, & 54 & dest_sz, dests, MPI_UNWEIGHTED, & 55 & ierr) 56 if (ierr .ne. MPI_SUCCESS) then 57 validate_dgraph = .false. 58 write(6,*) "MPI_Dist_graph_neighbors() fails!" 59 return 60 endif 61 62! Check if the neighbors returned from MPI are really 63! the nearest neighbors that within a ring. 64 call MPI_Comm_rank(MPI_COMM_WORLD, world_rank, ierr) 65 call MPI_Comm_size(MPI_COMM_WORLD, world_size, ierr) 66 67 do idx = 1, src_sz 68 nbr_sep = iabs(srcs(idx) - world_rank) 69 if (nbr_sep .ne. 1 .and. nbr_sep .ne. (world_size-1)) then 70 validate_dgraph = .false. 71 write(6,"('srcs[',I3,']=',I3, & 72 & ' is NOT a neighbor of my rank',I3)") & 73 & idx, srcs(idx), world_rank 74 return 75 endif 76 enddo 77 do idx = 1, dest_sz 78 nbr_sep = iabs(dests(idx) - world_rank) 79 if (nbr_sep .ne. 1 .and. nbr_sep .ne. (world_size-1)) then 80 validate_dgraph = .false. 81 write(6,"('dests[',I3,']=',I3, & 82 & ' is NOT a neighbor of my rank',I3)") & 83 & idx, dests(idx), world_rank 84 return 85 endif 86 enddo 87 88 validate_dgraph = .true. 89 return 90 end 91 92 integer function ring_rank(world_size, in_rank) 93 integer world_size, in_rank 94 if (in_rank .ge. 0 .and. in_rank .lt. world_size) then 95 ring_rank = in_rank 96 return 97 endif 98 if (in_rank .lt. 0 ) then 99 ring_rank = in_rank + world_size 100 return 101 endif 102 if (in_rank .ge. world_size) then 103 ring_rank = in_rank - world_size 104 return 105 endif 106 ring_rank = -99999 107 return 108 end 109 110 111 112 program dgraph_unwgt 113 use mpi_f08 114 115 integer ring_rank 116 external ring_rank 117 logical validate_dgraph 118 external validate_dgraph 119 integer errs, ierr 120 121 type(MPI_Comm) dgraph_comm 122 integer world_size, world_rank 123 integer src_sz, dest_sz 124 integer degs(1) 125 integer srcs(2), dests(2) 126 127 errs = 0 128 call MTEST_Init(ierr) 129 call MPI_Comm_rank(MPI_COMM_WORLD, world_rank, ierr) 130 call MPI_Comm_size(MPI_COMM_WORLD, world_size, ierr) 131 132 srcs(1) = world_rank 133 degs(1) = 2; 134 dests(1) = ring_rank(world_size, world_rank-1) 135 dests(2) = ring_rank(world_size, world_rank+1) 136 call MPI_Dist_graph_create(MPI_COMM_WORLD, 1, srcs, degs, dests, & 137 & MPI_UNWEIGHTED, MPI_INFO_NULL, & 138 & .true., dgraph_comm, ierr) 139 if (ierr .ne. MPI_SUCCESS) then 140 write(6,*) "MPI_Dist_graph_create() fails!" 141 call MPI_Abort(MPI_COMM_WORLD, 1, ierr) 142 stop 143 endif 144 if (.not. validate_dgraph(dgraph_comm)) then 145 write(6,*) "MPI_Dist_graph_create() does not create" & 146 & //"a bidirectional ring graph!" 147 call MPI_Abort(MPI_COMM_WORLD, 1, ierr) 148 stop 149 endif 150 call MPI_Comm_free(dgraph_comm, ierr) 151 152! now create one with MPI_WEIGHTS_EMPTY 153! NOTE that MPI_WEIGHTS_EMPTY was added in MPI-3 and does not 154! appear before then. Incluing this test means that this test cannot 155! be compiled if the MPI version is less than 3 (see the testlist file) 156 157 degs(1) = 0; 158 call MPI_Dist_graph_create(MPI_COMM_WORLD, 1, srcs, degs, dests, & 159 & MPI_WEIGHTS_EMPTY, MPI_INFO_NULL, & 160 & .true., dgraph_comm, ierr) 161 if (ierr .ne. MPI_SUCCESS) then 162 write(6,*) "MPI_Dist_graph_create() fails!" 163 call MPI_Abort(MPI_COMM_WORLD, 1, ierr) 164 stop 165 endif 166 call MPI_Comm_free(dgraph_comm, ierr) 167 168 src_sz = 2 169 srcs(1) = ring_rank(world_size, world_rank-1) 170 srcs(2) = ring_rank(world_size, world_rank+1) 171 dest_sz = 2 172 dests(1) = ring_rank(world_size, world_rank-1) 173 dests(2) = ring_rank(world_size, world_rank+1) 174 call MPI_Dist_graph_create_adjacent(MPI_COMM_WORLD, & 175 & src_sz, srcs, & 176 & MPI_UNWEIGHTED, & 177 & dest_sz, dests, & 178 & MPI_UNWEIGHTED, & 179 & MPI_INFO_NULL, .true., & 180 & dgraph_comm, ierr) 181 if (ierr .ne. MPI_SUCCESS) then 182 write(6,*) "MPI_Dist_graph_create_adjacent() fails!" 183 call MPI_Abort(MPI_COMM_WORLD, 1, ierr) 184 stop 185 endif 186 if (.not. validate_dgraph(dgraph_comm)) then 187 write(6,*) "MPI_Dist_graph_create_adjacent() does not create" & 188 & //"a bidirectional ring graph!" 189 call MPI_Abort(MPI_COMM_WORLD, 1, ierr) 190 stop 191 endif 192 call MPI_Comm_free(dgraph_comm, ierr) 193 194! now create one with MPI_WEIGHTS_EMPTY 195 src_sz = 0 196 dest_sz = 0 197 call MPI_Dist_graph_create_adjacent(MPI_COMM_WORLD, & 198 & src_sz, srcs, & 199 & MPI_WEIGHTS_EMPTY, & 200 & dest_sz, dests, & 201 & MPI_WEIGHTS_EMPTY, & 202 & MPI_INFO_NULL, .true., & 203 & dgraph_comm, ierr) 204 if (ierr .ne. MPI_SUCCESS) then 205 write(6,*) "MPI_Dist_graph_create_adjacent() fails!" 206 call MPI_Abort(MPI_COMM_WORLD, 1, ierr) 207 stop 208 endif 209 call MPI_Comm_free(dgraph_comm, ierr) 210 211 call MTEST_Finalize(errs) 212 end 213