- Timestamp:
- 01/10/20 12:20:51 (4 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
codes/icosagcm/trunk/src/parallel/transfert_mpi.f90
r965 r998 13 13 type t_local_submessage 14 14 integer :: src_ind_loc, dest_ind_loc ! index of local and remote domain 15 integer :: npoints ! Number of cells to transfer (dim12) 15 16 integer, allocatable :: displ_src(:) ! List of indexes to copy from domain src_ind_loc 16 17 integer, allocatable :: displ_dest(:) ! List of indexes to copy to domain dest_ind_loc … … 21 22 type t_submessage 22 23 integer :: ind_loc, remote_ind_glo ! index of local and remote domain 24 integer :: npoints ! Number of cells to transfer (dim12) 23 25 integer, allocatable :: displs(:) ! List of indexes to copy from field to buffer for each level 24 26 integer, allocatable :: sign(:) ! Sign change to be applied for vector requests … … 145 147 submessage_local%src_ind_loc = ind_loc 146 148 submessage_local%dest_ind_loc = domglo_loc_ind(remote_ind_glo) 149 submessage_local%npoints = submessage_out%npoints 147 150 submessage_local%displ_src = submessage_out%displs 148 151 submessage_local%displ_dest = submessage_in%displs … … 202 205 submessage%ind_loc = ind_loc 203 206 submessage%remote_ind_glo = remote_ind_glo 207 submessage%npoints = points%npoints 204 208 allocate( submessage%buff( points%npoints, dim3, dim4 ) ) 205 209 allocate( submessage%displs( points%npoints ) ) … … 406 410 do d3 = d3_begin, d3_end 407 411 !$acc loop 408 do k = 1, size(message%message_out(i)%displs)412 do k = 1, message%message_out(i)%npoints 409 413 message%message_out(i)%buff(k,d3,d4) = message%field(message%message_out(i)%ind_loc)%rval4d(message%message_out(i)%displs(k),d3, d4) 410 414 end do … … 423 427 ! Cannot collapse because size(displ_dest) varies with i 424 428 !$acc loop vector 425 do k = 1, size(message%message_local(i)%displ_dest)429 do k = 1, message%message_local(i)%npoints 426 430 message%field(message%message_local(i)%dest_ind_loc)%rval4d(message%message_local(i)%displ_dest(k),d3,d4) = & 427 431 message%message_local(i)%sign(k)*message%field(message%message_local(i)%src_ind_loc)%rval4d(message%message_local(i)%displ_src(k),d3,d4) … … 504 508 do d3 = d3_begin, d3_end 505 509 !$acc loop 506 do k = 1, size(message%message_in(i)%displs)510 do k = 1, message%message_in(i)%npoints 507 511 message%field(message%message_in(i)%ind_loc)%rval4d(message%message_in(i)%displs(k),d3,d4) = message%message_in(i)%sign(k)*message%message_in(i)%buff(k,d3,d4) 508 512 end do
Note: See TracChangeset
for help on using the changeset viewer.