[963] | 1 | ! This module contains the requests for transfert_mpi |
---|
| 2 | module transfert_request_mod |
---|
| 3 | use abort_mod |
---|
[965] | 4 | use field_mod, only : t_field, field_T, field_U, field_Z |
---|
[963] | 5 | use domain_mod, only : ndomain, ndomain_glo, domain, domain_glo, domloc_glo_ind, domglo_rank, domglo_loc_ind |
---|
| 6 | implicit none |
---|
| 7 | private |
---|
| 8 | |
---|
| 9 | ! Points for t_request : list of points coordinates in a domain |
---|
| 10 | type t_points |
---|
[965] | 11 | integer, allocatable :: i(:), j(:) |
---|
| 12 | integer, allocatable :: elt(:) ! Element : cell edge or cell point |
---|
[963] | 13 | integer :: npoints |
---|
| 14 | end type |
---|
| 15 | |
---|
| 16 | ! Describes how to exchange data from/to one domain to/from other domains |
---|
| 17 | type t_request |
---|
| 18 | integer :: field_type |
---|
| 19 | type(t_points), allocatable :: points_BtoH(:) ! req(local_ind)%points_BtoH(remote_ind) are the points of |
---|
| 20 | ! domain local_ind to unpack from mpi buffer recieved from remote_ind |
---|
| 21 | type(t_points), allocatable :: points_HtoB(:) ! req(local_ind)%points_HtoB(remote_ind) are the points of |
---|
| 22 | ! domain local_ind to pack to mpi buffer to send to remote_ind |
---|
| 23 | logical :: vector = .false. ! Is a vector request (true if sign change needed) |
---|
| 24 | end type t_request |
---|
| 25 | |
---|
| 26 | type(t_request),save,pointer :: req_i1(:) ! Halos for T fields |
---|
| 27 | type(t_request),save,pointer :: req_e1_scal(:) ! Halos for U fields |
---|
| 28 | type(t_request),save,pointer :: req_e1_vect(:) ! Halos for U fields (with sign change) |
---|
[965] | 29 | type(t_request),save,pointer :: req_z1_scal(:) ! Halos for Z fields |
---|
[963] | 30 | |
---|
| 31 | type(t_request),save,pointer :: req_i0(:) ! Duplicated cells for T fields |
---|
| 32 | type(t_request),save,pointer :: req_e0_scal(:) |
---|
| 33 | type(t_request),save,pointer :: req_e0_vect(:) |
---|
| 34 | |
---|
| 35 | public :: t_points, t_request, & |
---|
| 36 | req_i1, req_e1_scal, req_e1_vect, & |
---|
| 37 | req_i0, req_e0_scal, req_e0_vect, & |
---|
| 38 | req_z1_scal, & |
---|
| 39 | init_all_requests |
---|
| 40 | contains |
---|
| 41 | |
---|
| 42 | ! Initialize requests |
---|
| 43 | subroutine init_all_requests |
---|
| 44 | use dimensions, only : swap_dimensions, ii_begin, ii_end, jj_begin, jj_end |
---|
[965] | 45 | use metric, only : right, rdown, ldown, left, lup, rup, vdown, vlup, vrdown, vrup, vup |
---|
[963] | 46 | integer :: ind, i, j |
---|
| 47 | |
---|
| 48 | ! Create requests req_* see corresponding module variables for description |
---|
| 49 | call request_init(req_i0, field_T) |
---|
| 50 | DO ind=1,ndomain |
---|
| 51 | CALL swap_dimensions(ind) |
---|
| 52 | DO i=ii_begin,ii_end |
---|
| 53 | CALL request_add_point(req_i0,ind,i,jj_begin) |
---|
| 54 | ENDDO |
---|
| 55 | DO j=jj_begin,jj_end |
---|
| 56 | CALL request_add_point(req_i0,ind,ii_begin,j) |
---|
| 57 | CALL request_add_point(req_i0,ind,ii_end,j) |
---|
| 58 | ENDDO |
---|
| 59 | DO i=ii_begin,ii_end |
---|
| 60 | CALL request_add_point(req_i0,ind,i,jj_end) |
---|
| 61 | ENDDO |
---|
| 62 | ENDDO |
---|
| 63 | call request_exchange(req_i0) |
---|
| 64 | |
---|
| 65 | call request_init(req_i1, field_T) |
---|
| 66 | do ind=1,ndomain |
---|
| 67 | call swap_dimensions(ind) |
---|
| 68 | do i=ii_begin,ii_end+1 |
---|
| 69 | call request_add_point(req_i1,ind,i,jj_begin-1) |
---|
| 70 | enddo |
---|
| 71 | do j=jj_begin,jj_end |
---|
| 72 | call request_add_point(req_i1,ind,ii_begin-1,j) |
---|
| 73 | call request_add_point(req_i1,ind,ii_end+1,j) |
---|
| 74 | enddo |
---|
| 75 | call request_add_point(req_i1,ind,ii_begin-1,jj_end+1) |
---|
| 76 | do i=ii_begin,ii_end |
---|
| 77 | call request_add_point(req_i1,ind,i,jj_end+1) |
---|
| 78 | enddo |
---|
| 79 | enddo |
---|
| 80 | call request_exchange(req_i1) |
---|
| 81 | |
---|
| 82 | call request_init(req_e0_scal, field_U) |
---|
| 83 | DO ind=1,ndomain |
---|
| 84 | CALL swap_dimensions(ind) |
---|
| 85 | DO i=ii_begin+1,ii_end-1 |
---|
| 86 | CALL request_add_point(req_e0_scal,ind,i,jj_begin,right) |
---|
| 87 | CALL request_add_point(req_e0_scal,ind,i,jj_end,right) |
---|
| 88 | ENDDO |
---|
| 89 | DO j=jj_begin+1,jj_end-1 |
---|
| 90 | CALL request_add_point(req_e0_scal,ind,ii_begin,j,rup) |
---|
| 91 | CALL request_add_point(req_e0_scal,ind,ii_end,j,rup) |
---|
| 92 | ENDDO |
---|
| 93 | CALL request_add_point(req_e0_scal,ind,ii_begin+1,jj_begin,left) |
---|
| 94 | CALL request_add_point(req_e0_scal,ind,ii_begin,jj_begin+1,ldown) |
---|
| 95 | CALL request_add_point(req_e0_scal,ind,ii_begin+1,jj_end,left) |
---|
| 96 | CALL request_add_point(req_e0_scal,ind,ii_end,jj_begin+1,ldown) |
---|
| 97 | ENDDO |
---|
| 98 | call request_exchange(req_e0_scal) |
---|
| 99 | |
---|
| 100 | ! req_e0_vect is the same request than req_e0_scal but with the vector flag to true |
---|
| 101 | allocate( req_e0_vect(ndomain) ) |
---|
| 102 | req_e0_vect = req_e0_scal |
---|
| 103 | req_e0_vect%vector = .true. |
---|
| 104 | |
---|
| 105 | call request_init(req_e1_scal, field_U) |
---|
| 106 | DO ind=1,ndomain |
---|
| 107 | CALL swap_dimensions(ind) |
---|
| 108 | DO i=ii_begin,ii_end |
---|
| 109 | CALL request_add_point(req_e1_scal,ind,i,jj_begin-1,rup) |
---|
| 110 | CALL request_add_point(req_e1_scal,ind,i+1,jj_begin-1,lup) |
---|
| 111 | ENDDO |
---|
| 112 | DO j=jj_begin,jj_end |
---|
| 113 | CALL request_add_point(req_e1_scal,ind,ii_end+1,j,left) |
---|
| 114 | ENDDO |
---|
| 115 | DO j=jj_begin,jj_end |
---|
| 116 | CALL request_add_point(req_e1_scal,ind,ii_end+1,j-1,lup) |
---|
| 117 | ENDDO |
---|
| 118 | DO i=ii_begin,ii_end |
---|
| 119 | CALL request_add_point(req_e1_scal,ind,i,jj_end+1,ldown) |
---|
| 120 | CALL request_add_point(req_e1_scal,ind,i-1,jj_end+1,rdown) |
---|
| 121 | ENDDO |
---|
| 122 | DO j=jj_begin,jj_end |
---|
| 123 | CALL request_add_point(req_e1_scal,ind,ii_begin-1,j,right) |
---|
| 124 | ENDDO |
---|
| 125 | DO j=jj_begin,jj_end |
---|
| 126 | CALL request_add_point(req_e1_scal,ind,ii_begin-1,j+1,rdown) |
---|
| 127 | ENDDO |
---|
| 128 | ENDDO |
---|
| 129 | call request_exchange(req_e1_scal) |
---|
| 130 | |
---|
| 131 | allocate( req_e1_vect(ndomain) ) |
---|
| 132 | req_e1_vect = req_e1_scal |
---|
| 133 | req_e1_vect%vector = .true. |
---|
[965] | 134 | |
---|
| 135 | call request_init(req_z1_scal, field_z) |
---|
| 136 | do ind=1,ndomain |
---|
| 137 | call swap_dimensions(ind) |
---|
| 138 | do i=ii_begin,ii_end |
---|
| 139 | call request_add_point(req_z1_scal,ind,i,jj_begin-1,vrup) |
---|
| 140 | call request_add_point(req_z1_scal,ind,i+1,jj_begin-1,vup) |
---|
| 141 | enddo |
---|
| 142 | do j=jj_begin,jj_end |
---|
| 143 | call request_add_point(req_z1_scal,ind,ii_end+1,j,vlup) |
---|
| 144 | enddo |
---|
| 145 | do j=jj_begin,jj_end |
---|
| 146 | call request_add_point(req_z1_scal,ind,ii_end+1,j-1,vup) |
---|
| 147 | enddo |
---|
| 148 | do i=ii_begin,ii_end |
---|
| 149 | call request_add_point(req_z1_scal,ind,i,jj_end+1,vdown) |
---|
| 150 | call request_add_point(req_z1_scal,ind,i-1,jj_end+1,vrdown) |
---|
| 151 | enddo |
---|
| 152 | do j=jj_begin,jj_end |
---|
| 153 | call request_add_point(req_z1_scal,ind,ii_begin-1,j,vrup) |
---|
| 154 | enddo |
---|
| 155 | do j=jj_begin,jj_end |
---|
| 156 | call request_add_point(req_z1_scal,ind,ii_begin-1,j,vrdown) |
---|
| 157 | enddo |
---|
| 158 | enddo |
---|
| 159 | call request_exchange(req_z1_scal) |
---|
[963] | 160 | contains |
---|
| 161 | |
---|
| 162 | ! ---- Points ---- |
---|
| 163 | ! Initialize (allocate) points |
---|
[965] | 164 | subroutine points_init(points, has_elt) |
---|
[963] | 165 | type(t_points):: points |
---|
[965] | 166 | logical :: has_elt |
---|
[963] | 167 | integer, parameter :: INITIAL_ALLOC_SIZE = 10 |
---|
| 168 | |
---|
| 169 | allocate(points%i(INITIAL_ALLOC_SIZE)) |
---|
| 170 | allocate(points%j(INITIAL_ALLOC_SIZE)) |
---|
[965] | 171 | if(has_elt) allocate(points%elt(INITIAL_ALLOC_SIZE)) |
---|
[963] | 172 | points%npoints = 0 |
---|
| 173 | end subroutine |
---|
| 174 | |
---|
[965] | 175 | ! Append point (i, j [,elt]) to point list |
---|
| 176 | subroutine points_append(points, i, j, elt) |
---|
[963] | 177 | type(t_points), intent(inout):: points |
---|
| 178 | integer, intent(in) :: i, j |
---|
[965] | 179 | integer, intent(in), optional :: elt |
---|
[963] | 180 | integer :: begin_size |
---|
| 181 | |
---|
| 182 | begin_size=points%npoints |
---|
| 183 | call array_append( points%i, begin_size, i ) |
---|
| 184 | begin_size=points%npoints |
---|
| 185 | call array_append( points%j, begin_size, j ) |
---|
[965] | 186 | if(present(elt)) then |
---|
[963] | 187 | begin_size=points%npoints |
---|
[965] | 188 | call array_append( points%elt, begin_size, elt ) |
---|
[963] | 189 | end if |
---|
| 190 | |
---|
| 191 | points%npoints = points%npoints + 1 |
---|
| 192 | end subroutine |
---|
[965] | 193 | |
---|
[963] | 194 | ! Add element to array, and reallocate if necessary |
---|
| 195 | subroutine array_append( a, a_size, elt ) |
---|
| 196 | integer, allocatable, intent(inout) :: a(:) |
---|
| 197 | integer, intent(inout) :: a_size |
---|
| 198 | integer, intent(in) :: elt |
---|
| 199 | integer, allocatable :: a_tmp(:) |
---|
| 200 | integer, parameter :: GROW_FACTOR = 2 |
---|
[965] | 201 | |
---|
[963] | 202 | if( size( a ) <= a_size ) then |
---|
| 203 | allocate( a_tmp ( a_size * GROW_FACTOR ) ) |
---|
| 204 | a_tmp(1:a_size) = a(1:a_size) |
---|
| 205 | call move_alloc(a_tmp, a) |
---|
| 206 | end if |
---|
| 207 | a_size = a_size + 1 |
---|
| 208 | a(a_size) = elt; |
---|
| 209 | end subroutine |
---|
| 210 | |
---|
| 211 | ! Sort and remove duplicates in points |
---|
| 212 | ! dimensions%iim must be set to the right value with swap_dimensions before |
---|
| 213 | ! After this call, size(points%i) == points%npoints |
---|
| 214 | subroutine points_sort( points ) |
---|
| 215 | use dimensions, only : iim, jjm |
---|
| 216 | type(t_points):: points |
---|
| 217 | integer :: displs(points%npoints) |
---|
| 218 | integer :: k, i, last, i_min, tmp |
---|
| 219 | |
---|
[999] | 220 | displs = (points%i(1:points%npoints)-1) + (points%j(1:points%npoints)-1)*iim |
---|
| 221 | if(allocated(points%elt)) displs = displs + (points%elt(1:points%npoints)-1)*iim*jjm |
---|
[963] | 222 | |
---|
| 223 | last=0 |
---|
| 224 | do i=1, points%npoints |
---|
| 225 | i_min = minloc(displs(i:),1)+i-1 |
---|
| 226 | tmp = displs(i_min) |
---|
| 227 | displs(i_min) = displs(i) |
---|
[999] | 228 | if( last==0 ) then |
---|
| 229 | last=last+1 |
---|
| 230 | else if( displs(i_min) /= displs(last) ) then |
---|
| 231 | last=last+1 |
---|
| 232 | endif |
---|
[963] | 233 | displs(last) = tmp |
---|
| 234 | end do |
---|
| 235 | |
---|
| 236 | points%npoints = last |
---|
| 237 | deallocate(points%i); allocate(points%i(last)) |
---|
| 238 | deallocate(points%j); allocate(points%j(last)) |
---|
| 239 | do k=1, last |
---|
| 240 | points%i(k) = mod(displs(k), iim)+1 |
---|
| 241 | points%j(k) = mod((displs(k)-(points%i(k)-1))/iim, jjm)+1 |
---|
| 242 | end do |
---|
[965] | 243 | if(allocated(points%elt)) then |
---|
| 244 | deallocate(points%elt); allocate(points%elt(last)) |
---|
| 245 | points%elt(:) = displs(1:last)/(iim*jjm) + 1 |
---|
[963] | 246 | end if |
---|
| 247 | |
---|
| 248 | if( size(points%i) /= points%npoints ) call dynamico_abort( "Storage too big for points" ) |
---|
| 249 | end subroutine |
---|
| 250 | |
---|
| 251 | ! ---- Requests ---- |
---|
| 252 | ! Initialize request |
---|
| 253 | subroutine request_init( request, field_type ) |
---|
| 254 | type(t_request), pointer, intent(out) :: request(:) |
---|
| 255 | integer, intent(in) :: field_type |
---|
| 256 | integer :: ind, remote_ind_glo |
---|
| 257 | |
---|
| 258 | allocate( request(ndomain) ) |
---|
| 259 | do ind=1, ndomain |
---|
| 260 | request(ind)%field_type = field_type |
---|
| 261 | allocate(request(ind)%points_BtoH(ndomain_glo)) |
---|
| 262 | allocate(request(ind)%points_HtoB(ndomain_glo)) |
---|
| 263 | do remote_ind_glo = 1, ndomain_glo |
---|
[965] | 264 | call points_init(request(ind)%points_BtoH(remote_ind_glo), field_type /= field_T) |
---|
[963] | 265 | !call points_init(request(ind)%points_HtoB(remote_ind_glo)) ! Initialized before MPI communication in request_exchange |
---|
| 266 | end do |
---|
| 267 | end do |
---|
| 268 | end subroutine |
---|
| 269 | |
---|
| 270 | ! Append local point (ind_glo,i,j) to the point list of domain ind, where ind_glo is the domain owning the point |
---|
[965] | 271 | subroutine request_add_point(req, ind, i, j, elt) |
---|
[963] | 272 | type(t_request), intent(inout) :: req(:) |
---|
| 273 | integer :: ind, i, j |
---|
[965] | 274 | integer, optional :: elt |
---|
[963] | 275 | |
---|
[965] | 276 | if( req(ind)%field_type == field_U ) then |
---|
| 277 | if(domain(ind)%edge_assign_domain(elt-1,i,j) /= domloc_glo_ind(ind) ) & ! Add only if remote domain != local domain |
---|
| 278 | call points_append(req(ind)%points_BtoH(domain(ind)%edge_assign_domain(elt-1,i,j)), i, j, elt ) |
---|
| 279 | else if( req(ind)%field_type == field_Z ) then |
---|
| 280 | if(domain(ind)%vertex_assign_domain(elt-1,i,j) /= domloc_glo_ind(ind) ) & ! Add only if remote domain != local domain |
---|
| 281 | call points_append(req(ind)%points_BtoH(domain(ind)%vertex_assign_domain(elt-1,i,j)), i, j, elt ) |
---|
[963] | 282 | else ! T field |
---|
| 283 | if(domain(ind)%assign_domain(i,j) /= domloc_glo_ind(ind) ) & |
---|
| 284 | call points_append(req(ind)%points_BtoH(domain(ind)%assign_domain(i,j)), i, j ) |
---|
| 285 | endif |
---|
| 286 | end subroutine |
---|
| 287 | |
---|
| 288 | ! Each domain send the cell index it needs from other domains |
---|
| 289 | subroutine request_exchange(req) |
---|
| 290 | use mpi_mod |
---|
| 291 | use dimensions, only : swap_dimensions |
---|
| 292 | use mpipara, only : comm_icosa |
---|
| 293 | type(t_request), target, intent(in) :: req(:) |
---|
| 294 | type(t_points) :: points_send(ndomain,ndomain_glo) |
---|
| 295 | type(t_points), pointer :: points |
---|
| 296 | |
---|
| 297 | integer :: ind_loc, remote_ind_glo, k |
---|
[965] | 298 | integer :: i, j, elt |
---|
[963] | 299 | integer :: reqs(ndomain, ndomain_glo, 2), reqs2(ndomain, ndomain_glo, 2), reqs3(ndomain, ndomain_glo, 3), ierr |
---|
| 300 | |
---|
| 301 | ! Transform position in local domain to position in remote domain pos : (i,j) -> pos_send : (assign_i(i,j), assign_j(i,j)) |
---|
| 302 | do ind_loc = 1, ndomain |
---|
| 303 | call swap_dimensions(ind_loc) |
---|
| 304 | do remote_ind_glo = 1, ndomain_glo |
---|
| 305 | points => req(ind_loc)%points_BtoH(remote_ind_glo) |
---|
| 306 | call points_sort(points) |
---|
[965] | 307 | call points_init( points_send(ind_loc,remote_ind_glo), req(ind_loc)%field_type /= field_T ) |
---|
[963] | 308 | do k = 1, points%npoints |
---|
| 309 | i = points%i(k) |
---|
| 310 | j = points%j(k) |
---|
[965] | 311 | if( req(ind_loc)%field_type == field_U ) then |
---|
| 312 | elt = points%elt(k) - 1 ! WARNING : domain%edge_assign_* use index in 0:5 instead of 1:6 everywhere else |
---|
| 313 | call points_append(points_send(ind_loc,remote_ind_glo), domain(ind_loc)%edge_assign_i(elt,i,j), & |
---|
| 314 | domain(ind_loc)%edge_assign_j(elt,i,j), & |
---|
| 315 | domain(ind_loc)%edge_assign_pos(elt,i,j)+1) |
---|
| 316 | else if( req(ind_loc)%field_type == field_Z ) then |
---|
| 317 | elt = points%elt(k) - 1 ! WARNING : domain%vertex_assign_* use index in 0:5 instead of 1:6 everywhere else |
---|
| 318 | call points_append(points_send(ind_loc,remote_ind_glo), domain(ind_loc)%vertex_assign_i(elt,i,j), & |
---|
| 319 | domain(ind_loc)%vertex_assign_j(elt,i,j), & |
---|
| 320 | domain(ind_loc)%vertex_assign_pos(elt,i,j)+1) |
---|
| 321 | else ! T field |
---|
[963] | 322 | call points_append(points_send(ind_loc,remote_ind_glo), domain(ind_loc)%assign_i(i,j), & |
---|
| 323 | domain(ind_loc)%assign_j(i,j)) |
---|
[965] | 324 | endif |
---|
[963] | 325 | end do |
---|
| 326 | end do |
---|
| 327 | end do |
---|
| 328 | |
---|
| 329 | ! Exchange sizes |
---|
| 330 | do ind_loc = 1, ndomain |
---|
| 331 | do remote_ind_glo = 1, ndomain_glo |
---|
| 332 | call MPI_Isend( points_send(ind_loc,remote_ind_glo)%npoints, 1, MPI_INT, domglo_rank(remote_ind_glo), & |
---|
| 333 | remote_ind_glo+ndomain_glo*domloc_glo_ind(ind_loc), comm_icosa, reqs(ind_loc,remote_ind_glo,1), ierr ) |
---|
| 334 | call MPI_Irecv( req(ind_loc)%points_HtoB(remote_ind_glo)%npoints, 1, MPI_INT, domglo_rank(remote_ind_glo), & |
---|
| 335 | domloc_glo_ind(ind_loc)+ndomain_glo*remote_ind_glo, comm_icosa, reqs(ind_loc,remote_ind_glo,2), ierr ) |
---|
| 336 | end do |
---|
| 337 | end do |
---|
| 338 | call waitall(reqs,ndomain*ndomain_glo*2) |
---|
| 339 | |
---|
| 340 | ! Exchange points_send -> points_recv |
---|
| 341 | do ind_loc = 1, ndomain |
---|
| 342 | do remote_ind_glo = 1, ndomain_glo |
---|
| 343 | call MPI_Isend( points_send(ind_loc,remote_ind_glo)%i, points_send(ind_loc,remote_ind_glo)%npoints, MPI_INT, & |
---|
| 344 | domglo_rank(remote_ind_glo), remote_ind_glo+ndomain_glo*domloc_glo_ind(ind_loc), & |
---|
| 345 | comm_icosa, reqs(ind_loc,remote_ind_glo,1), ierr ) |
---|
| 346 | call MPI_Isend( points_send(ind_loc,remote_ind_glo)%j, points_send(ind_loc,remote_ind_glo)%npoints, MPI_INT, & |
---|
| 347 | domglo_rank(remote_ind_glo), remote_ind_glo+ndomain_glo*domloc_glo_ind(ind_loc), & |
---|
| 348 | comm_icosa, reqs2(ind_loc,remote_ind_glo,1), ierr ) |
---|
| 349 | |
---|
| 350 | points => req(ind_loc)%points_HtoB(remote_ind_glo) ! Points to recieve (HtoB) |
---|
| 351 | allocate( points%i( points%npoints ) ) |
---|
| 352 | allocate( points%j( points%npoints ) ) |
---|
| 353 | call MPI_Irecv( points%i, points%npoints, MPI_INT, domglo_rank(remote_ind_glo), & |
---|
| 354 | domloc_glo_ind(ind_loc)+ndomain_glo*remote_ind_glo, comm_icosa, reqs(ind_loc,remote_ind_glo,2), ierr) |
---|
| 355 | call MPI_Irecv( points%j, points%npoints, MPI_INT, domglo_rank(remote_ind_glo), & |
---|
| 356 | domloc_glo_ind(ind_loc)+ndomain_glo*remote_ind_glo, comm_icosa, reqs2(ind_loc,remote_ind_glo,2), ierr) |
---|
| 357 | |
---|
[965] | 358 | if(req(1)%field_type /= field_T) then |
---|
| 359 | call MPI_Isend( points_send(ind_loc,remote_ind_glo)%elt, points_send(ind_loc,remote_ind_glo)%npoints, MPI_INT, & |
---|
[963] | 360 | domglo_rank(remote_ind_glo), remote_ind_glo+ndomain_glo*domloc_glo_ind(ind_loc), & |
---|
| 361 | comm_icosa, reqs3(ind_loc,remote_ind_glo,1), ierr ) |
---|
[965] | 362 | allocate( points%elt( points%npoints ) ) |
---|
| 363 | call MPI_Irecv( points%elt, points%npoints, MPI_INT, & |
---|
[963] | 364 | domglo_rank(remote_ind_glo), domloc_glo_ind(ind_loc)+ndomain_glo*remote_ind_glo, & |
---|
| 365 | comm_icosa, reqs3(ind_loc,remote_ind_glo,2), ierr ) |
---|
| 366 | endif |
---|
| 367 | end do |
---|
| 368 | end do |
---|
| 369 | call waitall(reqs,ndomain*ndomain_glo*2) |
---|
| 370 | call waitall(reqs2,ndomain*ndomain_glo*2) |
---|
[965] | 371 | if(req(1)%field_type /= field_T) call waitall(reqs3,ndomain*ndomain_glo*2) |
---|
[963] | 372 | end subroutine |
---|
| 373 | |
---|
| 374 | ! MPI_Waitall, but with a N-Dim array of requests |
---|
| 375 | subroutine waitall(reqs, size) |
---|
| 376 | use mpi_mod |
---|
| 377 | integer :: size, reqs(size), ierr |
---|
| 378 | call MPI_Waitall( size, reqs, MPI_STATUSES_IGNORE, ierr) |
---|
| 379 | end subroutine |
---|
| 380 | end subroutine |
---|
[999] | 381 | end module |
---|