File: N:\mfix\model\dmp_modules\sendrecv_mod.f

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