Changeset 965 for codes/icosagcm/trunk/src/parallel/transfert_requests.F90
- Timestamp:
- 07/25/19 19:53:31 (5 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
codes/icosagcm/trunk/src/parallel/transfert_requests.F90
r963 r965 2 2 module transfert_request_mod 3 3 use abort_mod 4 use field_mod, only : t_field, field_T, field_U 4 use field_mod, only : t_field, field_T, field_U, field_Z 5 5 use domain_mod, only : ndomain, ndomain_glo, domain, domain_glo, domloc_glo_ind, domglo_rank, domglo_loc_ind 6 6 implicit none … … 9 9 ! Points for t_request : list of points coordinates in a domain 10 10 type t_points 11 integer, allocatable :: i(:), j(:), edge(:) 11 integer, allocatable :: i(:), j(:) 12 integer, allocatable :: elt(:) ! Element : cell edge or cell point 12 13 integer :: npoints 13 14 end type … … 26 27 type(t_request),save,pointer :: req_e1_scal(:) ! Halos for U fields 27 28 type(t_request),save,pointer :: req_e1_vect(:) ! Halos for U fields (with sign change) 28 type(t_request),save,pointer :: req_z1_scal(:) => null() ! DO NOT USE : not initialized29 type(t_request),save,pointer :: req_z1_scal(:) ! Halos for Z fields 29 30 30 31 type(t_request),save,pointer :: req_i0(:) ! Duplicated cells for T fields … … 42 43 subroutine init_all_requests 43 44 use dimensions, only : swap_dimensions, ii_begin, ii_end, jj_begin, jj_end 44 use metric, only : right, rdown, ldown, left, lup, rup 45 use metric, only : right, rdown, ldown, left, lup, rup, vdown, vlup, vrdown, vrup, vup 45 46 integer :: ind, i, j 46 47 … … 131 132 req_e1_vect = req_e1_scal 132 133 req_e1_vect%vector = .true. 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) 133 160 contains 134 161 135 162 ! ---- Points ---- 136 163 ! Initialize (allocate) points 137 subroutine points_init(points, has_e dges)164 subroutine points_init(points, has_elt) 138 165 type(t_points):: points 139 logical :: has_e dges166 logical :: has_elt 140 167 integer, parameter :: INITIAL_ALLOC_SIZE = 10 141 168 142 169 allocate(points%i(INITIAL_ALLOC_SIZE)) 143 170 allocate(points%j(INITIAL_ALLOC_SIZE)) 144 if(has_e dges) allocate(points%edge(INITIAL_ALLOC_SIZE))171 if(has_elt) allocate(points%elt(INITIAL_ALLOC_SIZE)) 145 172 points%npoints = 0 146 173 end subroutine 147 174 148 ! Append point (i, j [,e dge]) to point list149 subroutine points_append(points, i, j, e dge)175 ! Append point (i, j [,elt]) to point list 176 subroutine points_append(points, i, j, elt) 150 177 type(t_points), intent(inout):: points 151 178 integer, intent(in) :: i, j 152 integer, intent(in), optional :: e dge179 integer, intent(in), optional :: elt 153 180 integer :: begin_size 154 181 … … 157 184 begin_size=points%npoints 158 185 call array_append( points%j, begin_size, j ) 159 if(present(e dge)) then186 if(present(elt)) then 160 187 begin_size=points%npoints 161 call array_append( points%e dge, begin_size, edge)188 call array_append( points%elt, begin_size, elt ) 162 189 end if 163 190 164 191 points%npoints = points%npoints + 1 165 192 end subroutine 166 193 167 194 ! Add element to array, and reallocate if necessary 168 195 subroutine array_append( a, a_size, elt ) … … 172 199 integer, allocatable :: a_tmp(:) 173 200 integer, parameter :: GROW_FACTOR = 2 174 201 175 202 if( size( a ) <= a_size ) then 176 203 allocate( a_tmp ( a_size * GROW_FACTOR ) ) … … 192 219 193 220 displs = (points%i-1) + (points%j-1)*iim 194 if(allocated(points%e dge)) displs = displs + (points%edge-1)*iim*jjm221 if(allocated(points%elt)) displs = displs + (points%elt-1)*iim*jjm 195 222 196 223 last=0 … … 210 237 points%j(k) = mod((displs(k)-(points%i(k)-1))/iim, jjm)+1 211 238 end do 212 if(allocated(points%e dge)) then213 deallocate(points%e dge); allocate(points%edge(last))214 points%e dge(:) = displs(1:last)/(iim*jjm) + 1239 if(allocated(points%elt)) then 240 deallocate(points%elt); allocate(points%elt(last)) 241 points%elt(:) = displs(1:last)/(iim*jjm) + 1 215 242 end if 216 243 … … 231 258 allocate(request(ind)%points_HtoB(ndomain_glo)) 232 259 do remote_ind_glo = 1, ndomain_glo 233 call points_init(request(ind)%points_BtoH(remote_ind_glo), field_type == field_U)260 call points_init(request(ind)%points_BtoH(remote_ind_glo), field_type /= field_T) 234 261 !call points_init(request(ind)%points_HtoB(remote_ind_glo)) ! Initialized before MPI communication in request_exchange 235 262 end do … … 238 265 239 266 ! Append local point (ind_glo,i,j) to the point list of domain ind, where ind_glo is the domain owning the point 240 subroutine request_add_point(req, ind, i, j, e dge)267 subroutine request_add_point(req, ind, i, j, elt) 241 268 type(t_request), intent(inout) :: req(:) 242 269 integer :: ind, i, j 243 integer, optional :: edge 244 245 if( present(edge) ) then ! U field 246 if(domain(ind)%edge_assign_domain(edge-1,i,j) /= domloc_glo_ind(ind) ) & ! Add only if remote domain != local domain 247 call points_append(req(ind)%points_BtoH(domain(ind)%edge_assign_domain(edge-1,i,j)), i, j, edge ) 270 integer, optional :: elt 271 272 if( req(ind)%field_type == field_U ) then 273 if(domain(ind)%edge_assign_domain(elt-1,i,j) /= domloc_glo_ind(ind) ) & ! Add only if remote domain != local domain 274 call points_append(req(ind)%points_BtoH(domain(ind)%edge_assign_domain(elt-1,i,j)), i, j, elt ) 275 else if( req(ind)%field_type == field_Z ) then 276 if(domain(ind)%vertex_assign_domain(elt-1,i,j) /= domloc_glo_ind(ind) ) & ! Add only if remote domain != local domain 277 call points_append(req(ind)%points_BtoH(domain(ind)%vertex_assign_domain(elt-1,i,j)), i, j, elt ) 248 278 else ! T field 249 279 if(domain(ind)%assign_domain(i,j) /= domloc_glo_ind(ind) ) & … … 262 292 263 293 integer :: ind_loc, remote_ind_glo, k 264 integer :: i, j, e dge294 integer :: i, j, elt 265 295 integer :: reqs(ndomain, ndomain_glo, 2), reqs2(ndomain, ndomain_glo, 2), reqs3(ndomain, ndomain_glo, 3), ierr 266 296 … … 271 301 points => req(ind_loc)%points_BtoH(remote_ind_glo) 272 302 call points_sort(points) 273 call points_init( points_send(ind_loc,remote_ind_glo), req(ind_loc)%field_type == field_U)303 call points_init( points_send(ind_loc,remote_ind_glo), req(ind_loc)%field_type /= field_T ) 274 304 do k = 1, points%npoints 275 305 i = points%i(k) 276 306 j = points%j(k) 277 if(req(ind_loc)%field_type == field_T) then 307 if( req(ind_loc)%field_type == field_U ) then 308 elt = points%elt(k) - 1 ! WARNING : domain%edge_assign_* use index in 0:5 instead of 1:6 everywhere else 309 call points_append(points_send(ind_loc,remote_ind_glo), domain(ind_loc)%edge_assign_i(elt,i,j), & 310 domain(ind_loc)%edge_assign_j(elt,i,j), & 311 domain(ind_loc)%edge_assign_pos(elt,i,j)+1) 312 else if( req(ind_loc)%field_type == field_Z ) then 313 elt = points%elt(k) - 1 ! WARNING : domain%vertex_assign_* use index in 0:5 instead of 1:6 everywhere else 314 call points_append(points_send(ind_loc,remote_ind_glo), domain(ind_loc)%vertex_assign_i(elt,i,j), & 315 domain(ind_loc)%vertex_assign_j(elt,i,j), & 316 domain(ind_loc)%vertex_assign_pos(elt,i,j)+1) 317 else ! T field 278 318 call points_append(points_send(ind_loc,remote_ind_glo), domain(ind_loc)%assign_i(i,j), & 279 319 domain(ind_loc)%assign_j(i,j)) 280 else 281 edge = points%edge(k) - 1 ! WARNING : domain%edge_assign_* use index in 0:5 instead of 1:6 everywhere else 282 call points_append(points_send(ind_loc,remote_ind_glo), domain(ind_loc)%edge_assign_i(edge,i,j), & 283 domain(ind_loc)%edge_assign_j(edge,i,j), & 284 domain(ind_loc)%edge_assign_pos(edge,i,j)+1) 285 end if 320 endif 286 321 end do 287 322 end do … … 317 352 domloc_glo_ind(ind_loc)+ndomain_glo*remote_ind_glo, comm_icosa, reqs2(ind_loc,remote_ind_glo,2), ierr) 318 353 319 if(req(1)%field_type == field_U) then320 call MPI_Isend( points_send(ind_loc,remote_ind_glo)%e dge, points_send(ind_loc,remote_ind_glo)%npoints, MPI_INT, &354 if(req(1)%field_type /= field_T) then 355 call MPI_Isend( points_send(ind_loc,remote_ind_glo)%elt, points_send(ind_loc,remote_ind_glo)%npoints, MPI_INT, & 321 356 domglo_rank(remote_ind_glo), remote_ind_glo+ndomain_glo*domloc_glo_ind(ind_loc), & 322 357 comm_icosa, reqs3(ind_loc,remote_ind_glo,1), ierr ) 323 allocate( points%e dge( points%npoints ) )324 call MPI_Irecv( points%e dge, points%npoints, MPI_INT, &358 allocate( points%elt( points%npoints ) ) 359 call MPI_Irecv( points%elt, points%npoints, MPI_INT, & 325 360 domglo_rank(remote_ind_glo), domloc_glo_ind(ind_loc)+ndomain_glo*remote_ind_glo, & 326 361 comm_icosa, reqs3(ind_loc,remote_ind_glo,2), ierr ) … … 330 365 call waitall(reqs,ndomain*ndomain_glo*2) 331 366 call waitall(reqs2,ndomain*ndomain_glo*2) 332 if(req(1)%field_type == field_U) call waitall(reqs3,ndomain*ndomain_glo*2)367 if(req(1)%field_type /= field_T) call waitall(reqs3,ndomain*ndomain_glo*2) 333 368 end subroutine 334 369
Note: See TracChangeset
for help on using the changeset viewer.