Ignore:
Timestamp:
01/10/20 12:20:51 (4 years ago)
Author:
adurocher
Message:

transfert_mpi : Added 'npoints' in submessages

message_*%npoints is used instead of size(message_*) : PGI wasn't able to understand how many gangs should be used

File:
1 edited

Legend:

Unmodified
Added
Removed
  • codes/icosagcm/trunk/src/parallel/transfert_mpi.f90

    r965 r998  
    1313  type t_local_submessage 
    1414    integer :: src_ind_loc, dest_ind_loc ! index of local and remote domain 
     15    integer :: npoints ! Number of cells to transfer (dim12) 
    1516    integer, allocatable :: displ_src(:) ! List of indexes to copy from domain src_ind_loc 
    1617    integer, allocatable :: displ_dest(:) ! List of indexes to copy to domain dest_ind_loc 
     
    2122  type t_submessage 
    2223    integer :: ind_loc, remote_ind_glo ! index of local and remote domain 
     24    integer :: npoints ! Number of cells to transfer (dim12) 
    2325    integer, allocatable :: displs(:) ! List of indexes to copy from field to buffer for each level 
    2426    integer, allocatable :: sign(:) ! Sign change to be applied for vector requests 
     
    145147            submessage_local%src_ind_loc = ind_loc 
    146148            submessage_local%dest_ind_loc = domglo_loc_ind(remote_ind_glo) 
     149            submessage_local%npoints = submessage_out%npoints 
    147150            submessage_local%displ_src = submessage_out%displs 
    148151            submessage_local%displ_dest = submessage_in%displs 
     
    202205      submessage%ind_loc = ind_loc 
    203206      submessage%remote_ind_glo = remote_ind_glo 
     207      submessage%npoints = points%npoints 
    204208      allocate( submessage%buff( points%npoints, dim3, dim4 ) ) 
    205209      allocate( submessage%displs( points%npoints ) ) 
     
    406410          do d3 = d3_begin, d3_end 
    407411            !$acc loop 
    408              do k = 1, size(message%message_out(i)%displs) 
     412             do k = 1, message%message_out(i)%npoints 
    409413              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) 
    410414            end do 
     
    423427            ! Cannot collapse because size(displ_dest) varies with i 
    424428            !$acc loop vector 
    425             do k = 1, size(message%message_local(i)%displ_dest) 
     429            do k = 1, message%message_local(i)%npoints 
    426430              message%field(message%message_local(i)%dest_ind_loc)%rval4d(message%message_local(i)%displ_dest(k),d3,d4) = & 
    427431                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) 
     
    504508          do d3 = d3_begin, d3_end 
    505509            !$acc loop 
    506             do k = 1, size(message%message_in(i)%displs) 
     510            do k = 1, message%message_in(i)%npoints 
    507511              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) 
    508512            end do 
Note: See TracChangeset for help on using the changeset viewer.