File: /nfs/home/0/users/jenkins/mfix.git/model/dmp_modules/sendrecv_mod.f
1
2
3
4
5
6
7
8
9
10 module sendrecv
11
12
13
14
15 use parallel_mpi
16 use debug
17 use geometry
18 use compar
19 use indices
20 use functions
21 implicit none
22
23
24 logical,parameter :: localfunc=.false.
25 logical,parameter :: use_persistent_message=.true.
26
27
28 integer, pointer, dimension(:) :: &
29 recvproc1, recvtag1, xrecv1, recvijk1, &
30 sendproc1, sendtag1, xsend1, sendijk1, &
31 recvproc2, recvtag2, xrecv2, recvijk2, &
32 sendproc2, sendtag2, xsend2, sendijk2
33
34 integer,pointer, dimension(:) :: &
35 send_persistent_request, recv_persistent_request, &
36 send_persistent_request1, send_persistent_request2, &
37 recv_persistent_request1, recv_persistent_request2
38
39 integer :: nrecv1,nsend1, nrecv2,nsend2
40
41
42 double precision, dimension(:), pointer :: &
43 dsendbuffer, drecvbuffer
44 integer, dimension(:), pointer :: &
45 isendbuffer, irecvbuffer
46 character, dimension(:), pointer :: &
47 csendbuffer, crecvbuffer
48
49 integer :: nrecv,nsend
50 integer, pointer, dimension(:) :: &
51 recvrequest, sendrequest, &
52 xrecv,recvproc, recvijk, recvtag, &
53 xsend,sendproc, sendijk, sendtag
54
55 integer :: communicator
56
57
58
59
60 interface sendrecv_begin
61 module procedure &
62 sendrecv_begin_1d, &
63 sendrecv_begin_1i, &
64 sendrecv_begin_1c
65 end interface
66
67 interface sendrecv_end
68 module procedure &
69 sendrecv_end_1d, &
70 sendrecv_end_1i, &
71 sendrecv_end_1c
72 end interface
73
74 interface send_recv
75 module procedure &
76 send_recv_1d, send_recv_2d, send_recv_3d, &
77 send_recv_1i, &
78 send_recv_1c
79 end interface
80
81
82 contains
83
84
85
86
87 subroutine ijk_of( ijkp, i,j,k )
88
89
90
91
92 integer, intent(in) :: ijkp
93 integer, intent(out) :: i,j,k
94
95
96
97 integer :: k1,k2, j1,j2, i1,i2, &
98 ijk, isize,jsize,ksize, gijk
99
100 character(len=32), parameter :: name = "ijk_of"
101 logical :: isok_k, isok_j, isok_i, is_same, isok
102
103
104 = ijkp
105
106 i1 = istart3_all(myPE)
107 i2 = iend3_all(myPE)
108 j1 = jstart3_all(myPE)
109 j2 = jend3_all(myPE)
110 k1 = kstart3_all(myPE)
111 k2 = kend3_all(myPE)
112
113 ksize = (k2-k1+1)
114 jsize = (j2-j1+1)
115 isize = (i2-i1+1)
116
117
118 if (mod(ijk,isize*jsize).ne.0) then
119 k = int( ijk/(isize*jsize) ) + k1
120 else
121 k = int( ijk/(isize*jsize) ) + k1 -1
122 endif
123 ijk = ijk - (k-k1)*(isize*jsize)
124
125 if (mod(ijk,isize).ne.0) then
126 j = int( ijk/isize ) + j1
127 else
128 j = int( ijk/isize ) + j1 - 1
129 endif
130 ijk = ijk - (j-j1)*isize
131
132 i = (ijk-1) + i1
133
134
135 = (i1 <= i) .and. (i <= i2)
136 isok_j = (j1 <= j) .and. (j <= j2)
137 isok_k = (k1 <= k) .and. (k <= k2)
138 gijk = 1 + (i-i1) + (j-j1)*(i2-i1+1) + &
139 (k-k1)*(j2-j1+1)*(i2-i1+1)
140 is_same = (gijk .eq. ijkp)
141 isok = isok_i .and. isok_j .and. isok_k .and. is_same
142 if (.not.isok) then
143 call write_debug( name, 'i,j,k ', i,j,k )
144 call write_debug( name, 'ijkp, gijk ', ijkp, gijk )
145 endif
146
147
148 return
149 end subroutine ijk_of
150
151
152
153
154
155 subroutine ijk_of_gl( ijkp, i,j,k )
156
157
158
159
160 integer, intent(in) :: ijkp
161 integer, intent(out) :: i,j,k
162
163
164
165 integer :: k1,k2, j1,j2, i1,i2, &
166 ijk, isize,jsize,ksize, gijk
167
168 character(len=32), parameter :: name = "ijk_of_gl"
169 logical :: isok_k, isok_j, isok_i, is_same, isok
170
171
172 = ijkp
173
174 k1 = minval( kstart3_all(:) )
175 k2 = maxval( kend3_all(:) )
176 j1 = minval( jstart3_all(:) )
177 j2 = maxval( jend3_all(:) )
178 i1 = minval( istart3_all(:) )
179 i2 = maxval( iend3_all(:) )
180
181 ksize = (k2-k1+1)
182 jsize = (j2-j1+1)
183 isize = (i2-i1+1)
184
185 if (mod(ijk,isize*jsize).ne.0) then
186 k = int( ijk/(isize*jsize) ) + k1
187 else
188 k = int( ijk/(isize*jsize) ) + k1 -1
189 endif
190 ijk = ijk - (k-k1)*(isize*jsize)
191
192 if (mod(ijk,isize).ne.0) then
193 j = int( ijk/isize ) + j1
194 else
195 j = int( ijk/isize ) + j1 - 1
196 endif
197 ijk = ijk - (j-j1)*isize
198
199 i = (ijk-1) + i1
200
201
202 = (i1 <= i) .and. (i <= i2)
203 isok_j = (j1 <= j) .and. (j <= j2)
204 isok_k = (k1 <= k) .and. (k <= k2)
205 gijk = 1 + (i-i1) + (j-j1)*(i2-i1+1) + &
206 (k-k1)*(j2-j1+1)*(i2-i1+1)
207 is_same = (gijk .eq. ijkp)
208 isok = isok_i .and. isok_j .and. isok_k .and. is_same
209 if (.not.isok) then
210 call write_debug( name, 'i,j,k ', i,j,k )
211 call write_debug( name, 'ijkp, gijk ', ijkp, gijk )
212 endif
213
214 return
215 end subroutine ijk_of_gl
216
217
218
219
220
221
222
223 subroutine sendrecv_init(comm, &
224 cyclic_i,cyclic_j,cyclic_k, idebug )
225
226 implicit none
227
228
229
230 integer, intent(in) :: comm
231 logical,intent(in) :: cyclic_i,cyclic_j,cyclic_k
232 integer, intent(in), optional :: idebug
233
234
235
236 logical, parameter :: jfastest = .true.
237 integer, parameter :: message_tag_offset = 11
238
239 character(len=80), parameter :: name = 'sendrecv_init'
240
241 character(len=80), pointer, dimension(:) :: line
242 integer :: ip, lmax
243
244 integer :: layer,request, source, tag, datatype
245
246 integer :: lidebug
247 integer :: isize,jsize,ksize, ijksize
248 integer :: recvsize1, recvsize2, &
249 sendsize1, sendsize2
250
251 integer :: iter, i,j,k, ii, jj,kk, &
252 ntotal, icount,ipos, &
253 ilayer, i1,i2, j1,j2, k1,k2, &
254 ijk, ijk2, iproc, jproc, src,dest, &
255 ierror
256
257 logical :: isok, isvalid, ismine, is_halobc
258
259 integer, dimension(:,:,:), pointer :: ijk2proc
260 integer, pointer, dimension(:) :: &
261 istartx,iendx, jstartx,jendx, kstartx,kendx, &
262 ncount, &
263 recvproc, recvtag, xrecv, recvijk, &
264 sendproc, sendtag, xsend, sendijk
265
266
267
268 integer :: message_tag
269
270
271 (src,dest) = message_tag_offset + (1+src + dest*numPEs)
272
273 nullify( &
274 recvproc1, recvtag1, xrecv1, recvijk1, &
275 sendproc1, sendtag1, xsend1, sendijk1, &
276 recvproc2, recvtag2, xrecv2, recvijk2, &
277 sendproc2, sendtag2, xsend2, sendijk2)
278
279 nullify( &
280 send_persistent_request, recv_persistent_request, &
281 send_persistent_request1, send_persistent_request2, &
282 recv_persistent_request1, recv_persistent_request2 )
283
284 nullify( dsendbuffer, drecvbuffer )
285 nullify( isendbuffer, irecvbuffer )
286 nullify( csendbuffer, crecvbuffer )
287
288 nullify( &
289 recvrequest, sendrequest, &
290 xrecv,recvproc, recvijk, recvtag, &
291 xsend,sendproc, sendijk, sendtag )
292
293
294 = 0
295 if (present(idebug)) then
296 lidebug = idebug
297 endif
298
299 communicator = comm
300 call MPI_COMM_SIZE( comm, numPEs, ierror )
301 call MPI_Check( 'sendrecv_init:MPI_COMM_SIZE ', ierror )
302
303 call MPI_COMM_RANK( comm, myPE, ierror )
304 call MPI_Check( 'sendrecv_init:MPI_COMM_RANK ', ierror )
305
306
307
308 call assert( kmin1 .eq. minval( kstart1_all(:) ), &
309 '** sendrecv_init: invalid kmin1, ' // &
310 ' kmin1, minval(kstart1_all(:)) ', &
311 kmin1, minval(kstart1_all(:)) )
312
313 call assert( kmin2 .eq. minval( kstart2_all(:) ), &
314 '** sendrecv_init: invalid kmin2, ' // &
315 ' kmin2, minval(kstart2_all(:)) ', &
316 kmin2, minval(kstart2_all(:)) )
317
318 call assert( kmin3 .eq. minval( kstart3_all(:) ), &
319 '** sendrecv_init: invalid kmin3, ' // &
320 ' kmin3, minval(kstart3_all(:)) ', &
321 kmin3, minval(kstart3_all(:)) )
322
323 call assert( kmax1 .eq. maxval( kend1_all(:) ), &
324 '** sendrecv_init: invalid kmax1, ' // &
325 ' kmax1, maxval(kend1_all(:)) ', &
326 kmax1, maxval(kend1_all(:)) )
327
328 call assert( kmax2 .eq. maxval( kend2_all(:) ), &
329 '** sendrecv_init: invalid kmax2, ' // &
330 ' kmax2, maxval(kend2_all(:)) ', &
331 kmax2, maxval(kend2_all(:)) )
332
333 call assert( kmax3 .eq. maxval( kend3_all(:) ), &
334 '** sendrecv_init: invalid kmax3, ' // &
335 ' kmax3, maxval(kend3_all(:)) ', &
336 kmax3, maxval(kend3_all(:)) )
337
338
339 call assert( jmin1 .eq. minval( jstart1_all(:) ), &
340 '** sendrecv_init: invalid jmin1, ' // &
341 ' jmin1, minval(jstart1_all(:)) ', &
342 jmin1, minval(jstart1_all(:)) )
343
344 call assert( jmin2 .eq. minval( jstart2_all(:) ), &
345 '** sendrecv_init: invalid jmin2, ' // &
346 ' jmin2, minval(jstart2_all(:)) ', &
347 jmin2, minval(jstart2_all(:)) )
348
349 call assert( jmin3 .eq. minval( jstart3_all(:) ), &
350 '** sendrecv_init: invalid jmin3, ' // &
351 ' jmin3, minval(jstart3_all(:)) ', &
352 jmin3, minval(jstart3_all(:)) )
353
354 call assert( jmax1 .eq. maxval( jend1_all(:) ), &
355 '** sendrecv_init: invalid jmax1, ' // &
356 ' jmax1, maxval(jend1_all(:)) ', &
357 jmax1, maxval(jend1_all(:)) )
358
359 call assert( jmax2 .eq. maxval( jend2_all(:) ), &
360 '** sendrecv_init: invalid jmax2, ' // &
361 ' jmax2, maxval(jend2_all(:)) ', &
362 jmax2, maxval(jend2_all(:)) )
363
364 call assert( jmax3 .eq. maxval( jend3_all(:) ), &
365 '** sendrecv_init: invalid jmax3, ' // &
366 ' jmax3, maxval(jend3_all(:)) ', &
367 jmax3, maxval(jend3_all(:)) )
368
369
370 call assert( imin1 .eq. minval( istart1_all(:) ), &
371 '** sendrecv_init: invalid imin1, ' // &
372 ' imin1, minval(istart1_all(:)) ', &
373 imin1, minval(istart1_all(:)) )
374
375 call assert( imin2 .eq. minval( istart2_all(:) ), &
376 '** sendrecv_init: invalid imin2, ' // &
377 ' imin2, minval(istart2_all(:)) ', &
378 imin2, minval(istart2_all(:)) )
379
380 call assert( imin3 .eq. minval( istart3_all(:) ), &
381 '** sendrecv_init: invalid imin3, ' // &
382 ' imin3, minval(istart3_all(:)) ', &
383 imin3, minval(istart3_all(:)) )
384
385 call assert( imax1 .eq. maxval( iend1_all(:) ), &
386 '** sendrecv_init: invalid imax1, ' // &
387 ' imax1, maxval(iend1_all(:)) ', &
388 imax1, maxval(iend1_all(:)) )
389
390 call assert( imax2 .eq. maxval( iend2_all(:) ), &
391 '** sendrecv_init: invalid imax2, ' // &
392 ' imax2, maxval(iend2_all(:)) ', &
393 imax2, maxval(iend2_all(:)) )
394
395 call assert( imax3 .eq. maxval( iend3_all(:) ), &
396 '** sendrecv_init: invalid imax3, ' // &
397 ' imax3, maxval(iend3_all(:)) ', &
398 imax3, maxval(iend3_all(:)) )
399
400
401
402 call assert( jmin1 .le. jmax1, &
403 '** sendrecv_init: jmin1,jmax1 ', jmin1,jmax1 )
404 call assert( jmin2 .le. jmax2, &
405 '** sendrecv_init: jmin2,jmax2 ', jmin2,jmax2 )
406 call assert( jmin3 .le. jmax3, &
407 '** sendrecv_init: jmin3,jmax3 ', jmin3,jmax3 )
408
409 call assert( kmin1 .le. kmax1, &
410 '** sendrecv_init: kmin1,kmax1 ', kmin1,kmax1 )
411 call assert( kmin2 .le. kmax2, &
412 '** sendrecv_init: kmin2,kmax2 ', kmin2,kmax2 )
413 call assert( kmin3 .le. kmax3, &
414 '** sendrecv_init: kmin3,kmax3 ', kmin3,kmax3 )
415
416 call assert( imin1 .le. imax1, &
417 '** sendrecv_init: imin1,imax1 ', imin1,imax1 )
418 call assert( imin2 .le. imax2, &
419 '** sendrecv_init: imin2,imax2 ', imin2,imax2 )
420 call assert( imin3 .le. imax3, &
421 '** sendrecv_init: imin3,imax3 ', imin3,imax3 )
422
423
424
425 k1 = min( kmin1, min(kmin2, kmin3) )
426 k2 = max( kmax1, max(kmax2, kmax3) )
427 j1 = min( jmin1, min(jmin2, jmin3) )
428 j2 = max( jmax1, max(jmax2, jmax3) )
429 i1 = min( imin1, min(imin2, imin3) )
430 i2 = max( imax1, max(imax2, imax3) )
431
432 allocate( ijk2proc( i1:i2, j1:j2, k1:k2 ) )
433
434 if(localfunc) then
435
436 do k=kstart3_all(myPE),kend3_all(myPE)
437 do j=jstart3_all(myPE),jend3_all(myPE)
438 do i=istart3_all(myPE),iend3_all(myPE)
439 ijk = funijk(i,j,k)
440 call ijk_of(ijk, ii,jj,kk)
441 ijk2 = funijk( ii,jj,kk)
442
443 isvalid = (ii.eq.i).and.(jj.eq.j).and.(kk.eq.k).and.(ijk.eq.ijk2)
444 if (.not.isvalid) then
445 call write_debug( name, 'error with ijk_of ')
446
447 call write_debug( name, &
448 'istart3_all(myPE),iend3_all(myPE) ', &
449 istart3_all(myPE),iend3_all(myPE) )
450 call write_debug( name, &
451 'jstart3_all(myPE),jend3_all(myPE) ', &
452 jstart3_all(myPE),jend3_all(myPE) )
453 call write_debug( name, &
454 'kstart3_all(myPE),kend3_all(myPE) ', &
455 kstart3_all(myPE),kend3_all(myPE) )
456
457 call write_debug( name, 'i,j,k, ijk ', i,j,k, ijk )
458 call write_debug( name, 'ii,jj,kk, ijk2 ',&
459 ii,jj,kk,ijk2 )
460
461 endif
462 enddo
463 enddo
464 enddo
465 endif
466
467 if (lidebug.ge.1) then
468 call write_debug( name, 'imap ', imap )
469 call write_debug( name, 'jmap ', jmap )
470 call write_debug( name, 'kmap ', kmap )
471 endif
472
473
474
475
476
477
478
479
480
481 ( :,:,: ) = 0
482
483
484
485
486
487 do iproc=0,numPEs-1
488 i1 = istart1_all(iproc)
489 i2 = iend1_all(iproc)
490 j1 = jstart1_all(iproc)
491 j2 = jend1_all(iproc)
492 k1 = kstart1_all(iproc)
493 k2 = kend1_all(iproc)
494 if(istart3_all(iproc).eq.imin3) i1 = istart3_all(iproc)
495 if(iend3_all(iproc).eq.imax3) i2 = iend3_all(iproc)
496 if(jstart3_all(iproc).eq.jmin3) j1 = jstart3_all(iproc)
497 if(jend3_all(iproc).eq.jmax3) j2 = jend3_all(iproc)
498 if(kstart3_all(iproc).eq.kmin3) k1 = kstart3_all(iproc)
499 if(kend3_all(iproc).eq.kmax3) k2 = kend3_all(iproc)
500 do k=k1,k2
501 do j=j1,j2
502 do i=i1,i2
503 ijk2proc(i,j,k) = ijk2proc(i,j,k) + 1
504 enddo
505 enddo
506 enddo
507 enddo
508
509 do k=lbound(ijk2proc,3),ubound(ijk2proc,3)
510 do j=lbound(ijk2proc,2),ubound(ijk2proc,2)
511 do i=lbound(ijk2proc,1),ubound(ijk2proc,1)
512 isvalid = (ijk2proc(i,j,k) .eq. 1)
513 if (.not.isvalid) then
514 call write_debug(name, ' invalid decomposition ')
515 call write_debug(name, 'i,j,k ',i,j,k )
516 call write_debug(name, 'ijk2proc(i,j,k) ', &
517 ijk2proc(i,j,k))
518 call mfix_exit( myPE )
519 endif
520 enddo
521 enddo
522 enddo
523
524 ijk2proc(:,:,:) = -1
525 do iproc=0,numPEs-1
526 i1 = istart1_all(iproc)
527 i2 = iend1_all(iproc)
528 j1 = jstart1_all(iproc)
529 j2 = jend1_all(iproc)
530 k1 = kstart1_all(iproc)
531 k2 = kend1_all(iproc)
532 if(istart3_all(iproc).eq.imin3) i1 = istart3_all(iproc)
533 if(iend3_all(iproc).eq.imax3) i2 = iend3_all(iproc)
534 if(jstart3_all(iproc).eq.jmin3) j1 = jstart3_all(iproc)
535 if(jend3_all(iproc).eq.jmax3) j2 = jend3_all(iproc)
536 if(kstart3_all(iproc).eq.kmin3) k1 = kstart3_all(iproc)
537 if(kend3_all(iproc).eq.kmax3) k2 = kend3_all(iproc)
538 do k=k1,k2
539 do j=j1,j2
540 do i=i1,i2
541 ijk2proc(i,j,k) = iproc
542 enddo
543 enddo
544 enddo
545 enddo
546
547
548 allocate( ncount(0:numPEs-1) )
549 allocate( istartx(0:numPEs-1) )
550 allocate( jstartx(0:numPEs-1) )
551 allocate( kstartx(0:numPEs-1) )
552 allocate( iendx(0:numPEs-1) )
553 allocate( jendx(0:numPEs-1) )
554 allocate( kendx(0:numPEs-1) )
555
556
557 do ilayer=1,2
558 if (ilayer.eq.1) then
559 kstartx(:) = kstart2_all(:)
560 kendx(:) = kend2_all(:)
561 jstartx(:) = jstart2_all(:)
562 jendx(:) = jend2_all(:)
563 istartx(:) = istart2_all(:)
564 iendx(:) = iend2_all(:)
565 else
566 kstartx(:) = kstart3_all(:)
567 kendx(:) = kend3_all(:)
568 jstartx(:) = jstart3_all(:)
569 jendx(:) = jend3_all(:)
570 istartx(:) = istart3_all(:)
571 iendx(:) = iend3_all(:)
572 endif
573
574 if (lidebug.ge.1) then
575 call write_debug(name, 'determine send schedule ', myPE )
576 endif
577
578
579
580
581
582
583
584
585 (:) = 0
586 do iproc=0,numPEs-1
587 if (iproc.ne.myPE) then
588 k1 = lbound(ijk2proc,3)
589 k2 = ubound(ijk2proc,3)
590 j1 = lbound(ijk2proc,2)
591 j2 = ubound(ijk2proc,2)
592 i1 = lbound(ijk2proc,1)
593 i2 = ubound(ijk2proc,1)
594
595 do k=kstartx(iproc),kendx(iproc)
596 do j=jstartx(iproc),jendx(iproc)
597 do i=istartx(iproc),iendx(iproc)
598 ii = imap(i)
599 jj = jmap(j)
600 kk = kmap(k)
601
602 isvalid = (k1.le.kk).and.(kk.le.k2)
603 call assert( isvalid, '** sendrecv_init: invalid kk ', kk )
604 isvalid = (j1.le.jj).and.(jj.le.j2)
605 call assert( isvalid, '** sendrecv_init: invalid jj ', jj )
606 isvalid = (i1.le.ii).and.(ii.le.i2)
607 call assert( isvalid, '** sendrecv_init: invalid ii ', ii )
608 jproc = ijk2proc( ii,jj,kk )
609
610 ismine = (jproc .eq. myPE)
611 if (ismine) then
612 ncount(iproc) = ncount(iproc) + 1
613 endif
614 enddo
615 enddo
616 enddo
617 endif
618 enddo
619
620
621
622
623 = 0
624 nsend = 0
625 do iproc=0,numPEs-1
626 ntotal = ntotal + ncount(iproc)
627 if (ncount(iproc).ge.1) then
628 nsend = nsend + 1
629 endif
630 enddo
631
632 if (lidebug.ge.1) then
633 call write_debug( name, 'ncount = ', ncount )
634 call write_debug( name, 'nsend, ntotal ', nsend, ntotal )
635 endif
636
637 allocate( xsend(nsend+1) )
638 allocate( sendijk( max(1,ntotal) ) )
639 allocate( sendproc(max(1,nsend)) )
640
641 nsend = 0
642 do iproc=0,numPEs-1
643 if (ncount(iproc).ne.0) then
644 nsend = nsend + 1
645 sendproc(nsend) = iproc
646 endif
647 enddo
648
649 xsend(1) = 1
650 do i=1,nsend
651 iproc = sendproc(i)
652 xsend(i+1) = xsend(i) + ncount(iproc)
653 enddo
654
655 allocate( sendtag( max(1,nsend) ) )
656 do ii=1,nsend
657 iproc = sendproc(ii)
658 src = myPE
659 dest = iproc
660 sendtag(ii) = message_tag( src, dest )
661 enddo
662
663
664
665
666
667 = 1
668 do iter=1,nsend
669 iproc = sendproc(iter)
670 icount = 0
671 do k=kstartx(iproc),kendx(iproc)
672
673 if (jfastest) then
674 do i=istartx(iproc),iendx(iproc)
675 do j=jstartx(iproc),jendx(iproc)
676 ii = imap(i)
677 jj = jmap(j)
678 kk = kmap(k)
679 jproc = ijk2proc(ii,jj,kk)
680 ismine = (jproc.eq.myPE)
681 if (ismine) then
682 icount = icount + 1
683 ijk = funijk(ii,jj,kk)
684 ipos = xsend(iter)-1 + icount
685 sendijk( ipos ) = ijk
686 endif
687 enddo
688 enddo
689 else
690 do j=jstartx(iproc),jendx(iproc)
691 do i=istartx(iproc),iendx(iproc)
692 ii = imap(i)
693 jj = jmap(j)
694 kk = kmap(k)
695 jproc = ijk2proc(ii,jj,kk)
696 ismine = (jproc.eq.myPE)
697 if (ismine) then
698 icount = icount + 1
699 ijk = funijk(ii,jj,kk)
700 ipos = xsend(iter)-1 + icount
701 sendijk( ipos ) = ijk
702 endif
703 enddo
704 enddo
705 endif
706 enddo
707 = (icount .eq. ncount(iproc))
708 call assert( isvalid, &
709 '** sendrecv_init: icount != ncount(iproc) ', iproc)
710 enddo
711
712 if (lidebug.ge.1) then
713 call write_debug(name, 'determine recv schedule ', myPE )
714 endif
715
716
717
718
719
720
721
722
723
724
725
726
727
728 (:) = 0
729 k1 = lbound(ijk2proc,3)
730 k2 = ubound(ijk2proc,3)
731 j1 = lbound(ijk2proc,2)
732 j2 = ubound(ijk2proc,2)
733 i1 = lbound(ijk2proc,1)
734 i2 = ubound(ijk2proc,1)
735
736 do k=kstartx(myPE),kendx(myPE)
737 do j=jstartx(myPE),jendx(myPE)
738 do i=istartx(myPE),iendx(myPE)
739 ii = imap(i)
740 jj = jmap(j)
741 kk = kmap(k)
742
743 isvalid = (k1.le.kk).and.(kk.le.k2)
744 call assert( isvalid, '** sendrecv_init: invalid kk ', kk )
745
746 isvalid = (j1.le.jj).and.(jj.le.j2)
747 call assert( isvalid, '** sendrecv_init: invalid jj ', jj )
748
749 isvalid = (i1.le.ii).and.(ii.le.i2)
750 call assert( isvalid, '** sendrecv_init: invalid ii ', ii )
751
752
753 iproc = ijk2proc(ii,jj,kk)
754 is_halobc = (iproc.eq.-1)
755 ismine = (iproc.eq.myPE)
756 if (.not.ismine) then
757 isvalid = (0 .le. iproc) .and. &
758 (iproc.le.numPEs-1) .and. &
759 (iproc.ne.myPE)
760 call assert( isvalid, &
761 '** sendrecv_init: invalid iproc ',iproc)
762
763 ncount(iproc) = ncount(iproc) + 1
764 endif
765 enddo
766 enddo
767 enddo
768
769 ncount(myPE) = 0
770
771 ntotal = 0
772 do iproc=0,numPEs-1
773 ntotal = ntotal + ncount(iproc)
774 enddo
775
776 nrecv = count( ncount(:) .ne. 0)
777
778 allocate( recvproc( max(1,nrecv) ) )
779
780 nrecv = 0
781 do iproc=0,numPEs-1
782 if (ncount(iproc).ne.0) then
783 nrecv = nrecv + 1
784 recvproc(nrecv) = iproc
785 endif
786 enddo
787
788 allocate( xrecv(nrecv+1) )
789 allocate( recvijk(max(1,ntotal)) )
790
791 xrecv(1) = 1
792 do iter=1,nrecv
793 iproc = recvproc(iter)
794 xrecv(iter+1) = xrecv(iter) + ncount(iproc)
795 enddo
796
797 allocate( recvtag( max(1,nrecv) ) )
798
799 do iter=1,nrecv
800 iproc = recvproc(iter)
801 src = iproc
802 dest = myPE
803 recvtag(iter) = message_tag( src, dest )
804 enddo
805
806
807
808
809
810 if (lidebug.ge.1) then
811 call write_debug( name, 'recv second pass ', myPE )
812 endif
813
814 ipos = 1
815
816 do iter=1,nrecv
817 jproc = recvproc(iter)
818 do k=kstartx(myPE),kendx(myPE)
819
820 if (jfastest) then
821 do i=istartx(myPE),iendx(myPE)
822 do j=jstartx(myPE),jendx(myPE)
823 ii = imap(i)
824 jj = jmap(j)
825 kk = kmap(k)
826
827 iproc = ijk2proc(ii,jj,kk)
828 is_halobc = (iproc.eq.-1)
829 ismine = (iproc.eq.myPE)
830
831 if ((.not.ismine) .and. (iproc.eq.jproc)) then
832 ijk = funijk( i,j,k)
833 recvijk( ipos ) = ijk
834 ipos = ipos + 1
835 endif
836 enddo
837 enddo
838
839 else
840 do j=jstartx(myPE),jendx(myPE)
841 do i=istartx(myPE),iendx(myPE)
842 ii = imap(i)
843 jj = jmap(j)
844 kk = kmap(k)
845
846 iproc = ijk2proc(ii,jj,kk)
847 is_halobc = (iproc.eq.-1)
848 ismine = (iproc.eq.myPE)
849
850 if ((.not.ismine) .and. (iproc.eq.jproc)) then
851 ijk = funijk( i,j,k)
852 recvijk( ipos ) = ijk
853 ipos = ipos + 1
854 endif
855 enddo
856 enddo
857 endif
858 enddo
859 enddo
860
861
862 if (ilayer.eq.1) then
863 nsend1 = nsend
864 xsend1 => xsend
865 sendijk1 => sendijk
866 sendproc1 => sendproc
867 sendtag1 => sendtag
868
869 nrecv1 = nrecv
870 xrecv1 => xrecv
871 recvijk1 => recvijk
872 recvproc1 => recvproc
873 recvtag1 => recvtag
874 else
875 nsend2 = nsend
876 xsend2 => xsend
877 sendijk2 => sendijk
878 sendproc2 => sendproc
879 sendtag2 => sendtag
880
881 nrecv2 = nrecv
882 xrecv2 => xrecv
883 recvijk2 => recvijk
884 recvproc2 => recvproc
885 recvtag2 => recvtag
886 endif
887
888
889 nullify( xsend )
890 nullify( sendijk )
891 nullify( sendproc )
892 nullify( sendtag )
893 nullify( xrecv )
894 nullify( recvijk )
895 nullify( recvproc )
896 nullify( recvtag )
897
898 enddo
899
900
901 deallocate( ncount )
902 deallocate( ijk2proc )
903
904 deallocate( istartx )
905 deallocate( jstartx )
906 deallocate( kstartx )
907 deallocate( iendx )
908 deallocate( jendx )
909 deallocate( kendx )
910
911 nullify( ncount )
912 nullify( ijk2proc )
913
914 nullify( istartx )
915 nullify( jstartx )
916 nullify( kstartx )
917 nullify( iendx )
918 nullify( jendx )
919 nullify( kendx )
920
921
922
923 if (lidebug.ge.1) then
924 call write_debug( name, ' allocate message buffers ' )
925 call write_debug( name, 'nrecv1 ', nrecv1 )
926 call write_debug( name, 'recvproc1 ', recvproc1 )
927 call write_debug( name, 'recvtag1 ', recvtag1 )
928 call write_debug( name, 'xrecv1 ', xrecv1 )
929
930 lmax = size(recvijk1)
931 allocate( line(lmax) )
932 line(:) = " "
933 ip = 1
934 do ii=lbound(recvijk1,1),ubound(recvijk1,1)
935 ijk = recvijk1(ii)
936 if(localfunc) then
937 call ijk_of(ijk,i,j,k)
938 else
939 i = i_of(ijk)
940 j = j_of(ijk)
941 k = k_of(ijk)
942 endif
943
944 write(line(ip),9001) ii,ijk, i,j,k
945 9001 format('recvijk1( ', i6,') = ', &
946 i6, '( ', i6,',',i6,',',i6,') ')
947
948 ip = ip + 1
949 enddo
950 call write_error( name, line, lmax )
951 deallocate( line )
952 nullify( line )
953
954
955 lmax = size(recvijk2)
956 allocate( line(lmax) )
957 line(:) = " "
958 ip = 1
959 do ii=lbound(recvijk2,1),ubound(recvijk2,1)
960 ijk = recvijk2(ii)
961 if(localfunc) then
962 call ijk_of(ijk,i,j,k)
963 else
964 i = i_of(ijk)
965 j = j_of(ijk)
966 k = k_of(ijk)
967 endif
968
969 write(line(ip),9101) ii,ijk, i,j,k
970 9101 format('recvijk2( ', i6,') = ', &
971 i6, '( ', i6,',',i6,',',i6,') ')
972
973 ip = ip + 1
974 enddo
975 call write_error( name, line, lmax )
976 deallocate( line )
977 nullify( line )
978
979 call write_debug( name, ' allocate message buffers ' )
980 call write_debug( name, 'nsend1 ', nsend1 )
981 call write_debug( name, 'sendproc1 ', sendproc1 )
982 call write_debug( name, 'sendtag1 ', sendtag1 )
983 call write_debug( name, 'xsend1 ', xsend1 )
984
985 lmax = size(sendijk1)
986 allocate(line(lmax))
987 line(:) = " "
988 ip = 1
989 do ii=lbound(sendijk1,1),ubound(sendijk1,1)
990 ijk = sendijk1(ii)
991 if(localfunc) then
992 call ijk_of(ijk,i,j,k)
993 else
994 i = i_of(ijk)
995 j = j_of(ijk)
996 k = k_of(ijk)
997 endif
998
999 write(line(ip),9002) ii,ijk, i,j,k
1000 9002 format('sendijk1( ', i6,') = ', &
1001 i6, '( ', i6,',',i6,',',i6,') ')
1002
1003 ip = ip + 1
1004 enddo
1005 call write_error( name, line, lmax )
1006 deallocate( line )
1007 nullify( line )
1008
1009 lmax = size(sendijk2)
1010 allocate(line(lmax))
1011 line(:) = " "
1012 ip = 1
1013 do ii=lbound(sendijk2,1),ubound(sendijk2,1)
1014 ijk = sendijk2(ii)
1015 if(localfunc) then
1016 call ijk_of(ijk,i,j,k)
1017 else
1018 i = i_of(ijk)
1019 j = j_of(ijk)
1020 k = k_of(ijk)
1021 endif
1022
1023 write(line(ip),9102) ii,ijk, i,j,k
1024 9102 format('sendijk2( ', i6,') = ', &
1025 i6, '( ', i6,',',i6,',',i6,') ')
1026
1027 ip = ip + 1
1028 enddo
1029 call write_error( name, line, lmax )
1030 deallocate( line )
1031 nullify( line )
1032 endif
1033
1034
1035
1036
1037
1038 = max(1,max(nsend1,nsend2))
1039 allocate( sendrequest( isize ) )
1040 allocate( send_persistent_request1( isize ) )
1041 allocate( send_persistent_request2( isize ) )
1042
1043 isize = max(1,max(nrecv1,nrecv2))
1044 allocate( recvrequest( isize ) )
1045 allocate( recv_persistent_request1( isize ) )
1046 allocate( recv_persistent_request2( isize ) )
1047
1048
1049
1050 = xrecv1( nrecv1+1)-1
1051 recvsize2 = xrecv2( nrecv2+1)-1
1052
1053 isize = max(1,max(recvsize1,recvsize2))
1054 allocate( drecvbuffer( isize ) )
1055
1056 sendsize1 = xsend1( nsend1+1)-1
1057 sendsize2 = xsend2( nsend2+1)-1
1058
1059 isize = max(1,max(sendsize1,sendsize2))
1060 allocate( dsendbuffer( isize ) )
1061
1062
1063 if (use_persistent_message) then
1064 datatype = MPI_DOUBLE_PRECISION
1065
1066 do layer=1,2
1067 if (layer.eq.1) then
1068 nrecv = nrecv1
1069 recvtag =>recvtag1
1070 recvproc => recvproc1
1071 recvijk => recvijk1
1072 xrecv => xrecv1
1073
1074 nsend = nsend1
1075 sendtag => sendtag1
1076 sendproc => sendproc1
1077 sendijk => sendijk1
1078 xsend => xsend1
1079
1080 send_persistent_request => send_persistent_request1
1081 recv_persistent_request => recv_persistent_request1
1082
1083 else
1084 nrecv = nrecv2
1085 recvtag =>recvtag2
1086 recvproc => recvproc2
1087 recvijk => recvijk2
1088 xrecv => xrecv2
1089
1090 nsend = nsend2
1091 sendtag => sendtag2
1092 sendproc => sendproc2
1093 sendijk => sendijk2
1094 xsend => xsend2
1095
1096 send_persistent_request => send_persistent_request2
1097 recv_persistent_request => recv_persistent_request2
1098 endif
1099
1100 do ii=1,nrecv
1101 j1 = xrecv(ii)
1102 j2 = xrecv(ii+1)-1
1103 icount = j2-j1+1
1104 source = recvproc( ii )
1105 tag = recvtag( ii )
1106
1107 if (lidebug.ge.2) then
1108 call write_debug(name, 'mpi_recv_init: ii,j1,j2 ', &
1109 ii,j1,j2 )
1110 call write_debug(name, 'icount, source, tag ', &
1111 icount,source,tag )
1112 endif
1113
1114 call MPI_RECV_INIT( drecvbuffer(j1), icount, datatype, &
1115 source, tag, comm, request, ierror )
1116 call MPI_Check( 'sendrecv_begin_1d:MPI_IRECV ', ierror )
1117
1118 recv_persistent_request(ii) = request
1119 enddo
1120
1121 do ii=1,nsend
1122 j1 = xsend(ii)
1123 j2 = xsend(ii+1)-1
1124 dest = sendproc( ii )
1125 tag = sendtag( ii )
1126 icount = j2-j1+1
1127
1128 if (lidebug.ge.2) then
1129 call write_debug(name, 'mpi_send_init: ii,j1,j2 ', &
1130 ii,j1,j2)
1131 call write_debug(name, 'icount, dest, tag ', &
1132 icount,dest,tag )
1133 endif
1134
1135 call MPI_SEND_INIT( dsendbuffer(j1), icount, datatype, &
1136 dest, tag, comm, request, ierror )
1137 call MPI_Check( 'sendrecv_begin_1d:MPI_SEND_INIT ', &
1138 ierror)
1139
1140 send_persistent_request( ii ) = request
1141 enddo
1142 enddo
1143
1144 endif
1145
1146 if (lidebug.ge.1) then
1147 call write_debug(name, ' end of sendrecv_init ', myPE )
1148 endif
1149
1150 end subroutine sendrecv_init
1151
1152
1153
1154
1155
1156
1157 subroutine sendrecv_begin_1d( XX, ilayer, idebug )
1158
1159 implicit none
1160
1161
1162
1163 integer, intent(in),optional :: ilayer
1164 double precision, intent(inout), dimension(:) :: XX
1165 integer, intent(in), optional :: idebug
1166
1167
1168
1169 character(len=80), parameter :: name = 'sendrecv_begin_1d'
1170 integer :: lidebug
1171 integer :: layer, datatype, comm, recvsize, sendsize, &
1172 request, count, source,dest, tag, ierror
1173 integer :: ijk, jj, j1, j2, ii
1174
1175
1176
1177
1178
1179
1180
1181
1182
1183
1184
1185
1186
1187
1188
1189
1190
1191
1192 = 0
1193
1194 if (present(idebug)) then
1195 lidebug = idebug
1196 endif
1197
1198 layer = 1
1199 if (present(ilayer)) then
1200 layer = ilayer
1201 endif
1202
1203 if (layer.eq.1) then
1204 nrecv = nrecv1
1205 recvtag =>recvtag1
1206 recvproc => recvproc1
1207 recvijk => recvijk1
1208 xrecv => xrecv1
1209
1210 nsend = nsend1
1211 sendtag => sendtag1
1212 sendproc => sendproc1
1213 sendijk => sendijk1
1214 xsend => xsend1
1215
1216 send_persistent_request => send_persistent_request1
1217 recv_persistent_request => recv_persistent_request1
1218
1219 else
1220 nrecv = nrecv2
1221 recvtag =>recvtag2
1222 recvproc => recvproc2
1223 recvijk => recvijk2
1224 xrecv => xrecv2
1225
1226 nsend = nsend2
1227 sendtag => sendtag2
1228 sendproc => sendproc2
1229 sendijk => sendijk2
1230 xsend => xsend2
1231
1232 send_persistent_request => send_persistent_request2
1233 recv_persistent_request => recv_persistent_request2
1234
1235 endif
1236
1237
1238
1239
1240 if (lidebug.ge.1) then
1241 call write_debug(name, 'post asynchronous receives, nrecv = ', &
1242 nrecv)
1243 endif
1244
1245 if (nrecv.ge.1) then
1246 recvsize = xrecv( nrecv+1)-1
1247
1248 if (lidebug.ge.1) then
1249 call write_debug( name, 'recvsize, ubound(drecvbuffer,1) ',&
1250 recvsize, ubound(drecvbuffer,1) )
1251
1252 call write_debug( name, 'ubound(xrecv,1) ', &
1253 ubound(xrecv,1) )
1254 call write_debug( name, 'ubound(recvproc,1) ', &
1255 ubound(recvproc,1) )
1256 call write_debug( name, 'ubound(recvtag,1) ', &
1257 ubound(recvtag,1) )
1258 endif
1259
1260
1261 = MPI_DOUBLE_PRECISION
1262 comm = communicator
1263
1264 if (use_persistent_message) then
1265
1266 if (lidebug.ge.2) then
1267 call write_debug( name,'before startall for recv ',&
1268 recv_persistent_request)
1269 endif
1270
1271 call MPI_STARTALL( nrecv, recv_persistent_request, ierror )
1272
1273 if (lidebug.ge.2) then
1274 call write_debug( name,'after startall for recv, ierror',&
1275 ierror)
1276 endif
1277
1278 call MPI_Check( 'sendrecv_begin: MPI_STARTALL ', ierror )
1279
1280 else
1281
1282 do ii=1,nrecv
1283 j1 = xrecv(ii)
1284 j2 = xrecv(ii+1)-1
1285 count = j2-j1+1
1286 source = recvproc( ii )
1287 tag = recvtag( ii )
1288
1289 if (lidebug.ge.2) then
1290 call write_debug(name, 'mpi_irecv: ii,j1,j2 ', &
1291 ii, j1, j2)
1292 call write_debug(name, 'count, source, tag ', &
1293 count,source,tag )
1294 endif
1295
1296 call MPI_IRECV( drecvbuffer(j1), count, datatype, &
1297 source, tag, comm, request, ierror )
1298
1299 call MPI_Check( 'sendrecv_begin_1d:MPI_IRECV ', ierror )
1300
1301 recvrequest( ii ) = request
1302 enddo
1303 endif
1304 endif
1305
1306
1307
1308
1309 if (lidebug.ge.1) then
1310 call write_debug(name, 'post asynchronous sends ')
1311 endif
1312
1313 if (nsend.ge.1) then
1314 sendsize = xsend( nsend+1)-1
1315
1316 if (lidebug.ge.1) then
1317 call write_debug( name, &
1318 'sendsize, ubound(dsendbuffer,1) ', &
1319 sendsize, ubound(dsendbuffer,1) )
1320
1321 call write_debug( name, 'ubound(xsend,1) ', &
1322 ubound(xsend,1) )
1323 call write_debug( name, 'ubound(sendproc,1) ', &
1324 ubound(sendproc,1) )
1325 call write_debug( name, 'ubound(sendtag,1) ', &
1326 ubound(sendtag,1) )
1327 endif
1328
1329
1330 = MPI_DOUBLE_PRECISION
1331 comm = communicator
1332
1333 if (use_persistent_message) then
1334
1335
1336
1337 = xsend(1)
1338 j2 = xsend(nsend+1)-1
1339
1340 do jj=j1,j2
1341 ijk = sendijk( jj )
1342 dsendbuffer( jj ) = XX(ijk)
1343 enddo
1344
1345 if (lidebug.ge.2) then
1346 call write_debug(name,'before mpi_startall send ',&
1347 send_persistent_request )
1348 endif
1349
1350 call MPI_STARTALL( nsend, send_persistent_request, ierror )
1351
1352 if (lidebug .ge.2) then
1353 call write_debug(name,'after mpi_startall send ',&
1354 send_persistent_request )
1355 endif
1356
1357 call MPI_Check( 'sendrecv_begin_1d:MPI_STARTALL ', ierror )
1358
1359 else
1360
1361 do ii=1,nsend
1362
1363 = xsend(ii)
1364 j2 = xsend(ii+1)-1
1365 count = j2-j1+1
1366
1367 do jj=j1,j2
1368 ijk = sendijk( jj )
1369 dsendbuffer(jj) = XX(ijk)
1370 enddo
1371
1372 dest = sendproc( ii )
1373 tag = sendtag( ii )
1374
1375 if (lidebug.ge.2) then
1376 call write_debug(name, 'mpi_isend: ii,j1,j2 ', &
1377 ii,j1,j2)
1378 call write_debug(name, 'count, dest, tag ', &
1379 count,dest,tag)
1380 endif
1381
1382 call MPI_ISEND( dsendbuffer(j1), count, datatype, dest, &
1383 tag, comm, request, ierror )
1384 call MPI_Check( 'sendrecv_begin_1d:MPI_ISEND ', ierror )
1385
1386 sendrequest( ii ) = request
1387 enddo
1388 endif
1389 endif
1390
1391
1392 return
1393 end subroutine sendrecv_begin_1d
1394
1395
1396
1397
1398
1399
1400
1401 subroutine sendrecv_begin_1i( XX, ilayer, idebug )
1402
1403 implicit none
1404
1405
1406
1407 integer, intent(in),optional :: ilayer
1408 integer, intent(inout), dimension(:) :: XX
1409 integer, intent(in), optional :: idebug
1410
1411
1412
1413 character(len=80), parameter :: name = 'sendrecv_begin_1i'
1414 integer :: lidebug
1415 integer :: layer, datatype, comm, recvsize, sendsize, &
1416 request, count, source, dest, tag, ierror
1417 integer :: ijk, jj, j1, j2, ii
1418
1419
1420
1421
1422
1423
1424
1425
1426
1427
1428
1429
1430
1431
1432
1433
1434
1435
1436
1437 = 0
1438 if (present(idebug)) then
1439 lidebug = idebug
1440 endif
1441
1442 layer = 1
1443 if (present(ilayer)) then
1444 layer = ilayer
1445 endif
1446
1447 if (layer.eq.1) then
1448 nrecv = nrecv1
1449 recvtag =>recvtag1
1450 recvproc => recvproc1
1451 recvijk => recvijk1
1452 xrecv => xrecv1
1453
1454 nsend = nsend1
1455 sendtag => sendtag1
1456 sendproc => sendproc1
1457 sendijk => sendijk1
1458 xsend => xsend1
1459 else
1460 nrecv = nrecv2
1461 recvtag =>recvtag2
1462 recvproc => recvproc2
1463 recvijk => recvijk2
1464 xrecv => xrecv2
1465
1466 nsend = nsend2
1467 sendtag => sendtag2
1468 sendproc => sendproc2
1469 sendijk => sendijk2
1470 xsend => xsend2
1471 endif
1472
1473
1474
1475
1476
1477 if (lidebug.ge.1) then
1478 call write_debug(name, &
1479 'post asynchronous receives, nrecv = ', nrecv )
1480 endif
1481
1482 if (nrecv.ge.1) then
1483 recvsize = xrecv( nrecv+1)-1
1484 allocate( irecvbuffer( recvsize ) )
1485
1486 if (lidebug.ge.1) then
1487 call write_debug( name, &
1488 'recvsize, ubound(irecvbuffer,1) ', &
1489 recvsize, ubound(irecvbuffer,1) )
1490 call write_debug( name, 'ubound(xrecv,1) ', &
1491 ubound(xrecv,1) )
1492 call write_debug( name, 'ubound(recvproc,1) ', &
1493 ubound(recvproc,1) )
1494 call write_debug( name, 'ubound(recvtag,1) ', &
1495 ubound(recvtag,1) )
1496 endif
1497
1498
1499 = MPI_INTEGER
1500 comm = communicator
1501
1502 do ii=1,nrecv
1503 j1 = xrecv(ii)
1504 j2 = xrecv(ii+1)-1
1505 count = j2-j1+1
1506 source = recvproc( ii )
1507 tag = recvtag( ii )
1508
1509 if (lidebug.ge.2) then
1510 call write_debug(name, 'mpi_irecv: ii,j1,j2 ', ii,j1,j2 )
1511 call write_debug(name, 'count, source, tag ', &
1512 count,source,tag )
1513 endif
1514
1515 call MPI_IRECV( irecvbuffer(j1), count, datatype, &
1516 source, tag, comm, request, ierror )
1517 call MPI_Check( 'sendrecv_begin_1i:MPI_IRECV ', ierror )
1518
1519 recvrequest( ii ) = request
1520 enddo
1521 endif
1522
1523
1524
1525
1526
1527
1528 if (lidebug.ge.1) then
1529 call write_debug(name, 'post asynchronous sends ')
1530 endif
1531
1532 if (nsend.ge.1) then
1533 sendsize = xsend( nsend+1)-1
1534 allocate( isendbuffer( sendsize ) )
1535
1536 if (lidebug.ge.1) then
1537 call write_debug( name, 'sendsize, ubound(isendbuffer,1) ',&
1538 sendsize, ubound(isendbuffer,1) )
1539 call write_debug( name, 'ubound(xsend,1) ', &
1540 ubound(xsend,1) )
1541 call write_debug( name, 'ubound(sendproc,1) ', &
1542 ubound(sendproc,1) )
1543 call write_debug( name, 'ubound(sendtag,1) ', &
1544 ubound(sendtag,1) )
1545 endif
1546
1547
1548 = MPI_INTEGER
1549 comm = communicator
1550
1551 do ii=1,nsend
1552
1553 = xsend(ii)
1554 j2 = xsend(ii+1)-1
1555 count = j2-j1+1
1556
1557 do jj=j1,j2
1558 ijk = sendijk( jj )
1559 isendbuffer(jj) = XX(ijk)
1560 enddo
1561
1562 dest = sendproc( ii )
1563 tag = sendtag( ii )
1564
1565 if (lidebug.ge.2) then
1566 call write_debug(name, 'mpi_isend: ii,j1,j2 ', ii,j1,j2)
1567 call write_debug(name, 'count, dest, tag ', count, &
1568 dest, tag)
1569 endif
1570
1571 call MPI_ISEND( isendbuffer(j1), count, datatype, dest, &
1572 tag, comm, request, ierror )
1573 call MPI_Check( 'sendrecv_begin_1i:MPI_ISEND ', ierror )
1574
1575 sendrequest( ii ) = request
1576 enddo
1577 endif
1578
1579
1580 return
1581 end subroutine sendrecv_begin_1i
1582
1583
1584
1585
1586
1587
1588 subroutine sendrecv_begin_1c( XX, ilayer, idebug )
1589
1590 use functions
1591
1592 implicit none
1593
1594
1595
1596
1597 integer, intent(in),optional :: ilayer
1598 character(len=*), intent(inout), dimension(:) :: XX
1599 integer, intent(in), optional :: idebug
1600
1601
1602
1603 character(len=80), parameter :: name = 'sendrecv_begin_1c'
1604 integer :: lidebug
1605 integer :: layer, datatype, comm, recvsize, sendsize, &
1606 request, count, source, dest, tag, ierror
1607 integer :: ijk, jj, j1, j2, ii
1608 integer :: ic, clen, jpos
1609
1610
1611
1612
1613
1614
1615
1616
1617
1618
1619
1620
1621
1622
1623
1624
1625
1626
1627
1628 = 0
1629 if (present(idebug)) then
1630 lidebug = idebug
1631 endif
1632
1633 layer = 1
1634 if (present(ilayer)) then
1635 layer = ilayer
1636 endif
1637
1638 jpos = lbound(XX,1)
1639 clen = len( XX( jpos ) )
1640
1641 if (layer.eq.1) then
1642 nrecv = nrecv1
1643 recvtag =>recvtag1
1644 recvproc => recvproc1
1645 recvijk => recvijk1
1646 xrecv => xrecv1
1647
1648 nsend = nsend1
1649 sendtag => sendtag1
1650 sendproc => sendproc1
1651 sendijk => sendijk1
1652 xsend => xsend1
1653 else
1654 nrecv = nrecv2
1655 recvtag =>recvtag2
1656 recvproc => recvproc2
1657 recvijk => recvijk2
1658 xrecv => xrecv2
1659
1660 nsend = nsend2
1661 sendtag => sendtag2
1662 sendproc => sendproc2
1663 sendijk => sendijk2
1664 xsend => xsend2
1665 endif
1666
1667
1668
1669
1670 if (lidebug.ge.1) then
1671 call write_debug(name, 'post asynchronous receives, nrecv = ',&
1672 nrecv )
1673 endif
1674
1675 if (nrecv.ge.1) then
1676 recvsize = xrecv( nrecv+1)-1
1677
1678 allocate( crecvbuffer( recvsize*clen ) )
1679
1680 if (lidebug.ge.1) then
1681 call write_debug( name, 'recvsize, ubound(crecvbuffer,1) ', &
1682 recvsize, ubound(crecvbuffer,1) )
1683 call write_debug( name, 'ubound(xrecv,1) ', &
1684 ubound(xrecv,1) )
1685 call write_debug( name, 'ubound(recvproc,1) ', &
1686 ubound(recvproc,1) )
1687 call write_debug( name, 'ubound(recvtag,1) ', &
1688 ubound(recvtag,1) )
1689 endif
1690
1691
1692 = MPI_CHARACTER
1693 comm = communicator
1694
1695 do ii=1,nrecv
1696 j1 = xrecv(ii)
1697 j2 = xrecv(ii+1)-1
1698
1699 count = j2-j1+1
1700 count = count*clen
1701
1702 source = recvproc( ii )
1703 tag = recvtag( ii )
1704
1705 if (lidebug.ge.2) then
1706 call write_debug(name, 'mpi_irecv: ii,j1,j2 ', ii,j1,j2 )
1707 call write_debug(name, 'count, source, tag ', &
1708 count,source,tag )
1709 endif
1710
1711 jpos = 1 + (j1-1)*clen
1712 call MPI_IRECV( crecvbuffer(jpos), count, datatype, source, &
1713 tag, comm, request, ierror )
1714 call MPI_Check( 'sendrecv_begin_1c:MPI_IRECV ', ierror )
1715
1716 recvrequest( ii ) = request
1717 enddo
1718 endif
1719
1720
1721
1722
1723
1724 if (lidebug.ge.1) then
1725 call write_debug(name, 'post asynchronous sends ')
1726 endif
1727
1728 if (nsend.ge.1) then
1729 sendsize = xsend( nsend+1)-1
1730
1731 allocate( csendbuffer( sendsize*clen ) )
1732
1733 if (lidebug.ge.1) then
1734 call write_debug( name, 'sendsize, ubound(csendbuffer,1) ', &
1735 sendsize, ubound(csendbuffer,1) )
1736 call write_debug( name, 'ubound(xsend,1) ', &
1737 ubound(xsend,1) )
1738 call write_debug( name, 'ubound(sendproc,1) ', &
1739 ubound(sendproc,1) )
1740 call write_debug( name, 'ubound(sendtag,1) ', &
1741 ubound(sendtag,1) )
1742 endif
1743
1744
1745 = MPI_CHARACTER
1746 comm = communicator
1747
1748 do ii=1,nsend
1749
1750 = xsend(ii)
1751 j2 = xsend(ii+1)-1
1752
1753 count = j2-j1+1
1754 count = count*clen
1755
1756 do jj=j1,j2
1757 ijk = sendijk( jj )
1758 do ic=1,clen
1759 jpos = (jj-1)*clen + ic
1760 csendbuffer(jpos) = XX(ijk)(ic:ic)
1761 enddo
1762 enddo
1763
1764 dest = sendproc( ii )
1765 tag = sendtag( ii )
1766
1767 if (lidebug.ge.2) then
1768 call write_debug(name, 'mpi_isend: ii,j1,j2 ', ii,j1,j2)
1769 call write_debug(name, 'count, dest, tag ', count, &
1770 dest, tag )
1771 endif
1772
1773 jpos = (j1-1)*clen + 1
1774 call MPI_ISEND( csendbuffer(jpos), count, datatype, dest, &
1775 tag, comm, request, ierror )
1776 call MPI_Check( 'sendrecv_begin_1c:MPI_ISEND ', ierror )
1777 sendrequest( ii ) = request
1778 enddo
1779 endif
1780
1781
1782 return
1783 end subroutine sendrecv_begin_1c
1784
1785
1786
1787
1788
1789
1790 subroutine sendrecv_end_1d( XX, idebug )
1791
1792 use functions
1793
1794 implicit none
1795
1796
1797
1798 double precision, intent(inout), dimension(:) :: XX
1799 integer, intent(in), optional :: idebug
1800
1801
1802 interface
1803 subroutine MPI_WAITANY(count, array_of_requests, jindex, &
1804 status, ierror)
1805 use mpi, only: MPI_STATUS_SIZE
1806 integer count
1807 integer array_of_requests(*)
1808 integer jindex
1809 integer status(MPI_STATUS_SIZE)
1810 integer ierror
1811 end subroutine MPI_WAITANY
1812
1813 subroutine MPI_WAITALL( count, array_of_requests, &
1814 array_of_status, ierror )
1815 use mpi, only: MPI_STATUS_SIZE
1816 integer count
1817 integer array_of_requests(*)
1818 integer array_of_status( MPI_STATUS_SIZE,*)
1819 integer ierror
1820 end subroutine MPI_WAITALL
1821 end interface
1822
1823
1824
1825
1826 character(len=80), parameter :: name = 'sendrecv_end_1d'
1827 logical, parameter :: use_waitany = .false.
1828 integer :: lidebug
1829 integer :: jj, ijk, jindex, ii, j1, j2, ierror
1830 integer, dimension(MPI_STATUS_SIZE) :: recv_status_any
1831 integer, dimension(:,:), pointer :: recv_status
1832 integer, dimension(:,:), pointer :: send_status
1833
1834
1835
1836
1837 = 0
1838 if (present(idebug)) then
1839 lidebug = idebug
1840 endif
1841
1842 if (nsend.ge.1) then
1843 if (lidebug.ge.1) then
1844 call write_debug(name, &
1845 'waiting for sends to complete, nsend = ', nsend )
1846 endif
1847
1848 allocate( send_status(MPI_STATUS_SIZE,nsend))
1849
1850 if (use_persistent_message) then
1851 call MPI_WAITALL( nsend, send_persistent_request, &
1852 send_status, ierror )
1853 else
1854 call MPI_WAITALL( nsend, sendrequest, send_status, ierror )
1855 endif
1856
1857 call MPI_Check( 'sendrecv_end_1d:MPI_WAITALL ', ierror )
1858
1859 deallocate( send_status )
1860 nullify( send_status )
1861 endif
1862
1863
1864
1865 if (nrecv.ge.1) then
1866 if (lidebug.ge.1) then
1867 call write_debug( name, &
1868 'waiting for receives to complete, nrecv = ', nrecv )
1869 endif
1870
1871 if (use_waitany) then
1872 do ii=1,nrecv
1873 if (use_persistent_message) then
1874 call MPI_WAITANY( nrecv, recv_persistent_request, &
1875 jindex, recv_status_any, ierror )
1876 else
1877 call MPI_WAITANY( nrecv, recvrequest, &
1878 jindex, recv_status_any, ierror )
1879 endif
1880
1881 call MPI_Check( 'sendrecv_end_1d:MPI_WAITANY ', ierror )
1882
1883 j1 = xrecv( jindex )
1884 j2 = xrecv( jindex + 1)-1
1885
1886 if (lidebug.ge.2) then
1887 call write_debug(name, 'jindex, j1,j2 ', jindex,j1,j2)
1888 endif
1889
1890 do jj=j1,j2
1891 ijk = recvijk( jj )
1892 XX(ijk) = drecvbuffer(jj)
1893 enddo
1894 enddo
1895
1896 else
1897
1898 allocate( recv_status(MPI_STATUS_SIZE,nrecv))
1899 if (use_persistent_message) then
1900 call MPI_WAITALL( nrecv, recv_persistent_request, &
1901 recv_status, ierror )
1902 else
1903 call MPI_WAITALL( nrecv, recvrequest, recv_status,&
1904 ierror )
1905 endif
1906
1907 call MPI_Check( 'sendrecv_end_1d:MPI_WAITALL recv ', &
1908 ierror )
1909 deallocate( recv_status )
1910 nullify( recv_status )
1911
1912 j1 = xrecv(1)
1913 j2 = xrecv( nrecv +1)-1
1914 do jj=j1,j2
1915 ijk = recvijk( jj )
1916 XX(ijk) = drecvbuffer(jj)
1917 enddo
1918 endif
1919 endif
1920
1921 return
1922 end subroutine sendrecv_end_1d
1923
1924
1925
1926
1927
1928
1929 subroutine sendrecv_end_1c( XX, idebug )
1930
1931 use functions
1932
1933 implicit none
1934
1935
1936
1937
1938 character(len=*), intent(inout), dimension(:) :: XX
1939 integer, intent(in), optional :: idebug
1940
1941
1942 interface
1943 subroutine MPI_WAITANY(count, array_of_requests, jindex, &
1944 status, ierror)
1945 use mpi, only: MPI_STATUS_SIZE
1946 integer count
1947 integer array_of_requests(*)
1948 integer jindex
1949 integer status(MPI_STATUS_SIZE)
1950 integer ierror
1951 end subroutine MPI_WAITANY
1952
1953 subroutine MPI_WAITALL( count, array_of_requests, &
1954 array_of_status, ierror )
1955 use mpi, only: MPI_STATUS_SIZE
1956 integer count
1957 integer array_of_requests(*)
1958 integer array_of_status( MPI_STATUS_SIZE,*)
1959 integer ierror
1960 end subroutine MPI_WAITALL
1961 end interface
1962
1963
1964
1965
1966 character(len=80), parameter :: name = 'sendrecv_end_1c'
1967 integer :: ic, clen, jpos
1968 logical, parameter :: use_waitany = .false.
1969 integer :: lidebug
1970 integer :: jj, ijk, jindex, ii, j1, j2, ierror
1971 integer, dimension(MPI_STATUS_SIZE) :: recv_status_any
1972 integer, dimension(:,:), pointer :: recv_status
1973 integer, dimension(:,:), pointer :: send_status
1974
1975
1976
1977
1978 = 0
1979 if (present(idebug)) then
1980 lidebug = idebug
1981 endif
1982
1983 jpos = lbound(XX,1)
1984 clen = len(XX(jpos))
1985
1986 if (nsend.ge.1) then
1987 if (lidebug.ge.1) then
1988 call write_debug(name, &
1989 'waiting for sends to complete, nsend = ', nsend )
1990 endif
1991
1992 allocate( send_status(MPI_STATUS_SIZE,nsend))
1993
1994 call MPI_WAITALL( nsend, sendrequest, send_status, ierror )
1995 call MPI_Check( 'sendrecv_end_1c:MPI_WAITALL ', ierror )
1996
1997 deallocate( send_status )
1998 nullify( send_status )
1999
2000 deallocate( csendbuffer )
2001 nullify( csendbuffer )
2002
2003 endif
2004
2005
2006
2007 if (nrecv.ge.1) then
2008 if (lidebug.ge.1) then
2009 call write_debug( name, &
2010 'waiting for receives to complete, nrecv = ', nrecv )
2011 endif
2012
2013 if (use_waitany) then
2014 do ii=1,nrecv
2015 call MPI_WAITANY( nrecv, recvrequest, jindex, &
2016 recv_status_any, ierror )
2017 call MPI_Check( 'sendrecv_end_1c:MPI_WAITANY ', ierror )
2018
2019 j1 = xrecv( jindex )
2020 j2 = xrecv( jindex + 1)-1
2021
2022 if (lidebug.ge.2) then
2023 call write_debug(name, 'jindex, j1,j2 ', jindex,j1,j2 )
2024 endif
2025
2026 do jj=j1,j2
2027 ijk = recvijk( jj )
2028
2029 do ic=1,clen
2030 jpos = (jj-1)*clen + ic
2031 XX(ijk)(ic:ic) = crecvbuffer(jpos)
2032 enddo
2033 enddo
2034 enddo
2035
2036 else
2037
2038 allocate( recv_status(MPI_STATUS_SIZE,nrecv))
2039 call MPI_WAITALL( nrecv, recvrequest, recv_status, ierror )
2040 call MPI_Check( 'sendrecv_end_1c:MPI_WAITALL recv ', ierror )
2041
2042 deallocate( recv_status )
2043 nullify( recv_status )
2044
2045 j1 = xrecv(1)
2046 j2 = xrecv( nrecv +1)-1
2047 do jj=j1,j2
2048 ijk = recvijk( jj )
2049 do ic=1,clen
2050 jpos = (jj-1)*clen + ic
2051 XX(ijk)(ic:ic) = crecvbuffer(jpos)
2052 enddo
2053 enddo
2054 endif
2055
2056 deallocate( crecvbuffer )
2057 nullify( crecvbuffer )
2058
2059 endif
2060
2061 return
2062 end subroutine sendrecv_end_1c
2063
2064
2065
2066
2067
2068 subroutine sendrecv_end_1i( XX, idebug )
2069
2070 use functions
2071
2072 implicit none
2073
2074
2075
2076 integer, intent(inout), dimension(:) :: XX
2077 integer, intent(in), optional :: idebug
2078
2079 interface
2080 subroutine MPI_WAITANY(count, array_of_requests, jindex, &
2081 status, ierror)
2082 use mpi, only: MPI_STATUS_SIZE
2083 integer count
2084 integer array_of_requests(*)
2085 integer jindex
2086 integer status(MPI_STATUS_SIZE)
2087 integer ierror
2088 end subroutine MPI_WAITANY
2089
2090 subroutine MPI_WAITALL( count, array_of_requests, &
2091 array_of_status, ierror )
2092 use mpi, only: MPI_STATUS_SIZE
2093 integer count
2094 integer array_of_requests(*)
2095 integer array_of_status( MPI_STATUS_SIZE,*)
2096 integer ierror
2097 end subroutine MPI_WAITALL
2098 end interface
2099
2100
2101
2102
2103 character(len=80), parameter :: name = 'sendrecv_end_1i'
2104 logical, parameter :: use_waitany = .false.
2105 integer :: lidebug
2106 integer :: jj, ijk, jindex, ii, j1, j2, ierror
2107 integer, dimension(MPI_STATUS_SIZE) :: recv_status_any
2108 integer, dimension(:,:), pointer :: recv_status
2109 integer, dimension(:,:), pointer :: send_status
2110
2111
2112
2113
2114 = 0
2115 if (present(idebug)) then
2116 lidebug = idebug
2117 endif
2118
2119 if (nsend.ge.1) then
2120 if (lidebug.ge.1) then
2121 call write_debug(name, &
2122 'waiting for sends to complete, nsend = ', nsend )
2123 endif
2124
2125 allocate( send_status(MPI_STATUS_SIZE,nsend))
2126
2127 call MPI_WAITALL( nsend, sendrequest, send_status, ierror )
2128 call MPI_Check( 'sendrecv_end_1i:MPI_WAITALL ', ierror )
2129
2130 deallocate( send_status )
2131 nullify( send_status )
2132
2133 deallocate( isendbuffer )
2134 nullify( isendbuffer )
2135 endif
2136
2137
2138
2139
2140 if (nrecv.ge.1) then
2141 if (lidebug.ge.1) then
2142 call write_debug( name, &
2143 'waiting for receives to complete, nrecv = ', nrecv )
2144 endif
2145
2146 if (use_waitany) then
2147 do ii=1,nrecv
2148 call MPI_WAITANY( nrecv, recvrequest, jindex, &
2149 recv_status_any, ierror )
2150 call MPI_Check( 'sendrecv_end_1i:MPI_WAITANY ', ierror )
2151
2152 j1 = xrecv( jindex )
2153 j2 = xrecv( jindex + 1)-1
2154
2155 if (lidebug.ge.2) then
2156 call write_debug(name, 'jindex, j1,j2 ', jindex,j1,j2 )
2157 endif
2158
2159 do jj=j1,j2
2160 ijk = recvijk( jj )
2161 XX(ijk) = irecvbuffer(jj)
2162 enddo
2163 enddo
2164 else
2165 allocate( recv_status(MPI_STATUS_SIZE,nrecv))
2166 call MPI_WAITALL( nrecv, recvrequest, recv_status, ierror )
2167 call MPI_Check( 'sendrecv_end_1i:MPI_WAITALL recv ', ierror )
2168 deallocate( recv_status )
2169 nullify( recv_status )
2170
2171 j1 = xrecv(1)
2172 j2 = xrecv( nrecv +1)-1
2173 do jj=j1,j2
2174 ijk = recvijk( jj )
2175 XX(ijk) = irecvbuffer(jj)
2176 enddo
2177 endif
2178
2179 deallocate( irecvbuffer )
2180 nullify( irecvbuffer )
2181
2182 endif
2183
2184 return
2185 end subroutine sendrecv_end_1i
2186
2187
2188
2189
2190
2191
2192
2193 subroutine send_recv_1c( XX, ilayer, idebug )
2194 implicit none
2195
2196
2197
2198 character(len=*), dimension(:), intent(inout) :: XX
2199 integer, intent(in), optional :: ilayer,idebug
2200
2201
2202
2203 integer :: lidebug, layer
2204
2205
2206 = 0
2207 if (present(idebug)) then
2208 lidebug = idebug
2209 endif
2210
2211 layer = 1
2212 if (present(ilayer)) then
2213 layer = ilayer
2214 endif
2215
2216 call sendrecv_begin(XX,layer,lidebug)
2217 call sendrecv_end( XX, lidebug )
2218
2219 return
2220 end subroutine send_recv_1c
2221
2222
2223
2224
2225
2226
2227 subroutine send_recv_1d( XX, ilayer, idebug )
2228 implicit none
2229
2230
2231
2232 double precision, dimension(:), intent(inout) :: XX
2233 integer, intent(in), optional :: ilayer,idebug
2234
2235
2236
2237 integer :: lidebug, layer
2238
2239
2240 = 0
2241 if (present(idebug)) then
2242 lidebug = idebug
2243 endif
2244
2245 layer = 1
2246 if (present(ilayer)) then
2247 layer = ilayer
2248 endif
2249
2250 call sendrecv_begin(XX,layer,lidebug)
2251 call sendrecv_end( XX, lidebug )
2252
2253 return
2254 end subroutine send_recv_1d
2255
2256
2257
2258
2259
2260
2261 subroutine send_recv_2d( XX, ilayer, idebug )
2262 implicit none
2263
2264
2265
2266 double precision, dimension(:,:), intent(inout) :: XX
2267 integer, intent(in), optional :: ilayer,idebug
2268
2269
2270
2271 integer :: lidebug, layer
2272 integer :: j
2273
2274
2275 = 0
2276 if (present(idebug)) then
2277 lidebug = idebug
2278 endif
2279
2280 layer = 1
2281 if (present(ilayer)) then
2282 layer = ilayer
2283 endif
2284
2285 do j=lbound(XX,2),ubound(XX,2)
2286 call sendrecv_begin(XX(:,j),layer,lidebug)
2287 call sendrecv_end( XX(:,j), lidebug )
2288 enddo
2289
2290 return
2291 end subroutine send_recv_2d
2292
2293
2294
2295
2296
2297
2298 subroutine send_recv_3d( XX, ilayer, idebug )
2299 implicit none
2300
2301
2302
2303 double precision, dimension(:,:,:), intent(inout) :: XX
2304 integer, intent(in), optional :: ilayer,idebug
2305
2306
2307
2308 integer :: lidebug, layer
2309 integer :: j,k
2310
2311 = 0
2312 if (present(idebug)) then
2313 lidebug = idebug
2314 endif
2315
2316 layer = 1
2317 if (present(ilayer)) then
2318 layer = ilayer
2319 endif
2320
2321 do k=lbound(XX,3),ubound(XX,3)
2322 do j=lbound(XX,2),ubound(XX,2)
2323 call sendrecv_begin(XX(:,j,k),layer,lidebug)
2324 call sendrecv_end( XX(:,j,k), lidebug )
2325 enddo
2326 enddo
2327
2328 return
2329 end subroutine send_recv_3d
2330
2331
2332
2333
2334
2335
2336 subroutine send_recv_1i( XX, ilayer, idebug )
2337 implicit none
2338
2339
2340
2341 integer, dimension(:), intent(inout) :: XX
2342 integer, intent(in), optional :: ilayer,idebug
2343
2344
2345
2346 integer :: lidebug, layer
2347
2348 = 0
2349 if (present(idebug)) then
2350 lidebug = idebug
2351 endif
2352
2353 layer = 1
2354 if (present(ilayer)) then
2355 layer = ilayer
2356 endif
2357
2358 call sendrecv_begin(XX,layer,lidebug)
2359 call sendrecv_end( XX, lidebug )
2360
2361 return
2362 end subroutine send_recv_1i
2363
2364
2365
2366
2367 subroutine sendrecv_re_init_after_re_indexing(comm, idebug )
2368
2369 use functions
2370
2371 implicit none
2372
2373 integer, intent(in) :: comm
2374
2375 integer, intent(in), optional :: idebug
2376
2377
2378
2379
2380
2381
2382
2383
2384
2385 character(len=80), parameter :: name = 'sendrecv_init'
2386
2387 character(len=80), pointer, dimension(:) :: line
2388 integer :: ip, lmax
2389
2390 integer :: layer,request, source, tag, datatype
2391
2392 integer :: lidebug
2393 integer :: isize,jsize,ksize, ijksize
2394 integer :: recvsize1, recvsize2, &
2395 sendsize1, sendsize2
2396
2397 integer :: iter, i,j,k, ii, jj,kk, &
2398 ntotal, icount,ipos, &
2399 ilayer, i1,i2, j1,j2, k1,k2, &
2400 ijk, ijk2, iproc, jproc, src,dest, &
2401 ierror
2402
2403 logical :: isok, isvalid, ismine, is_halobc
2404
2405 integer, dimension(:,:,:), pointer :: ijk2proc
2406 integer, pointer, dimension(:) :: &
2407 istartx,iendx, jstartx,jendx, kstartx,kendx, &
2408 ncount, &
2409 recvproc, recvtag, xrecv, recvijk, &
2410 sendproc, sendtag, xsend, sendijk
2411
2412 logical, parameter :: jfastest = .true.
2413
2414
2415 integer, parameter :: message_tag_offset = 11
2416
2417
2418
2419
2420
2421 integer :: message_tag
2422
2423
2424 if (use_persistent_message) then
2425
2426 datatype = MPI_DOUBLE_PRECISION
2427
2428 do layer=1,2
2429
2430
2431 if (layer.eq.1) then
2432 nrecv = nrecv1
2433 recvtag =>recvtag1
2434 recvproc => recvproc1
2435 recvijk => recvijk1
2436 xrecv => xrecv1
2437
2438 nsend = nsend1
2439 sendtag => sendtag1
2440 sendproc => sendproc1
2441 sendijk => sendijk1
2442 xsend => xsend1
2443
2444 send_persistent_request => send_persistent_request1
2445 recv_persistent_request => recv_persistent_request1
2446
2447 else
2448 nrecv = nrecv2
2449 recvtag =>recvtag2
2450 recvproc => recvproc2
2451 recvijk => recvijk2
2452 xrecv => xrecv2
2453
2454 nsend = nsend2
2455 sendtag => sendtag2
2456 sendproc => sendproc2
2457 sendijk => sendijk2
2458 xsend => xsend2
2459
2460 send_persistent_request => send_persistent_request2
2461 recv_persistent_request => recv_persistent_request2
2462
2463 endif
2464
2465
2466
2467 do ii=1,nrecv
2468 j1 = xrecv(ii)
2469 j2 = xrecv(ii+1)-1
2470 icount = j2-j1+1
2471 source = recvproc( ii )
2472 tag = recvtag( ii )
2473
2474
2475
2476 if (lidebug.ge.2) then
2477
2478
2479
2480
2481
2482 endif
2483
2484
2485 call MPI_RECV_INIT( drecvbuffer(j1), icount, datatype, &
2486 source, tag, comm, request, ierror )
2487 call MPI_Check( 'sendrecv_begin_1d:MPI_IRECV ', ierror )
2488
2489 recv_persistent_request(ii) = request
2490 enddo
2491
2492
2493 do ii=1,nsend
2494 j1 = xsend(ii)
2495 j2 = xsend(ii+1)-1
2496 dest = sendproc( ii )
2497 tag = sendtag( ii )
2498 icount = j2-j1+1
2499
2500 if (lidebug.ge.2) then
2501
2502
2503
2504
2505
2506 endif
2507
2508
2509 call MPI_SEND_INIT( dsendbuffer(j1), icount, datatype, &
2510 dest, tag, &
2511 comm, request, ierror )
2512 call MPI_Check( 'sendrecv_begin_1d:MPI_SEND_INIT ', ierror )
2513
2514 send_persistent_request( ii ) = request
2515 enddo
2516
2517 enddo
2518
2519 endif
2520
2521
2522 return
2523 end subroutine sendrecv_re_init_after_re_indexing
2524
2525
2526 end module sendrecv
2527