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