[26] | 1 | MODULE transfert_mpi_mod |
---|
| 2 | USE genmod |
---|
| 3 | |
---|
| 4 | TYPE array |
---|
| 5 | INTEGER,POINTER :: value(:) |
---|
[146] | 6 | INTEGER,POINTER :: sign(:) |
---|
[26] | 7 | INTEGER :: domain |
---|
| 8 | INTEGER :: rank |
---|
| 9 | INTEGER :: size |
---|
| 10 | INTEGER,POINTER :: buffer(:) |
---|
| 11 | REAL,POINTER :: buffer_r2(:) |
---|
| 12 | REAL,POINTER :: buffer_r3(:,:) |
---|
| 13 | REAL,POINTER :: buffer_r4(:,:,:) |
---|
| 14 | END TYPE array |
---|
| 15 | |
---|
| 16 | TYPE t_request |
---|
| 17 | INTEGER :: type_field |
---|
| 18 | INTEGER :: max_size |
---|
| 19 | INTEGER :: size |
---|
[146] | 20 | LOGICAL :: vector |
---|
[26] | 21 | INTEGER,POINTER :: src_domain(:) |
---|
| 22 | INTEGER,POINTER :: src_i(:) |
---|
| 23 | INTEGER,POINTER :: src_j(:) |
---|
| 24 | INTEGER,POINTER :: src_ind(:) |
---|
| 25 | INTEGER,POINTER :: target_ind(:) |
---|
| 26 | INTEGER,POINTER :: target_i(:) |
---|
| 27 | INTEGER,POINTER :: target_j(:) |
---|
[146] | 28 | INTEGER,POINTER :: target_sign(:) |
---|
[26] | 29 | INTEGER :: nrecv |
---|
| 30 | TYPE(ARRAY),POINTER :: recv(:) |
---|
| 31 | INTEGER :: nsend |
---|
| 32 | TYPE(ARRAY),POINTER :: send(:) |
---|
| 33 | END TYPE t_request |
---|
| 34 | |
---|
| 35 | TYPE(t_request),POINTER :: req_i1(:) |
---|
[146] | 36 | TYPE(t_request),POINTER :: req_e1_scal(:) |
---|
| 37 | TYPE(t_request),POINTER :: req_e1_vect(:) |
---|
[26] | 38 | |
---|
| 39 | |
---|
| 40 | CONTAINS |
---|
| 41 | |
---|
| 42 | SUBROUTINE init_transfert |
---|
| 43 | USE domain_mod |
---|
| 44 | USE dimensions |
---|
| 45 | USE field_mod |
---|
| 46 | USE metric |
---|
| 47 | USE mpipara |
---|
| 48 | IMPLICIT NONE |
---|
| 49 | INTEGER :: ind,i,j |
---|
| 50 | |
---|
| 51 | CALL create_request(field_t,req_i1) |
---|
| 52 | |
---|
| 53 | DO ind=1,ndomain |
---|
| 54 | CALL swap_dimensions(ind) |
---|
| 55 | DO i=ii_begin,ii_end+1 |
---|
| 56 | CALL request_add_point(ind,i,jj_begin-1,req_i1) |
---|
| 57 | ENDDO |
---|
| 58 | |
---|
| 59 | DO j=jj_begin,jj_end |
---|
| 60 | CALL request_add_point(ind,ii_end+1,j,req_i1) |
---|
| 61 | ENDDO |
---|
| 62 | DO i=ii_begin,ii_end |
---|
| 63 | CALL request_add_point(ind,i,jj_end+1,req_i1) |
---|
| 64 | ENDDO |
---|
| 65 | |
---|
| 66 | DO j=jj_begin,jj_end+1 |
---|
| 67 | CALL request_add_point(ind,ii_begin-1,j,req_i1) |
---|
| 68 | ENDDO |
---|
| 69 | |
---|
| 70 | DO i=ii_begin,ii_end |
---|
| 71 | CALL request_add_point(ind,i,jj_begin,req_i1) |
---|
| 72 | ENDDO |
---|
| 73 | |
---|
| 74 | DO j=jj_begin,jj_end |
---|
| 75 | CALL request_add_point(ind,ii_end,j,req_i1) |
---|
| 76 | ENDDO |
---|
| 77 | |
---|
| 78 | DO i=ii_begin,ii_end |
---|
| 79 | CALL request_add_point(ind,i,jj_end,req_i1) |
---|
| 80 | ENDDO |
---|
| 81 | |
---|
| 82 | DO j=jj_begin,jj_end |
---|
| 83 | CALL request_add_point(ind,ii_begin,j,req_i1) |
---|
| 84 | ENDDO |
---|
| 85 | |
---|
| 86 | ENDDO |
---|
| 87 | |
---|
| 88 | CALL finalize_request(req_i1) |
---|
| 89 | |
---|
[146] | 90 | CALL create_request(field_u,req_e1_scal) |
---|
[26] | 91 | DO ind=1,ndomain |
---|
| 92 | CALL swap_dimensions(ind) |
---|
| 93 | DO i=ii_begin,ii_end |
---|
[146] | 94 | CALL request_add_point(ind,i,jj_begin-1,req_e1_scal,rup) |
---|
| 95 | CALL request_add_point(ind,i+1,jj_begin-1,req_e1_scal,lup) |
---|
[26] | 96 | ENDDO |
---|
| 97 | |
---|
| 98 | DO j=jj_begin,jj_end |
---|
[146] | 99 | CALL request_add_point(ind,ii_end+1,j,req_e1_scal,left) |
---|
| 100 | CALL request_add_point(ind,ii_end+1,j-1,req_e1_scal,lup) |
---|
[26] | 101 | ENDDO |
---|
| 102 | |
---|
| 103 | DO i=ii_begin,ii_end |
---|
[146] | 104 | CALL request_add_point(ind,i,jj_end+1,req_e1_scal,ldown) |
---|
| 105 | CALL request_add_point(ind,i-1,jj_end+1,req_e1_scal,rdown) |
---|
[26] | 106 | ENDDO |
---|
| 107 | |
---|
| 108 | DO j=jj_begin,jj_end |
---|
[146] | 109 | CALL request_add_point(ind,ii_begin-1,j,req_e1_scal,right) |
---|
| 110 | CALL request_add_point(ind,ii_begin-1,j+1,req_e1_scal,rdown) |
---|
[26] | 111 | ENDDO |
---|
| 112 | |
---|
| 113 | DO i=ii_begin+1,ii_end-1 |
---|
[146] | 114 | CALL request_add_point(ind,i,jj_begin,req_e1_scal,right) |
---|
| 115 | CALL request_add_point(ind,i,jj_end,req_e1_scal,right) |
---|
[26] | 116 | ENDDO |
---|
| 117 | |
---|
| 118 | DO j=jj_begin+1,jj_end-1 |
---|
[146] | 119 | CALL request_add_point(ind,ii_begin,j,req_e1_scal,rup) |
---|
| 120 | CALL request_add_point(ind,ii_end,j,req_e1_scal,rup) |
---|
[26] | 121 | ENDDO |
---|
| 122 | |
---|
[146] | 123 | CALL request_add_point(ind,ii_begin+1,jj_begin,req_e1_scal,left) |
---|
| 124 | CALL request_add_point(ind,ii_begin,jj_begin+1,req_e1_scal,ldown) |
---|
| 125 | CALL request_add_point(ind,ii_begin+1,jj_end,req_e1_scal,left) |
---|
| 126 | CALL request_add_point(ind,ii_end,jj_begin+1,req_e1_scal,ldown) |
---|
[26] | 127 | |
---|
| 128 | ENDDO |
---|
[146] | 129 | |
---|
| 130 | CALL finalize_request(req_e1_scal) |
---|
| 131 | |
---|
| 132 | CALL create_request(field_u,req_e1_vect,.TRUE.) |
---|
| 133 | DO ind=1,ndomain |
---|
| 134 | CALL swap_dimensions(ind) |
---|
| 135 | DO i=ii_begin,ii_end |
---|
| 136 | CALL request_add_point(ind,i,jj_begin-1,req_e1_vect,rup) |
---|
| 137 | CALL request_add_point(ind,i+1,jj_begin-1,req_e1_vect,lup) |
---|
| 138 | ENDDO |
---|
| 139 | |
---|
| 140 | DO j=jj_begin,jj_end |
---|
| 141 | CALL request_add_point(ind,ii_end+1,j,req_e1_vect,left) |
---|
| 142 | CALL request_add_point(ind,ii_end+1,j-1,req_e1_vect,lup) |
---|
| 143 | ENDDO |
---|
| 144 | |
---|
| 145 | DO i=ii_begin,ii_end |
---|
| 146 | CALL request_add_point(ind,i,jj_end+1,req_e1_vect,ldown) |
---|
| 147 | CALL request_add_point(ind,i-1,jj_end+1,req_e1_vect,rdown) |
---|
| 148 | ENDDO |
---|
| 149 | |
---|
| 150 | DO j=jj_begin,jj_end |
---|
| 151 | CALL request_add_point(ind,ii_begin-1,j,req_e1_vect,right) |
---|
| 152 | CALL request_add_point(ind,ii_begin-1,j+1,req_e1_vect,rdown) |
---|
| 153 | ENDDO |
---|
| 154 | |
---|
| 155 | DO i=ii_begin+1,ii_end-1 |
---|
| 156 | CALL request_add_point(ind,i,jj_begin,req_e1_vect,right) |
---|
| 157 | CALL request_add_point(ind,i,jj_end,req_e1_vect,right) |
---|
| 158 | ENDDO |
---|
| 159 | |
---|
| 160 | DO j=jj_begin+1,jj_end-1 |
---|
| 161 | CALL request_add_point(ind,ii_begin,j,req_e1_vect,rup) |
---|
| 162 | CALL request_add_point(ind,ii_end,j,req_e1_vect,rup) |
---|
| 163 | ENDDO |
---|
| 164 | |
---|
| 165 | CALL request_add_point(ind,ii_begin+1,jj_begin,req_e1_vect,left) |
---|
| 166 | CALL request_add_point(ind,ii_begin,jj_begin+1,req_e1_vect,ldown) |
---|
| 167 | CALL request_add_point(ind,ii_begin+1,jj_end,req_e1_vect,left) |
---|
| 168 | CALL request_add_point(ind,ii_end,jj_begin+1,req_e1_vect,ldown) |
---|
| 169 | |
---|
[26] | 170 | |
---|
[146] | 171 | ENDDO |
---|
| 172 | |
---|
| 173 | CALL finalize_request(req_e1_vect) |
---|
| 174 | |
---|
[26] | 175 | END SUBROUTINE init_transfert |
---|
| 176 | |
---|
[146] | 177 | SUBROUTINE create_request(type_field,request,vector) |
---|
[26] | 178 | USE domain_mod |
---|
| 179 | USE field_mod |
---|
| 180 | IMPLICIT NONE |
---|
| 181 | INTEGER :: type_field |
---|
| 182 | TYPE(t_request),POINTER :: request(:) |
---|
[146] | 183 | LOGICAL,OPTIONAL :: vector |
---|
| 184 | |
---|
[26] | 185 | TYPE(t_request),POINTER :: req |
---|
| 186 | TYPE(t_domain),POINTER :: d |
---|
| 187 | INTEGER :: ind |
---|
| 188 | INTEGER :: max_size |
---|
[146] | 189 | |
---|
[26] | 190 | ALLOCATE(request(ndomain)) |
---|
| 191 | |
---|
| 192 | DO ind=1,ndomain |
---|
| 193 | req=>request(ind) |
---|
| 194 | d=>domain(ind) |
---|
| 195 | IF (type_field==field_t) THEN |
---|
| 196 | Max_size=2*(d%iim+2)+2*(d%jjm+2) |
---|
| 197 | ELSE IF (type_field==field_u) THEN |
---|
| 198 | Max_size=3*(2*(d%iim+2)+2*(d%jjm+2)) |
---|
| 199 | ELSE IF (type_field==field_z) THEN |
---|
| 200 | Max_size=2*(2*(d%iim+2)+2*(d%jjm+2)) |
---|
| 201 | ENDIF |
---|
| 202 | |
---|
| 203 | req%type_field=type_field |
---|
| 204 | req%max_size=max_size*2 |
---|
| 205 | req%size=0 |
---|
[146] | 206 | req%vector=.FALSE. |
---|
| 207 | IF (PRESENT(vector)) req%vector=vector |
---|
[26] | 208 | ALLOCATE(req%src_domain(req%max_size)) |
---|
| 209 | ALLOCATE(req%src_ind(req%max_size)) |
---|
| 210 | ALLOCATE(req%target_ind(req%max_size)) |
---|
| 211 | ALLOCATE(req%src_i(req%max_size)) |
---|
| 212 | ALLOCATE(req%src_j(req%max_size)) |
---|
| 213 | ALLOCATE(req%target_i(req%max_size)) |
---|
| 214 | ALLOCATE(req%target_j(req%max_size)) |
---|
[146] | 215 | ALLOCATE(req%target_sign(req%max_size)) |
---|
[26] | 216 | ENDDO |
---|
| 217 | |
---|
| 218 | END SUBROUTINE create_request |
---|
| 219 | |
---|
| 220 | SUBROUTINE reallocate_request(req) |
---|
| 221 | IMPLICIT NONE |
---|
| 222 | TYPE(t_request),POINTER :: req |
---|
| 223 | |
---|
| 224 | INTEGER,POINTER :: src_domain(:) |
---|
| 225 | INTEGER,POINTER :: src_ind(:) |
---|
| 226 | INTEGER,POINTER :: target_ind(:) |
---|
| 227 | INTEGER,POINTER :: src_i(:) |
---|
| 228 | INTEGER,POINTER :: src_j(:) |
---|
| 229 | INTEGER,POINTER :: target_i(:) |
---|
| 230 | INTEGER,POINTER :: target_j(:) |
---|
[146] | 231 | INTEGER,POINTER :: target_sign(:) |
---|
[26] | 232 | |
---|
| 233 | PRINT *,"REALLOCATE_REQUEST" |
---|
| 234 | src_domain=>req%src_domain |
---|
| 235 | src_ind=>req%src_ind |
---|
| 236 | target_ind=>req%target_ind |
---|
| 237 | src_i=>req%src_i |
---|
| 238 | src_j=>req%src_j |
---|
| 239 | target_i=>req%target_i |
---|
| 240 | target_j=>req%target_j |
---|
[146] | 241 | target_sign=>req%target_sign |
---|
[26] | 242 | ! req%max_size=req%max_size*2 |
---|
| 243 | ALLOCATE(req%src_domain(req%max_size*2)) |
---|
| 244 | ALLOCATE(req%src_ind(req%max_size*2)) |
---|
| 245 | ALLOCATE(req%target_ind(req%max_size*2)) |
---|
| 246 | ALLOCATE(req%src_i(req%max_size*2)) |
---|
| 247 | ALLOCATE(req%src_j(req%max_size*2)) |
---|
| 248 | ALLOCATE(req%target_i(req%max_size*2)) |
---|
| 249 | ALLOCATE(req%target_j(req%max_size*2)) |
---|
[146] | 250 | ALLOCATE(req%target_sign(req%max_size*2)) |
---|
[26] | 251 | |
---|
| 252 | req%src_domain(1:req%max_size)=src_domain(:) |
---|
| 253 | req%src_ind(1:req%max_size)=src_ind(:) |
---|
| 254 | req%target_ind(1:req%max_size)=target_ind(:) |
---|
| 255 | req%src_i(1:req%max_size)=src_i(:) |
---|
| 256 | req%src_j(1:req%max_size)=src_j(:) |
---|
| 257 | req%target_i(1:req%max_size)=target_i(:) |
---|
| 258 | req%target_j(1:req%max_size)=target_j(:) |
---|
[146] | 259 | req%target_sign(1:req%max_size)=target_sign(:) |
---|
[26] | 260 | |
---|
| 261 | req%max_size=req%max_size*2 |
---|
| 262 | |
---|
| 263 | DEALLOCATE(src_domain) |
---|
| 264 | DEALLOCATE(src_ind) |
---|
| 265 | DEALLOCATE(target_ind) |
---|
| 266 | DEALLOCATE(src_i) |
---|
| 267 | DEALLOCATE(src_j) |
---|
| 268 | DEALLOCATE(target_i) |
---|
| 269 | DEALLOCATE(target_j) |
---|
[146] | 270 | DEALLOCATE(target_sign) |
---|
[26] | 271 | |
---|
| 272 | END SUBROUTINE reallocate_request |
---|
| 273 | |
---|
| 274 | |
---|
| 275 | SUBROUTINE request_add_point(ind,i,j,request,pos) |
---|
| 276 | USE domain_mod |
---|
| 277 | USE field_mod |
---|
| 278 | IMPLICIT NONE |
---|
| 279 | INTEGER,INTENT(IN) :: ind |
---|
| 280 | INTEGER,INTENT(IN) :: i |
---|
| 281 | INTEGER,INTENT(IN) :: j |
---|
| 282 | TYPE(t_request),POINTER :: request(:) |
---|
| 283 | INTEGER,INTENT(IN),OPTIONAL :: pos |
---|
| 284 | |
---|
| 285 | INTEGER :: src_domain |
---|
| 286 | INTEGER :: src_iim,src_i,src_j,src_n,src_pos,src_delta |
---|
| 287 | TYPE(t_request),POINTER :: req |
---|
| 288 | TYPE(t_domain),POINTER :: d |
---|
| 289 | |
---|
| 290 | req=>request(ind) |
---|
| 291 | d=>domain(ind) |
---|
| 292 | |
---|
| 293 | IF (req%max_size==req%size) CALL reallocate_request(req) |
---|
| 294 | req%size=req%size+1 |
---|
| 295 | IF (req%type_field==field_t) THEN |
---|
| 296 | src_domain=domain(ind)%assign_domain(i,j) |
---|
| 297 | src_iim=domain_glo(src_domain)%iim |
---|
| 298 | src_i=domain(ind)%assign_i(i,j) |
---|
| 299 | src_j=domain(ind)%assign_j(i,j) |
---|
| 300 | |
---|
| 301 | req%target_ind(req%size)=(j-1)*d%iim+i |
---|
[146] | 302 | req%target_sign(req%size)=1 |
---|
[26] | 303 | req%src_domain(req%size)=src_domain |
---|
| 304 | req%src_ind(req%size)=(src_j-1)*src_iim+src_i |
---|
| 305 | ELSE IF (req%type_field==field_u) THEN |
---|
| 306 | IF (.NOT. PRESENT(pos)) STOP 'argument request_add_point non conforme' |
---|
| 307 | |
---|
| 308 | src_domain=domain(ind)%edge_assign_domain(pos-1,i,j) |
---|
| 309 | src_iim=domain_glo(src_domain)%iim |
---|
| 310 | src_i=domain(ind)%edge_assign_i(pos-1,i,j) |
---|
| 311 | src_j=domain(ind)%edge_assign_j(pos-1,i,j) |
---|
| 312 | src_n=(src_j-1)*src_iim+src_i |
---|
| 313 | src_delta=domain(ind)%delta(i,j) |
---|
| 314 | src_pos=domain(ind)%edge_assign_pos(pos-1,i,j)+1 |
---|
| 315 | |
---|
| 316 | req%target_ind(req%size)=(j-1)*d%iim+i+d%u_pos(pos) |
---|
[146] | 317 | |
---|
| 318 | req%target_sign(req%size)= 1 |
---|
| 319 | IF (req%vector) req%target_sign(req%size)= domain(ind)%edge_assign_sign(pos-1,i,j) |
---|
| 320 | |
---|
[26] | 321 | req%src_domain(req%size)=src_domain |
---|
| 322 | req%src_ind(req%size)=src_n+domain_glo(src_domain)%u_pos(src_pos) |
---|
| 323 | |
---|
| 324 | ELSE IF (req%type_field==field_z) THEN |
---|
| 325 | IF (.NOT. PRESENT(pos)) STOP 'argument request_add_point non conforme' |
---|
| 326 | |
---|
| 327 | src_domain=domain(ind)%assign_domain(i,j) |
---|
| 328 | src_iim=domain_glo(src_domain)%iim |
---|
| 329 | src_i=domain(ind)%assign_i(i,j) |
---|
| 330 | src_j=domain(ind)%assign_j(i,j) |
---|
| 331 | src_n=(src_j-1)*src_iim+src_i |
---|
| 332 | src_delta=domain(ind)%delta(i,j) |
---|
| 333 | |
---|
| 334 | src_pos=MOD(pos-1+src_delta+6,6)+1 |
---|
| 335 | |
---|
| 336 | req%target_ind(req%size)=(j-1)*d%iim+i+d%z_pos(pos) |
---|
[146] | 337 | req%target_sign(req%size)=1 |
---|
[26] | 338 | req%src_domain(req%size)=src_domain |
---|
| 339 | req%src_ind(req%size)=src_n+domain_glo(src_domain)%z_pos(src_pos) |
---|
| 340 | ENDIF |
---|
| 341 | END SUBROUTINE request_add_point |
---|
| 342 | |
---|
| 343 | |
---|
| 344 | SUBROUTINE Finalize_request(request) |
---|
| 345 | USE mpipara |
---|
| 346 | USE domain_mod |
---|
| 347 | USE mpi_mod |
---|
| 348 | IMPLICIT NONE |
---|
| 349 | TYPE(t_request),POINTER :: request(:) |
---|
| 350 | TYPE(t_request),POINTER :: req |
---|
| 351 | INTEGER :: nb_domain_recv(0:mpi_size-1) |
---|
| 352 | INTEGER :: nb_domain_send(0:mpi_size-1) |
---|
| 353 | INTEGER :: nb_data_domain_recv(ndomain_glo) |
---|
| 354 | INTEGER :: list_domain_recv(ndomain_glo) |
---|
| 355 | INTEGER,ALLOCATABLE :: list_domain_send(:) |
---|
| 356 | INTEGER :: list_domain(ndomain) |
---|
| 357 | |
---|
| 358 | INTEGER :: rank,i,j |
---|
| 359 | INTEGER :: size,ind_glo,ind_loc |
---|
| 360 | INTEGER :: isend, irecv, ireq, nreq |
---|
| 361 | INTEGER, ALLOCATABLE :: mpi_req(:) |
---|
| 362 | INTEGER, ALLOCATABLE :: status(:,:) |
---|
| 363 | |
---|
| 364 | IF (.NOT. using_mpi) RETURN |
---|
| 365 | |
---|
| 366 | DO ind_loc=1,ndomain |
---|
| 367 | req=>request(ind_loc) |
---|
| 368 | |
---|
| 369 | nb_data_domain_recv(:) = 0 |
---|
| 370 | nb_domain_recv(:) = 0 |
---|
| 371 | |
---|
| 372 | DO i=1,req%size |
---|
| 373 | ind_glo=req%src_domain(i) |
---|
| 374 | nb_data_domain_recv(ind_glo)=nb_data_domain_recv(ind_glo)+1 |
---|
| 375 | ENDDO |
---|
| 376 | |
---|
| 377 | DO ind_glo=1,ndomain_glo |
---|
| 378 | IF ( nb_data_domain_recv(ind_glo) > 0 ) nb_domain_recv(domglo_rank(ind_glo))=nb_domain_recv(domglo_rank(ind_glo))+1 |
---|
| 379 | ENDDO |
---|
| 380 | |
---|
| 381 | req%nrecv=sum(nb_domain_recv(:)) |
---|
| 382 | ALLOCATE(req%recv(req%nrecv)) |
---|
| 383 | |
---|
| 384 | irecv=0 |
---|
| 385 | DO ind_glo=1,ndomain_glo |
---|
| 386 | IF (nb_data_domain_recv(ind_glo)>0) THEN |
---|
| 387 | irecv=irecv+1 |
---|
| 388 | list_domain_recv(ind_glo)=irecv |
---|
| 389 | req%recv(irecv)%rank=domglo_rank(ind_glo) |
---|
| 390 | req%recv(irecv)%size=nb_data_domain_recv(ind_glo) |
---|
| 391 | req%recv(irecv)%domain=domglo_loc_ind(ind_glo) |
---|
| 392 | ALLOCATE(req%recv(irecv)%value(req%recv(irecv)%size)) |
---|
[146] | 393 | ALLOCATE(req%recv(irecv)%sign(req%recv(irecv)%size)) |
---|
[26] | 394 | ALLOCATE(req%recv(irecv)%buffer(req%recv(irecv)%size)) |
---|
| 395 | ENDIF |
---|
| 396 | ENDDO |
---|
| 397 | |
---|
| 398 | req%recv(:)%size=0 |
---|
| 399 | irecv=0 |
---|
| 400 | DO i=1,req%size |
---|
| 401 | irecv=list_domain_recv(req%src_domain(i)) |
---|
| 402 | req%recv(irecv)%size=req%recv(irecv)%size+1 |
---|
| 403 | size=req%recv(irecv)%size |
---|
| 404 | req%recv(irecv)%value(size)=req%src_ind(i) |
---|
| 405 | req%recv(irecv)%buffer(size)=req%target_ind(i) |
---|
[146] | 406 | req%recv(irecv)%sign(size)=req%target_sign(i) |
---|
[26] | 407 | ENDDO |
---|
| 408 | ENDDO |
---|
| 409 | |
---|
| 410 | nb_domain_recv(:) = 0 |
---|
| 411 | DO ind_loc=1,ndomain |
---|
| 412 | req=>request(ind_loc) |
---|
| 413 | |
---|
| 414 | DO irecv=1,req%nrecv |
---|
| 415 | rank=req%recv(irecv)%rank |
---|
| 416 | nb_domain_recv(rank)=nb_domain_recv(rank)+1 |
---|
| 417 | ENDDO |
---|
| 418 | ENDDO |
---|
| 419 | |
---|
| 420 | CALL MPI_Alltoall(nb_domain_recv,1,MPI_INTEGER,nb_domain_send,1,MPI_INTEGER,comm_icosa,ierr) |
---|
| 421 | |
---|
| 422 | |
---|
| 423 | ALLOCATE(list_domain_send(sum(nb_domain_send))) |
---|
| 424 | |
---|
| 425 | nreq=sum(nb_domain_recv(:))+sum(nb_domain_send(:)) |
---|
| 426 | ALLOCATE(mpi_req(nreq)) |
---|
| 427 | ALLOCATE(status(MPI_STATUS_SIZE,nreq)) |
---|
| 428 | |
---|
| 429 | ireq=0 |
---|
| 430 | DO ind_loc=1,ndomain |
---|
| 431 | req=>request(ind_loc) |
---|
| 432 | DO irecv=1,req%nrecv |
---|
| 433 | ireq=ireq+1 |
---|
| 434 | CALL MPI_ISEND(req%recv(irecv)%domain,1,MPI_INTEGER,req%recv(irecv)%rank,0,comm_icosa, mpi_req(ireq),ierr) |
---|
| 435 | ENDDO |
---|
| 436 | ENDDO |
---|
| 437 | |
---|
| 438 | j=0 |
---|
| 439 | DO rank=0,mpi_size-1 |
---|
| 440 | DO i=1,nb_domain_send(rank) |
---|
| 441 | j=j+1 |
---|
| 442 | ireq=ireq+1 |
---|
| 443 | CALL MPI_IRECV(list_domain_send(j),1,MPI_INTEGER,rank,0,comm_icosa, mpi_req(ireq),ierr) |
---|
| 444 | ENDDO |
---|
| 445 | ENDDO |
---|
| 446 | |
---|
| 447 | CALL MPI_WAITALL(nreq,mpi_req,status,ierr) |
---|
| 448 | |
---|
| 449 | list_domain(:)=0 |
---|
| 450 | DO i=1,sum(nb_domain_send) |
---|
| 451 | ind_loc=list_domain_send(i) |
---|
| 452 | list_domain(ind_loc)=list_domain(ind_loc)+1 |
---|
| 453 | ENDDO |
---|
| 454 | |
---|
| 455 | DO ind_loc=1,ndomain |
---|
| 456 | req=>request(ind_loc) |
---|
| 457 | req%nsend=list_domain(ind_loc) |
---|
| 458 | ALLOCATE(req%send(req%nsend)) |
---|
| 459 | ENDDO |
---|
| 460 | |
---|
| 461 | ireq=0 |
---|
| 462 | DO ind_loc=1,ndomain |
---|
| 463 | req=>request(ind_loc) |
---|
| 464 | |
---|
| 465 | DO irecv=1,req%nrecv |
---|
| 466 | ireq=ireq+1 |
---|
| 467 | CALL MPI_ISEND(mpi_rank,1,MPI_INTEGER,req%recv(irecv)%rank,req%recv(irecv)%domain,comm_icosa, mpi_req(ireq),ierr) |
---|
| 468 | ENDDO |
---|
| 469 | |
---|
| 470 | DO isend=1,req%nsend |
---|
| 471 | ireq=ireq+1 |
---|
| 472 | CALL MPI_IRECV(req%send(isend)%rank,1,MPI_INTEGER,MPI_ANY_SOURCE,ind_loc,comm_icosa, mpi_req(ireq),ierr) |
---|
| 473 | ENDDO |
---|
| 474 | ENDDO |
---|
| 475 | |
---|
| 476 | CALL MPI_WAITALL(nreq,mpi_req,status,ierr) |
---|
| 477 | CALL MPI_BARRIER(comm_icosa,ierr) |
---|
| 478 | |
---|
| 479 | ireq=0 |
---|
| 480 | DO ind_loc=1,ndomain |
---|
| 481 | req=>request(ind_loc) |
---|
| 482 | |
---|
| 483 | DO irecv=1,req%nrecv |
---|
| 484 | ireq=ireq+1 |
---|
| 485 | CALL MPI_ISEND(req%recv(irecv)%size,1,MPI_INTEGER,req%recv(irecv)%rank,req%recv(irecv)%domain,comm_icosa, mpi_req(ireq),ierr) |
---|
| 486 | ENDDO |
---|
| 487 | |
---|
| 488 | DO isend=1,req%nsend |
---|
| 489 | ireq=ireq+1 |
---|
| 490 | CALL MPI_IRECV(req%send(isend)%size,1,MPI_INTEGER,req%send(isend)%rank,ind_loc,comm_icosa, mpi_req(ireq),ierr) |
---|
| 491 | ENDDO |
---|
| 492 | ENDDO |
---|
| 493 | |
---|
| 494 | CALL MPI_WAITALL(nreq,mpi_req,status,ierr) |
---|
| 495 | |
---|
| 496 | ireq=0 |
---|
| 497 | DO ind_loc=1,ndomain |
---|
| 498 | req=>request(ind_loc) |
---|
| 499 | |
---|
| 500 | DO irecv=1,req%nrecv |
---|
| 501 | ireq=ireq+1 |
---|
[44] | 502 | CALL MPI_ISEND(req%recv(irecv)%value,req%recv(irecv)%size,MPI_INTEGER,& |
---|
| 503 | req%recv(irecv)%rank,req%recv(irecv)%domain,comm_icosa, mpi_req(ireq),ierr) |
---|
[26] | 504 | ENDDO |
---|
| 505 | |
---|
| 506 | DO isend=1,req%nsend |
---|
| 507 | ireq=ireq+1 |
---|
| 508 | ALLOCATE(req%send(isend)%value(req%send(isend)%size)) |
---|
[44] | 509 | CALL MPI_IRECV(req%send(isend)%value,req%send(isend)%size,MPI_INTEGER,& |
---|
| 510 | req%send(isend)%rank,ind_loc,comm_icosa, mpi_req(ireq),ierr) |
---|
[26] | 511 | ENDDO |
---|
| 512 | ENDDO |
---|
| 513 | |
---|
| 514 | CALL MPI_WAITALL(nreq,mpi_req,status,ierr) |
---|
| 515 | |
---|
| 516 | DO ind_loc=1,ndomain |
---|
| 517 | req=>request(ind_loc) |
---|
| 518 | |
---|
| 519 | DO irecv=1,req%nrecv |
---|
| 520 | req%recv(irecv)%value(:)=req%recv(irecv)%buffer(:) |
---|
[146] | 521 | req%recv(irecv)%sign(:) =req%recv(irecv)%sign(:) |
---|
[26] | 522 | DEALLOCATE(req%recv(irecv)%buffer) |
---|
| 523 | ENDDO |
---|
| 524 | ENDDO |
---|
| 525 | |
---|
| 526 | |
---|
[146] | 527 | END SUBROUTINE Finalize_request |
---|
[26] | 528 | |
---|
| 529 | |
---|
| 530 | SUBROUTINE transfert_request_mpi(field,request) |
---|
| 531 | USE field_mod |
---|
| 532 | USE domain_mod |
---|
| 533 | USE mpi_mod |
---|
| 534 | USE mpipara |
---|
| 535 | IMPLICIT NONE |
---|
| 536 | TYPE(t_field),POINTER :: field(:) |
---|
| 537 | TYPE(t_request),POINTER :: request(:) |
---|
| 538 | REAL(rstd),POINTER :: rval2d(:) |
---|
| 539 | REAL(rstd),POINTER :: rval3d(:,:) |
---|
| 540 | REAL(rstd),POINTER :: rval4d(:,:,:) |
---|
| 541 | REAL(rstd),POINTER :: buffer_r2(:) |
---|
| 542 | REAL(rstd),POINTER :: buffer_r3(:,:) |
---|
| 543 | REAL(rstd),POINTER :: buffer_r4(:,:,:) |
---|
| 544 | INTEGER,POINTER :: value(:) |
---|
[146] | 545 | INTEGER,POINTER :: sgn(:) |
---|
[26] | 546 | TYPE(ARRAY),POINTER :: recv,send |
---|
| 547 | TYPE(t_request),POINTER :: req |
---|
| 548 | INTEGER, ALLOCATABLE :: mpi_req(:) |
---|
| 549 | INTEGER, ALLOCATABLE :: status(:,:) |
---|
| 550 | INTEGER :: irecv,isend |
---|
| 551 | INTEGER :: ireq,nreq |
---|
| 552 | INTEGER :: ind,n |
---|
[74] | 553 | INTEGER :: dim3,dim4 |
---|
[26] | 554 | |
---|
| 555 | IF (field(1)%data_type==type_real) THEN |
---|
| 556 | IF (field(1)%ndim==2) THEN |
---|
| 557 | |
---|
| 558 | nreq=sum(request(:)%nsend)+sum(request(:)%nrecv) |
---|
| 559 | ALLOCATE(mpi_req(nreq)) |
---|
| 560 | ALLOCATE(status(MPI_STATUS_SIZE,nreq)) |
---|
| 561 | |
---|
| 562 | ireq=0 |
---|
| 563 | DO ind=1,ndomain |
---|
| 564 | rval2d=>field(ind)%rval2d |
---|
| 565 | |
---|
| 566 | req=>request(ind) |
---|
| 567 | DO isend=1,req%nsend |
---|
| 568 | send=>req%send(isend) |
---|
| 569 | |
---|
| 570 | ALLOCATE(send%buffer_r2(send%size)) |
---|
| 571 | buffer_r2=>send%buffer_r2 |
---|
| 572 | value=>send%value |
---|
| 573 | DO n=1,send%size |
---|
| 574 | buffer_r2(n)=rval2d(value(n)) |
---|
| 575 | ENDDO |
---|
| 576 | |
---|
| 577 | ireq=ireq+1 |
---|
| 578 | CALL MPI_ISEND(buffer_r2,send%size,MPI_REAL8,send%rank,ind,comm_icosa, mpi_req(ireq),ierr) |
---|
| 579 | ENDDO |
---|
| 580 | |
---|
| 581 | DO irecv=1,req%nrecv |
---|
| 582 | recv=>req%recv(irecv) |
---|
| 583 | ALLOCATE(recv%buffer_r2(recv%size)) |
---|
| 584 | |
---|
| 585 | ireq=ireq+1 |
---|
| 586 | CALL MPI_IRECV(recv%buffer_r2,recv%size,MPI_REAL8,recv%rank,recv%domain,comm_icosa, mpi_req(ireq),ierr) |
---|
| 587 | ENDDO |
---|
| 588 | |
---|
| 589 | ENDDO |
---|
| 590 | |
---|
| 591 | CALL MPI_WAITALL(nreq,mpi_req,status,ierr) |
---|
| 592 | ! DO ind=1,ndomain |
---|
| 593 | ! field(ind)%rval2d(:)=0. |
---|
| 594 | ! ENDDO |
---|
| 595 | |
---|
| 596 | DO ind=1,ndomain |
---|
| 597 | rval2d=>field(ind)%rval2d |
---|
| 598 | |
---|
| 599 | req=>request(ind) |
---|
| 600 | DO isend=1,req%nsend |
---|
| 601 | send=>req%send(isend) |
---|
| 602 | DEALLOCATE(send%buffer_r2) |
---|
| 603 | ENDDO |
---|
| 604 | |
---|
| 605 | DO irecv=1,req%nrecv |
---|
| 606 | recv=>req%recv(irecv) |
---|
| 607 | buffer_r2=>recv%buffer_r2 |
---|
| 608 | value=>recv%value |
---|
[146] | 609 | sgn=>recv%sign |
---|
[26] | 610 | DO n=1,recv%size |
---|
[146] | 611 | rval2d(value(n))=buffer_r2(n)*sgn(n) |
---|
[26] | 612 | ENDDO |
---|
| 613 | DEALLOCATE(recv%buffer_r2) |
---|
| 614 | ENDDO |
---|
| 615 | |
---|
| 616 | ENDDO |
---|
| 617 | |
---|
| 618 | |
---|
| 619 | ELSE IF (field(1)%ndim==3) THEN |
---|
| 620 | |
---|
| 621 | nreq=sum(request(:)%nsend)+sum(request(:)%nrecv) |
---|
| 622 | ALLOCATE(mpi_req(nreq)) |
---|
| 623 | ALLOCATE(status(MPI_STATUS_SIZE,nreq)) |
---|
| 624 | |
---|
| 625 | ireq=0 |
---|
| 626 | DO ind=1,ndomain |
---|
| 627 | dim3=size(field(ind)%rval3d,2) |
---|
| 628 | rval3d=>field(ind)%rval3d |
---|
| 629 | |
---|
| 630 | req=>request(ind) |
---|
| 631 | DO isend=1,req%nsend |
---|
| 632 | send=>req%send(isend) |
---|
| 633 | |
---|
| 634 | ALLOCATE(send%buffer_r3(send%size,dim3)) |
---|
| 635 | buffer_r3=>send%buffer_r3 |
---|
| 636 | value=>send%value |
---|
| 637 | DO n=1,send%size |
---|
| 638 | buffer_r3(n,:)=rval3d(value(n),:) |
---|
| 639 | ENDDO |
---|
| 640 | |
---|
| 641 | ireq=ireq+1 |
---|
| 642 | CALL MPI_ISEND(buffer_r3,send%size*dim3,MPI_REAL8,send%rank,ind,comm_icosa, mpi_req(ireq),ierr) |
---|
| 643 | ENDDO |
---|
| 644 | |
---|
| 645 | DO irecv=1,req%nrecv |
---|
| 646 | recv=>req%recv(irecv) |
---|
| 647 | ALLOCATE(recv%buffer_r3(recv%size,dim3)) |
---|
| 648 | |
---|
| 649 | ireq=ireq+1 |
---|
| 650 | CALL MPI_IRECV(recv%buffer_r3,recv%size*dim3,MPI_REAL8,recv%rank,recv%domain,comm_icosa, mpi_req(ireq),ierr) |
---|
| 651 | ENDDO |
---|
| 652 | |
---|
| 653 | ENDDO |
---|
| 654 | |
---|
| 655 | CALL MPI_WAITALL(nreq,mpi_req,status,ierr) |
---|
| 656 | ! DO ind=1,ndomain |
---|
| 657 | ! field(ind)%rval2d(:)=0. |
---|
| 658 | ! ENDDO |
---|
| 659 | |
---|
| 660 | DO ind=1,ndomain |
---|
| 661 | rval3d=>field(ind)%rval3d |
---|
| 662 | |
---|
| 663 | req=>request(ind) |
---|
| 664 | DO isend=1,req%nsend |
---|
| 665 | send=>req%send(isend) |
---|
| 666 | DEALLOCATE(send%buffer_r3) |
---|
| 667 | ENDDO |
---|
| 668 | |
---|
| 669 | DO irecv=1,req%nrecv |
---|
| 670 | recv=>req%recv(irecv) |
---|
| 671 | buffer_r3=>recv%buffer_r3 |
---|
| 672 | value=>recv%value |
---|
[146] | 673 | sgn=>recv%sign |
---|
[26] | 674 | DO n=1,recv%size |
---|
[146] | 675 | rval3d(value(n),:)=buffer_r3(n,:)*sgn(n) |
---|
[26] | 676 | ENDDO |
---|
| 677 | DEALLOCATE(recv%buffer_r3) |
---|
| 678 | ENDDO |
---|
| 679 | |
---|
| 680 | ENDDO |
---|
[74] | 681 | |
---|
| 682 | ELSE IF (field(1)%ndim==4) THEN |
---|
[26] | 683 | |
---|
[74] | 684 | nreq=sum(request(:)%nsend)+sum(request(:)%nrecv) |
---|
| 685 | ALLOCATE(mpi_req(nreq)) |
---|
| 686 | ALLOCATE(status(MPI_STATUS_SIZE,nreq)) |
---|
| 687 | |
---|
| 688 | ireq=0 |
---|
| 689 | DO ind=1,ndomain |
---|
| 690 | dim3=size(field(ind)%rval4d,2) |
---|
| 691 | dim4=size(field(ind)%rval4d,3) |
---|
| 692 | rval4d=>field(ind)%rval4d |
---|
| 693 | |
---|
| 694 | req=>request(ind) |
---|
| 695 | DO isend=1,req%nsend |
---|
| 696 | send=>req%send(isend) |
---|
| 697 | |
---|
| 698 | ALLOCATE(send%buffer_r4(send%size,dim3,dim4)) |
---|
| 699 | buffer_r4=>send%buffer_r4 |
---|
| 700 | value=>send%value |
---|
| 701 | DO n=1,send%size |
---|
| 702 | buffer_r4(n,:,:)=rval4d(value(n),:,:) |
---|
| 703 | ENDDO |
---|
| 704 | |
---|
| 705 | ireq=ireq+1 |
---|
| 706 | CALL MPI_ISEND(buffer_r4,send%size*dim3*dim4,MPI_REAL8,send%rank,ind,comm_icosa, mpi_req(ireq),ierr) |
---|
| 707 | ENDDO |
---|
| 708 | |
---|
| 709 | DO irecv=1,req%nrecv |
---|
| 710 | recv=>req%recv(irecv) |
---|
| 711 | ALLOCATE(recv%buffer_r4(recv%size,dim3,dim4)) |
---|
| 712 | |
---|
| 713 | ireq=ireq+1 |
---|
| 714 | CALL MPI_IRECV(recv%buffer_r4,recv%size*dim3*dim4,MPI_REAL8,recv%rank,recv%domain,comm_icosa, mpi_req(ireq),ierr) |
---|
| 715 | ENDDO |
---|
| 716 | |
---|
| 717 | ENDDO |
---|
| 718 | |
---|
| 719 | CALL MPI_WAITALL(nreq,mpi_req,status,ierr) |
---|
| 720 | ! DO ind=1,ndomain |
---|
| 721 | ! field(ind)%rval2d(:)=0. |
---|
| 722 | ! ENDDO |
---|
| 723 | |
---|
| 724 | DO ind=1,ndomain |
---|
| 725 | rval4d=>field(ind)%rval4d |
---|
| 726 | |
---|
| 727 | req=>request(ind) |
---|
| 728 | DO isend=1,req%nsend |
---|
| 729 | send=>req%send(isend) |
---|
| 730 | DEALLOCATE(send%buffer_r4) |
---|
| 731 | ENDDO |
---|
| 732 | |
---|
| 733 | DO irecv=1,req%nrecv |
---|
| 734 | recv=>req%recv(irecv) |
---|
| 735 | buffer_r4=>recv%buffer_r4 |
---|
| 736 | value=>recv%value |
---|
[146] | 737 | sgn=>recv%sign |
---|
[74] | 738 | DO n=1,recv%size |
---|
[146] | 739 | rval4d(value(n),:,:)=buffer_r4(n,:,:)*sgn(n) |
---|
[74] | 740 | ENDDO |
---|
| 741 | DEALLOCATE(recv%buffer_r4) |
---|
| 742 | ENDDO |
---|
| 743 | |
---|
| 744 | ENDDO |
---|
[26] | 745 | |
---|
[74] | 746 | ENDIF |
---|
[26] | 747 | |
---|
| 748 | ENDIF |
---|
| 749 | |
---|
| 750 | END SUBROUTINE transfert_request_mpi |
---|
| 751 | |
---|
| 752 | SUBROUTINE transfert_request(field,request) |
---|
| 753 | USE field_mod |
---|
| 754 | USE domain_mod |
---|
| 755 | IMPLICIT NONE |
---|
| 756 | TYPE(t_field),POINTER :: field(:) |
---|
| 757 | TYPE(t_request),POINTER :: request(:) |
---|
| 758 | REAL(rstd),POINTER :: rval2d(:) |
---|
| 759 | REAL(rstd),POINTER :: rval3d(:,:) |
---|
| 760 | REAL(rstd),POINTER :: rval4d(:,:,:) |
---|
| 761 | INTEGER :: ind |
---|
| 762 | TYPE(t_request),POINTER :: req |
---|
| 763 | INTEGER :: n |
---|
| 764 | REAL(rstd) :: var1,var2 |
---|
| 765 | |
---|
| 766 | DO ind=1,ndomain |
---|
| 767 | req=>request(ind) |
---|
| 768 | rval2d=>field(ind)%rval2d |
---|
| 769 | rval3d=>field(ind)%rval3d |
---|
| 770 | rval4d=>field(ind)%rval4d |
---|
| 771 | |
---|
| 772 | IF (field(ind)%data_type==type_real) THEN |
---|
| 773 | IF (field(ind)%ndim==2) THEN |
---|
| 774 | DO n=1,req%size |
---|
[146] | 775 | rval2d(req%target_ind(n))=field(req%src_domain(n))%rval2d(req%src_ind(n))*req%target_sign(n) |
---|
[26] | 776 | ENDDO |
---|
| 777 | ELSE IF (field(ind)%ndim==3) THEN |
---|
| 778 | DO n=1,req%size |
---|
[146] | 779 | rval3d(req%target_ind(n),:)=field(req%src_domain(n))%rval3d(req%src_ind(n),:)*req%target_sign(n) |
---|
[26] | 780 | ENDDO |
---|
| 781 | ELSE IF (field(ind)%ndim==4) THEN |
---|
| 782 | DO n=1,req%size |
---|
[146] | 783 | rval4d(req%target_ind(n),:,:)=field(req%src_domain(n))%rval4d(req%src_ind(n),:,:)*req%target_sign(n) |
---|
[26] | 784 | ENDDO |
---|
| 785 | ENDIF |
---|
| 786 | ENDIF |
---|
| 787 | |
---|
| 788 | ENDDO |
---|
| 789 | |
---|
| 790 | END SUBROUTINE transfert_request |
---|
| 791 | |
---|
| 792 | |
---|
| 793 | SUBROUTINE gather_field(field_loc,field_glo) |
---|
| 794 | USE field_mod |
---|
| 795 | USE domain_mod |
---|
| 796 | USE mpi_mod |
---|
| 797 | USE mpipara |
---|
| 798 | IMPLICIT NONE |
---|
| 799 | TYPE(t_field),POINTER :: field_loc(:) |
---|
| 800 | TYPE(t_field),POINTER :: field_glo(:) |
---|
| 801 | INTEGER, ALLOCATABLE :: mpi_req(:) |
---|
| 802 | INTEGER, ALLOCATABLE :: status(:,:) |
---|
| 803 | INTEGER :: ireq,nreq |
---|
| 804 | INTEGER :: ind_glo,ind_loc |
---|
| 805 | |
---|
| 806 | IF (.NOT. using_mpi) THEN |
---|
| 807 | |
---|
| 808 | DO ind_loc=1,ndomain |
---|
| 809 | IF (field_loc(ind_loc)%ndim==2) field_glo(ind_loc)%rval2d=field_loc(ind_loc)%rval2d |
---|
| 810 | IF (field_loc(ind_loc)%ndim==3) field_glo(ind_loc)%rval3d=field_loc(ind_loc)%rval3d |
---|
| 811 | IF (field_loc(ind_loc)%ndim==4) field_glo(ind_loc)%rval4d=field_loc(ind_loc)%rval4d |
---|
| 812 | ENDDO |
---|
| 813 | |
---|
| 814 | ELSE |
---|
| 815 | |
---|
| 816 | nreq=ndomain |
---|
| 817 | IF (mpi_rank==0) nreq=nreq+ndomain_glo |
---|
| 818 | ALLOCATE(mpi_req(nreq)) |
---|
| 819 | ALLOCATE(status(MPI_STATUS_SIZE,nreq)) |
---|
| 820 | |
---|
| 821 | |
---|
| 822 | ireq=0 |
---|
| 823 | IF (mpi_rank==0) THEN |
---|
| 824 | DO ind_glo=1,ndomain_glo |
---|
| 825 | ireq=ireq+1 |
---|
| 826 | |
---|
| 827 | IF (field_glo(ind_glo)%ndim==2) THEN |
---|
| 828 | CALL MPI_IRECV(field_glo(ind_glo)%rval2d,size(field_glo(ind_glo)%rval2d) , MPI_REAL8 , & |
---|
| 829 | domglo_rank(ind_glo),domglo_loc_ind(ind_glo), comm_icosa, mpi_req(ireq), ierr) |
---|
| 830 | |
---|
| 831 | ELSE IF (field_glo(ind_glo)%ndim==3) THEN |
---|
| 832 | CALL MPI_IRECV(field_glo(ind_glo)%rval3d,size(field_glo(ind_glo)%rval3d) , MPI_REAL8 , & |
---|
| 833 | domglo_rank(ind_glo),domglo_loc_ind(ind_glo), comm_icosa, mpi_req(ireq), ierr) |
---|
| 834 | |
---|
| 835 | ELSE IF (field_glo(ind_glo)%ndim==4) THEN |
---|
[31] | 836 | CALL MPI_IRECV(field_glo(ind_glo)%rval4d,size(field_glo(ind_glo)%rval4d) , MPI_REAL8 , & |
---|
[26] | 837 | domglo_rank(ind_glo),domglo_loc_ind(ind_glo), comm_icosa, mpi_req(ireq), ierr) |
---|
| 838 | ENDIF |
---|
| 839 | |
---|
| 840 | ENDDO |
---|
| 841 | ENDIF |
---|
| 842 | |
---|
| 843 | DO ind_loc=1,ndomain |
---|
| 844 | ireq=ireq+1 |
---|
| 845 | |
---|
| 846 | IF (field_loc(ind_loc)%ndim==2) THEN |
---|
| 847 | CALL MPI_ISEND(field_loc(ind_loc)%rval2d,size(field_loc(ind_loc)%rval2d) , MPI_REAL8 , & |
---|
| 848 | 0, ind_loc, comm_icosa, mpi_req(ireq), ierr) |
---|
| 849 | ELSE IF (field_loc(ind_loc)%ndim==3) THEN |
---|
| 850 | CALL MPI_ISEND(field_loc(ind_loc)%rval3d,size(field_loc(ind_loc)%rval3d) , MPI_REAL8 , & |
---|
| 851 | 0, ind_loc, comm_icosa, mpi_req(ireq), ierr) |
---|
| 852 | ELSE IF (field_loc(ind_loc)%ndim==4) THEN |
---|
| 853 | CALL MPI_ISEND(field_loc(ind_loc)%rval4d,size(field_loc(ind_loc)%rval4d) , MPI_REAL8 , & |
---|
| 854 | 0, ind_loc, comm_icosa, mpi_req(ireq), ierr) |
---|
| 855 | ENDIF |
---|
| 856 | |
---|
| 857 | ENDDO |
---|
| 858 | |
---|
| 859 | CALL MPI_WAITALL(nreq,mpi_req,status,ierr) |
---|
| 860 | |
---|
| 861 | ENDIF |
---|
| 862 | |
---|
| 863 | END SUBROUTINE gather_field |
---|
| 864 | |
---|
| 865 | |
---|
| 866 | END MODULE transfert_mpi_mod |
---|
| 867 | |
---|
| 868 | |
---|
| 869 | |
---|
| 870 | |
---|
| 871 | |
---|