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

1     !       MPI Modules written at ORNL by Ed and Sreekanth for MFIX
2     !       under joint effort with FETC - 06/08/99.
3     
4     #include "version.inc"
5     
6     module mpi_utility
7     
8       !       module to perform most of the mpi functionalities like scatter,
9       !       gather, bcast, globalsum and so on.
10     
11       use geometry
12       use compar
13       use parallel_mpi
14       use debug
15       use indices
16       implicit none
17     
18       !       Object-oriented approach to direct to the correct procedure
19       !       depending on the argument type. i stands for integer, r for real
20       !       and d for double precision. 0 for scalar, 1 for vector, 2 for
21       !       2-D array and similarly 3.
22     
23       !==============================================================================
24       !  JFD: Interfaces used for vtk file writting (Cartesian grid):
25       !==============================================================================
26     
27       interface allgather
28          module  procedure allgather_1i
29       end interface allgather
30     
31       interface gatherv
32          module  procedure gatherv_1d
33       end interface gatherv
34     
35       !==============================================================================
36       !  JFD: End of Interfaces used for vtk file writting (Cartesian grid):
37       !==============================================================================
38     
39       interface scatter
40          module  procedure scatter_1i, scatter_2i, scatter_3i, &
41               scatter_1r, scatter_2r, scatter_3r, &
42               scatter_1d, scatter_2d, scatter_3d, &
43               scatter_1c,scatter_1l
44       end interface scatter
45     
46       interface gather
47          module  procedure gather_1i, gather_2i, gather_3i, &
48               gather_1r, gather_2r, gather_3r, &
49               gather_1d, gather_2d, gather_3d, &
50               gather_1c, gather_1l
51       end interface gather
52     
53       interface bcast
54          module  procedure bcast_0i, bcast_1i, bcast_2i, bcast_3i, &
55               bcast_0r, bcast_1r, bcast_2r, bcast_3r, &
56               bcast_0d, bcast_1d, bcast_2d, bcast_3d, &
57               bcast_0l, bcast_1l, bcast_0c, bcast_1c
58       end interface bcast
59     
60       interface global_sum
61          module  procedure global_sum_0i, global_sum_1i, global_sum_2i, global_sum_3i, &
62               global_sum_0r, global_sum_1r, global_sum_2r, global_sum_3r, &
63               global_sum_0d, global_sum_1d, global_sum_2d, global_sum_3d
64       end interface global_sum
65     
66       interface global_all_sum
67          module  procedure &
68               global_all_sum_0i, global_all_sum_1i, &
69               global_all_sum_2i, global_all_sum_3i, &
70               global_all_sum_0r, global_all_sum_1r, &
71               global_all_sum_2r, global_all_sum_3r, &
72               global_all_sum_0d, global_all_sum_1d, &
73               global_all_sum_2d, global_all_sum_3d, &
74               global_all_sum_onevar_0i, global_all_sum_onevar_1i, &
75               global_all_sum_onevar_2i, global_all_sum_onevar_3i, &
76               global_all_sum_onevar_0r, global_all_sum_onevar_1r, &
77               global_all_sum_onevar_2r, global_all_sum_onevar_3r, &
78               global_all_sum_onevar_0d, global_all_sum_onevar_1d, &
79               global_all_sum_onevar_2d, global_all_sum_onevar_3d
80       end interface global_all_sum
81     
82       interface global_min
83          module  procedure global_min_0i, global_min_1i, global_min_2i, global_min_3i, &
84               global_min_0r, global_min_1r, global_min_2r, global_min_3r, &
85               global_min_0d, global_min_1d, global_min_2d, global_min_3d
86       end interface global_min
87     
88       interface global_all_min
89          module  procedure &
90               global_all_min_0i, global_all_min_1i, &
91               global_all_min_2i, global_all_min_3i, &
92               global_all_min_0r, global_all_min_1r, &
93               global_all_min_2r, global_all_min_3r, &
94               global_all_min_0d, global_all_min_1d, &
95               global_all_min_2d, global_all_min_3d, &
96               global_all_min_onevar_0i, global_all_min_onevar_1i, &
97               global_all_min_onevar_2i, global_all_min_onevar_3i, &
98               global_all_min_onevar_0r, global_all_min_onevar_1r, &
99               global_all_min_onevar_2r, global_all_min_onevar_3r, &
100               global_all_min_onevar_0d, global_all_min_onevar_1d, &
101               global_all_min_onevar_2d, global_all_min_onevar_3d
102       end interface global_all_min
103     
104       interface global_max
105          module  procedure global_max_0i, global_max_1i, global_max_2i, global_max_3i, &
106               global_max_0r, global_max_1r, global_max_2r, global_max_3r, &
107               global_max_0d, global_max_1d, global_max_2d, global_max_3d
108       end interface global_max
109     
110       interface global_all_max
111          module  procedure &
112               global_all_max_0i, global_all_max_1i, &
113               global_all_max_2i, global_all_max_3i, &
114               global_all_max_0r, global_all_max_1r, &
115               global_all_max_2r, global_all_max_3r, &
116               global_all_max_0d, global_all_max_1d, &
117               global_all_max_2d, global_all_max_3d, &
118               global_all_max_onevar_0i, global_all_max_onevar_1i, &
119               global_all_max_onevar_2i, global_all_max_onevar_3i, &
120               global_all_max_onevar_0r, global_all_max_onevar_1r, &
121               global_all_max_onevar_2r, global_all_max_onevar_3r, &
122               global_all_max_onevar_0d, global_all_max_onevar_1d, &
123               global_all_max_onevar_2d, global_all_max_onevar_3d
124       end interface global_all_max
125     
126       interface global_all_and
127          module procedure &
128               global_all_and_0d, global_all_and_1d, &
129               global_all_and_onevar_0d, global_all_and_onevar_1d
130       end interface global_all_and
131     
132       interface global_all_or
133          module procedure &
134               global_all_or_0d, global_all_or_1d, &
135               global_all_or_onevar_0d, global_all_or_onevar_1d
136       end interface global_all_or
137     
138     contains
139     
140     
141       !==============================================================================
142       !  JFD: Subroutines used for vtk file writting (Cartesian grid):
143       !==============================================================================
144     
145       subroutine allgather_1i( lbuf, gbuf, idebug )
146         integer, intent(in) :: lbuf
147         integer, intent(out), dimension(:) :: gbuf
148         integer, optional, intent(in) ::  idebug
149     
150     #ifdef MPI
151         integer :: sendtype,recvtype,sendcnt,recvcnt,ierr,lidebug,mpierr
152     
153         if (.not. present(idebug)) then
154            lidebug = 0
155         else
156            lidebug = idebug
157         endif
158     
159         recvtype = MPI_INTEGER
160         sendtype = recvtype
161     
162         sendcnt = 1
163         recvcnt = sendcnt
164     
165         CALL MPI_ALLGATHER(lbuf,sendcnt,sendtype,  &
166              gbuf,recvcnt,recvtype,MPI_COMM_WORLD, IERR)
167     #else
168         gbuf = 0
169     #endif
170     
171         return
172       end subroutine allgather_1i
173     
174       subroutine allgather_1d( lbuf, gbuf, idebug )
175         double precision, intent(in) :: lbuf
176         double precision, intent(out), dimension(:) :: gbuf
177         integer, optional, intent(in) ::  idebug
178     
179     #ifdef MPI
180         integer :: sendtype,recvtype,sendcnt,recvcnt,ierr,lidebug,mpierr
181     
182         if (.not. present(idebug)) then
183            lidebug = 0
184         else
185            lidebug = idebug
186         endif
187     
188         recvtype = MPI_DOUBLE_PRECISION
189         sendtype = MPI_DOUBLE_PRECISION
190     
191         sendcnt = 1
192         recvcnt = sendcnt
193     
194         CALL MPI_ALLGATHER(lbuf,sendcnt,sendtype,  &
195              gbuf,recvcnt,recvtype,MPI_COMM_WORLD, IERR)
196     #else
197         gbuf = 0
198     #endif
199     
200         return
201       end subroutine allgather_1d
202     
203       subroutine gatherv_1i( lbuf, sendcnt, gbuf, rcount, disp, mroot, idebug )
204         integer, intent(in), dimension(:) :: lbuf
205         integer, intent(in), dimension(:) :: rcount
206         integer, intent(in), dimension(:) :: disp
207         integer, intent(out), dimension(:) :: gbuf
208         integer, optional, intent(in) :: mroot, idebug
209         integer :: sendtype,recvtype,sendcnt,recvcnt,lroot,ierr,lidebug
210     
211     #ifdef MPI
212         !       check to see whether there is root
213     
214         if (.not. present(mroot)) then
215            lroot = 0
216         else
217            lroot = mroot
218         endif
219     
220         if (.not. present(idebug)) then
221            lidebug = 0
222         else
223            lidebug = idebug
224         endif
225     
226         recvtype = MPI_INTEGER
227         sendtype = MPI_INTEGER
228     
229         CALL MPI_GATHERV(lbuf,sendcnt,sendtype,  &
230              gbuf,rcount,disp,recvtype, &
231              lroot,MPI_COMM_WORLD, IERR)
232     #else
233         gbuf = lbuf
234     #endif
235     
236         return
237       end subroutine gatherv_1i
238     
239       subroutine gatherv_1d( lbuf, sendcnt, gbuf, rcount, disp, mroot, idebug )
240         double precision, intent(in), dimension(:) :: lbuf
241         integer, intent(in), dimension(:) :: rcount
242         integer, intent(in), dimension(:) :: disp
243         double precision, intent(out), dimension(:) :: gbuf
244         integer, optional, intent(in) :: mroot, idebug
245         integer :: sendtype,recvtype,sendcnt,recvcnt,lroot,ierr,lidebug
246     
247     #ifdef MPI
248         !       check to see whether there is root
249     
250         if (.not. present(mroot)) then
251            lroot = 0
252         else
253            lroot = mroot
254         endif
255     
256         if (.not. present(idebug)) then
257            lidebug = 0
258         else
259            lidebug = idebug
260         endif
261     
262         recvtype = MPI_DOUBLE_PRECISION
263         sendtype = MPI_DOUBLE_PRECISION
264     
265         CALL MPI_GATHERV(lbuf,sendcnt,sendtype,  &
266              gbuf,rcount,disp,recvtype, &
267              lroot,MPI_COMM_WORLD, IERR)
268     #else
269         gbuf = lbuf
270     #endif
271     
272         return
273       end subroutine gatherv_1d
274     
275       !==============================================================================
276       !  JFD: End of Subroutines used for vtk file writting (Cartesian grid):
277       !==============================================================================
278     
279     
280       !       Routine to scatter gbuf available on root to all the processors
281     
282       subroutine scatter_1i( lbuf, gbuf, mroot, idebug )
283     
284         use functions
285     
286         implicit none
287     
288         integer, intent(in), dimension(:) :: gbuf
289         integer, intent(out), dimension(:) :: lbuf
290         integer, optional, intent(in) :: mroot, idebug
291     
292         integer, allocatable, dimension(:) :: gbuf_pack
293     
294         integer :: sendtype, recvtype, ijk1, ijk2, recvcnt, ierr,lroot, lidebug
295         integer :: i,j,k,ibuffer,iproc, ioffset
296         integer :: ijk
297     
298     #ifdef MPI
299         !       check to see whether there is root
300     
301         if (.not. present(mroot)) then
302            lroot = 0
303         else
304            lroot = mroot
305         endif
306     
307         if (.not. present(idebug)) then
308            lidebug = 0
309         else
310            lidebug = idebug
311         endif
312     
313         if(myPE.eq.lroot) then
314            allocate(gbuf_pack(sum(ijksize3_all(:))))
315         else
316            allocate(gbuf_pack(10))
317         endif
318     
319         if( myPE.eq.lroot) then
320            ioffset = 0
321            do iproc = 0,numPEs-1
322               ibuffer = 0
323               do k = kstart3_all(iproc), kend3_all(iproc)
324                  do j = jstart3_all(iproc), jend3_all(iproc)
325                     do i = istart3_all(iproc), iend3_all(iproc)
326     
327                        ibuffer = funijk_proc(i,j,k,iproc) + ioffset
328                        gbuf_pack(ibuffer) = gbuf(funijk_gl(i,j,k))
329     
330                     enddo
331                  enddo
332               enddo
333               ioffset = ibuffer
334            enddo
335         endif
336     
337         sendtype = MPI_INTEGER
338         recvtype = sendtype
339     
340         ijk1 = ijkstart3
341         ijk2 = ijkend3
342     
343         recvcnt = ijk2-ijk1+1
344     
345         !       Call MPI routines
346     
347         call MPI_Scatterv( gbuf_pack, ijksize3_all, displs, sendtype, &
348              lbuf, recvcnt, recvtype,  &
349              lroot, MPI_COMM_WORLD, ierr )
350         call MPI_Check( 'scatter_1i:MPI_Scatterv', ierr )
351     
352         deallocate(gbuf_pack)
353     #else
354         lbuf = gbuf
355     #endif
356     
357         return
358       end subroutine scatter_1i
359     
360       subroutine scatter_2i( lbuf, gbuf, mroot, idebug )
361         integer, intent(in), dimension(:,:) :: gbuf
362         integer, intent(out), dimension(:,:) :: lbuf
363         integer, optional, intent(in) :: mroot, idebug
364     
365     #ifdef MPI
366         integer :: i,j,lroot, lidebug
367     
368         if (.not. present(mroot)) then
369            lroot = 0
370         else
371            lroot = mroot
372         endif
373     
374         if (.not. present(idebug)) then
375            lidebug = 0
376         else
377            lidebug = idebug
378         endif
379     
380         if(myPE.eq.lroot) then
381            call assert( size(lbuf,2).eq.size(gbuf,2),  &
382                 '** scatter_2i: size(lbuf,2).ne.size(gbuf,2) ', &
383                 size(lbuf,2), size(gbuf,2) )
384         endif
385     
386         do j=lbound(lbuf,2),ubound(lbuf,2)
387            call scatter_1i( lbuf(:,j), gbuf(:,j), lroot, lidebug )
388         enddo
389     #else
390         lbuf = gbuf
391     #endif
392     
393         return
394       end subroutine scatter_2i
395     
396       subroutine scatter_3i( lbuf, gbuf, mroot, idebug )
397         integer, intent(in), dimension(:,:,:) :: gbuf
398         integer, intent(out), dimension(:,:,:) :: lbuf
399         integer, optional, intent(in) :: mroot, idebug
400     
401     #ifdef MPI
402         integer :: j,k,lroot, lidebug
403     
404         if (.not. present(mroot)) then
405            lroot = 0
406         else
407            lroot = mroot
408         endif
409     
410         if (.not. present(idebug)) then
411            lidebug = 0
412         else
413            lidebug = idebug
414         endif
415     
416         if(myPE.eq.lroot) then
417            call assert( size(lbuf,2).eq.size(gbuf,2),  &
418                 '** scatter_3i: size(lbuf,2).ne.size(gbuf,2) ', &
419                 size(lbuf,2), size(gbuf,2) )
420     
421            call assert( size(lbuf,3).eq.size(gbuf,3),  &
422                 '** scatter_3i: size(lbuf,3).ne.size(gbuf,3) ', &
423                 size(lbuf,3), size(gbuf,3) )
424         endif
425     
426         do k=lbound(lbuf,3),ubound(lbuf,3)
427            do j=lbound(lbuf,2),ubound(lbuf,2)
428               call scatter_1i( lbuf(:,j,k), gbuf(:,j,k), lroot, lidebug )
429            enddo
430         enddo
431     #else
432         lbuf = gbuf
433     #endif
434     
435         return
436       end subroutine scatter_3i
437     
438       subroutine scatter_1r( lbuf, gbuf, mroot, idebug )
439     
440         use functions
441     
442         implicit none
443     
444         real, intent(in), dimension(:) :: gbuf
445         real, intent(out), dimension(:) :: lbuf
446         integer, optional, intent(in) :: mroot, idebug
447     
448     #ifdef MPI
449         real, allocatable, dimension(:) :: gbuf_pack
450     
451         integer :: sendtype, recvtype, ijk1, ijk2, recvcnt, ierr,lroot, lidebug
452         integer :: i,j,k,ibuffer,iproc, ioffset
453         integer :: ijk
454     
455         if (.not. present(mroot)) then
456            lroot = 0
457         else
458            lroot = mroot
459         endif
460     
461         if (.not. present(idebug)) then
462            lidebug = 0
463         else
464            lidebug = idebug
465         endif
466     
467         if(myPE.eq.lroot) then
468            allocate(gbuf_pack(sum(ijksize3_all(:))))
469         else
470            allocate(gbuf_pack(10))
471         endif
472     
473         if( myPE.eq.lroot) then
474            ioffset = 0
475            do iproc = 0,numPEs-1
476               ibuffer = 0
477               do k = kstart3_all(iproc), kend3_all(iproc)
478                  do j = jstart3_all(iproc), jend3_all(iproc)
479                     do i = istart3_all(iproc), iend3_all(iproc)
480     
481                        ibuffer = funijk_proc(i,j,k,iproc) + ioffset
482                        gbuf_pack(ibuffer) = gbuf(funijk_gl(i,j,k))
483     
484                     enddo
485                  enddo
486               enddo
487               ioffset = ibuffer
488            enddo
489         endif
490     
491         sendtype = MPI_REAL
492         recvtype = sendtype
493     
494         ijk1 = ijkstart3
495         ijk2 = ijkend3
496     
497         recvcnt = ijk2-ijk1+1
498     
499         call MPI_Scatterv( gbuf_pack, ijksize3_all, displs, sendtype, &
500              lbuf, recvcnt, recvtype,  &
501              lroot, MPI_COMM_WORLD, ierr )
502         call MPI_Check( 'scatter_1r:MPI_Scatterv', ierr )
503     
504         deallocate(gbuf_pack)
505     #else
506         lbuf = gbuf
507     #endif
508     
509         return
510       end subroutine scatter_1r
511     
512     
513       subroutine scatter_2r( lbuf, gbuf, mroot, idebug )
514         real, intent(in), dimension(:,:) :: gbuf
515         real, intent(out), dimension(:,:) :: lbuf
516         integer, optional, intent(in) :: mroot, idebug
517     
518     #ifdef MPI
519         integer :: i,j,lroot, lidebug
520     
521         if (.not. present(mroot)) then
522            lroot = 0
523         else
524            lroot = mroot
525         endif
526     
527         if (.not. present(idebug)) then
528            lidebug = 0
529         else
530            lidebug = idebug
531         endif
532     
533         if(myPE.eq.lroot) then
534            call assert( size(lbuf,2).eq.size(gbuf,2),  &
535                 '** scatter_2r: size(lbuf,2).ne.size(gbuf,2) ', &
536                 size(lbuf,2), size(gbuf,2) )
537         endif
538     
539         do j=lbound(lbuf,2),ubound(lbuf,2)
540            call scatter_1r( lbuf(:,j), gbuf(:,j), lroot, lidebug )
541         enddo
542     #else
543         lbuf = gbuf
544     #endif
545     
546         return
547       end subroutine scatter_2r
548     
549       subroutine scatter_3r( lbuf, gbuf, mroot, idebug )
550         real, intent(in), dimension(:,:,:) :: gbuf
551         real, intent(out), dimension(:,:,:) :: lbuf
552         integer, optional, intent(in) :: mroot, idebug
553     
554     #ifdef MPI
555         integer :: j,k,lroot, lidebug
556     
557         if (.not. present(mroot)) then
558            lroot = 0
559         else
560            lroot = mroot
561         endif
562     
563         if (.not. present(idebug)) then
564            lidebug = 0
565         else
566            lidebug = idebug
567         endif
568     
569         if(myPE.eq.lroot) then
570            call assert( size(lbuf,2).eq.size(gbuf,2),  &
571                 '** scatter_3r: size(lbuf,2).ne.size(gbuf,2) ', &
572                 size(lbuf,2), size(gbuf,2) )
573     
574            call assert( size(lbuf,3).eq.size(gbuf,3),  &
575                 '** scatter_3r: size(lbuf,3).ne.size(gbuf,3) ', &
576                 size(lbuf,3), size(gbuf,3) )
577         endif
578     
579         do k=lbound(lbuf,3),ubound(lbuf,3)
580            do j=lbound(lbuf,2),ubound(lbuf,2)
581               call scatter_1r( lbuf(:,j,k), gbuf(:,j,k), lroot, lidebug )
582            enddo
583         enddo
584     #else
585         lbuf = gbuf
586     #endif
587     
588         return
589       end subroutine scatter_3r
590     
591     
592       subroutine scatter_1d( lbuf, gbuf, mroot, idebug )
593     
594         use functions
595         implicit none
596     
597         double precision, intent(in), dimension(:) :: gbuf
598         double precision, intent(out), dimension(:) :: lbuf
599         integer, optional, intent(in) :: mroot, idebug
600     
601     #ifdef MPI
602         double precision, allocatable, dimension(:) :: gbuf_pack
603     
604         integer :: sendtype, recvtype, ijk1,ijk2,recvcnt, ierr,lroot, lidebug
605         integer :: i,j,k,ibuffer,iproc, ioffset
606         integer :: ijk
607     
608         if (.not. present(mroot)) then
609            lroot = 0
610         else
611            lroot = mroot
612         endif
613     
614         if (.not. present(idebug)) then
615            lidebug = 0
616         else
617            lidebug = idebug
618         endif
619     
620         if(myPE.eq.lroot) then
621            allocate(gbuf_pack(sum(ijksize3_all(:))))
622         else
623            allocate(gbuf_pack(10))
624         endif
625     
626         if( myPE.eq.lroot) then
627            ioffset = 0
628            do iproc = 0,numPEs-1
629               ibuffer = 0
630               do k = kstart3_all(iproc), kend3_all(iproc)
631                  do j = jstart3_all(iproc), jend3_all(iproc)
632                     do i = istart3_all(iproc), iend3_all(iproc)
633     
634                        ibuffer = funijk_proc(i,j,k,iproc) + ioffset
635                        gbuf_pack(ibuffer) = gbuf(funijk_gl(i,j,k))
636     
637                     enddo
638                  enddo
639               enddo
640               ioffset = ibuffer
641            enddo
642         endif
643     
644         sendtype = MPI_DOUBLE_PRECISION
645         recvtype = sendtype
646     
647         ijk1 = ijkstart3
648         ijk2 = ijkend3
649     
650         recvcnt = ijk2-ijk1+1
651     
652         call MPI_Scatterv( gbuf_pack, ijksize3_all, displs, sendtype, &
653              lbuf, recvcnt, recvtype,  &
654              lroot, MPI_COMM_WORLD, ierr )
655         call MPI_Check( 'scatter_1d:MPI_Scatterv', ierr )
656     
657         deallocate(gbuf_pack)
658     #else
659         lbuf = gbuf
660     #endif
661     
662         return
663       end subroutine scatter_1d
664     
665     
666       subroutine scatter_2d( lbuf, gbuf, mroot, idebug )
667         double precision, intent(in), dimension(:,:) :: gbuf
668         double precision, intent(out), dimension(:,:) :: lbuf
669         integer, optional, intent(in) :: mroot, idebug
670     
671     #ifdef MPI
672         integer :: i,j,lroot, lidebug
673     
674         if (.not. present(mroot)) then
675            lroot = 0
676         else
677            lroot = mroot
678         endif
679     
680         if (.not. present(idebug)) then
681            lidebug = 0
682         else
683            lidebug = idebug
684         endif
685     
686         if(myPE.eq.lroot) then
687            call assert( size(lbuf,2).eq.size(gbuf,2),  &
688                 '** scatter_2d: size(lbuf,2).ne.size(gbuf,2) ', &
689                 size(lbuf,2), size(gbuf,2) )
690         endif
691     
692         do j=lbound(lbuf,2),ubound(lbuf,2)
693            call scatter_1d( lbuf(:,j), gbuf(:,j), lroot, lidebug )
694         enddo
695     #else
696         lbuf = gbuf
697     #endif
698     
699         return
700       end subroutine scatter_2d
701     
702       subroutine scatter_3d( lbuf, gbuf, mroot, idebug )
703         double precision, intent(in), dimension(:,:,:) :: gbuf
704         double precision, intent(out), dimension(:,:,:) :: lbuf
705         integer, optional, intent(in) :: mroot, idebug
706     
707     #ifdef MPI
708         integer :: j,k,lroot, lidebug
709     
710         if (.not. present(mroot)) then
711            lroot = 0
712         else
713            lroot = mroot
714         endif
715     
716         if (.not. present(idebug)) then
717            lidebug = 0
718         else
719            lidebug = idebug
720         endif
721     
722         if(myPE.eq.lroot) then
723            call assert( size(lbuf,2).eq.size(gbuf,2),  &
724                 '** scatter_3d: size(lbuf,2).ne.size(gbuf,2) ', &
725                 size(lbuf,2), size(gbuf,2) )
726     
727            call assert( size(lbuf,3).eq.size(gbuf,3),  &
728                 '** scatter_3d: size(lbuf,3).ne.size(gbuf,3) ', &
729                 size(lbuf,3), size(gbuf,3) )
730         endif
731     
732         do k=lbound(lbuf,3),ubound(lbuf,3)
733            do j=lbound(lbuf,2),ubound(lbuf,2)
734               call scatter_1d( lbuf(:,j,k), gbuf(:,j,k), lroot, lidebug )
735            enddo
736         enddo
737     #else
738         lbuf = gbuf
739     #endif
740     
741         return
742       end subroutine scatter_3d
743     
744       subroutine scatter_1c( lbuf, gbuf, mroot, idebug )
745     
746         use functions
747         implicit none
748     
749         character(len=*), intent(in), dimension(:) :: gbuf
750         character(len=*), intent(out), dimension(:) :: lbuf
751         integer, optional, intent(in) :: mroot, idebug
752     
753     #ifdef MPI
754         integer, allocatable, dimension(:,:) :: gbuf_pack,lbuf1
755         character(len=80) :: string
756     
757         integer :: sendtype, recvtype, ijk1, ijk2, recvcnt, ierr,lroot, lidebug
758         integer :: i,j,k,ibuffer,iproc, ioffset
759         integer :: ijk
760         integer :: lenchar, icount
761     
762         !       check to see whether there is root
763     
764         if (.not. present(mroot)) then
765            lroot = 0
766         else
767            lroot = mroot
768         endif
769     
770         if (.not. present(idebug)) then
771            lidebug = 0
772         else
773            lidebug = idebug
774         endif
775     
776         ijk1 = ijkstart3
777         ijk2 = ijkend3
778     
779         lenchar = len(gbuf(1))
780     
781         if(myPE.eq.lroot) then
782            allocate(gbuf_pack(ijkmax3,lenchar))
783         else
784            allocate(gbuf_pack(10,lenchar))
785         endif
786     
787         allocate(lbuf1(ijk1:ijk2,lenchar))
788     
789         if(myPE.eq.lroot) then
790            do i = 1,ijkmax3
791               do j = 1,lenchar
792     
793                  string = gbuf(i)(1:lenchar)
794                  gbuf_pack(i,j) = ichar(string(j:j))
795     
796               enddo
797            enddo
798         endif
799     
800         call scatter_2i(lbuf1,gbuf_pack)
801     
802         do i = ijk1, ijk2
803            do j = 1,lenchar
804     
805               lbuf(i)(j:j) = char(lbuf1(i,j))
806     
807            enddo
808         enddo
809     
810         deallocate(gbuf_pack)
811         deallocate(lbuf1)
812     #else
813         lbuf = gbuf
814     #endif
815     
816         return
817       end subroutine scatter_1c
818     
819     
820       subroutine scatter_1l( lbuf, gbuf, mroot, idebug )
821     
822         use functions
823         implicit none
824     
825         logical, intent(in), dimension(:) :: gbuf
826         logical, intent(out), dimension(:) :: lbuf
827         integer, optional, intent(in) :: mroot, idebug
828     
829     #ifdef MPI
830         logical, allocatable, dimension(:) :: gbuf_pack
831     
832         integer :: sendtype, recvtype, ijk1, ijk2, recvcnt, ierr,lroot, lidebug
833         integer :: i,j,k,ibuffer,iproc, ioffset
834         integer :: ijk
835     
836         !       check to see whether there is root
837     
838         if (.not. present(mroot)) then
839            lroot = 0
840         else
841            lroot = mroot
842         endif
843     
844         if (.not. present(idebug)) then
845            lidebug = 0
846         else
847            lidebug = idebug
848         endif
849     
850         if(myPE.eq.lroot) then
851            allocate(gbuf_pack(sum(ijksize3_all(:))))
852         else
853            allocate(gbuf_pack(10))
854         endif
855     
856         if( myPE.eq.lroot) then
857            ioffset = 0
858            do iproc = 0,numPEs-1
859               ibuffer = 0
860               do k = kstart3_all(iproc), kend3_all(iproc)
861                  do j = jstart3_all(iproc), jend3_all(iproc)
862                     do i = istart3_all(iproc), iend3_all(iproc)
863     
864                        ibuffer = funijk_proc(i,j,k,iproc) + ioffset
865                        gbuf_pack(ibuffer) = gbuf(funijk_gl(i,j,k))
866     
867                     enddo
868                  enddo
869               enddo
870               ioffset = ibuffer
871            enddo
872         endif
873     
874         sendtype = MPI_LOGICAL
875         recvtype = sendtype
876     
877         ijk1 = ijkstart3
878         ijk2 = ijkend3
879     
880         recvcnt = ijk2-ijk1+1
881     
882         !       Call MPI routines
883     
884         call MPI_Scatterv( gbuf_pack, ijksize3_all, displs, sendtype, &
885              lbuf, recvcnt, recvtype,  &
886              lroot, MPI_COMM_WORLD, ierr )
887         call MPI_Check( 'scatter_1l:MPI_Scatterv', ierr )
888     
889         deallocate(gbuf_pack)
890     #else
891         lbuf = gbuf
892     #endif
893     
894         return
895       end subroutine scatter_1l
896     
897     
898       !       Routines to gather lbuf from individual processors and put it on
899       !       processor root in gbuf
900       !       Logic is similar to the scatter routines above.
901     
902       subroutine gather_1i( lbuf, gbuf, mroot, idebug )
903     
904         use functions
905         implicit none
906     
907         integer, intent(in), dimension(:) :: lbuf
908         integer, intent(out), dimension(:) :: gbuf
909         integer, optional, intent(in) :: mroot, idebug
910     
911     #ifdef MPI
912         integer, allocatable, dimension(:) :: gbuf_pack
913     
914         integer :: recvtype, sendtype, ijk1,ijk2,sendcnt, ierr,lroot, lidebug
915         integer :: i,j,k,ibuffer,iproc, ioffset
916         integer :: ijk, ijk_gl
917         integer :: istartl, iendl, jstartl, jendl, kstartl, kendl
918         logical :: isok_k,isok_j,isok_i, isinterior
919         logical :: isbc_k,isbc_j,isbc_i, isboundary, need_copy
920     
921         !       check to see whether there is root
922     
923         if (.not. present(mroot)) then
924            lroot = 0
925         else
926            lroot = mroot
927         endif
928     
929         if (.not. present(idebug)) then
930            lidebug = 0
931         else
932            lidebug = idebug
933         endif
934     
935         if(myPE.eq.lroot) then
936            allocate(gbuf_pack(sum(ijksize3_all(:))))
937         else
938            allocate(gbuf_pack(10))
939         endif
940     
941         recvtype = MPI_INTEGER
942         sendtype = recvtype
943     
944         ijk1 = ijkstart3
945         !        ijk2 = ijkend3
946         ijk2 = max(ijkend3,BACKGROUND_IJKEND3)   !  For cell re-indexing
947     
948         sendcnt = ijk2-ijk1+1
949     
950         call MPI_Gatherv( lbuf, sendcnt, sendtype,  &
951              gbuf_pack, ijksize3_all, displs, recvtype, &
952              lroot, MPI_COMM_WORLD, ierr )
953         call MPI_Check( 'gather_1i:MPI_Gatherv', ierr )
954     
955         if( myPE.eq.lroot) then
956            ioffset = 0
957            do iproc = 0,numPEs-1
958               ibuffer = 0
959               istartl = istart1_all(iproc)
960               iendl = iend1_all(iproc)
961               jstartl = jstart1_all(iproc)
962               jendl = jend1_all(iproc)
963               kstartl = kstart1_all(iproc)
964               kendl = kend1_all(iproc)
965     
966               if(istart3_all(iproc).eq.imin3) istartl = istart3_all(iproc)
967               if(iend3_all(iproc).eq.imax3) iendl = iend3_all(iproc)
968               if(jstart3_all(iproc).eq.jmin3) jstartl = jstart3_all(iproc)
969               if(jend3_all(iproc).eq.jmax3) jendl = jend3_all(iproc)
970               if(kstart3_all(iproc).eq.kmin3) kstartl = kstart3_all(iproc)
971               if(kend3_all(iproc).eq.kmax3) kendl = kend3_all(iproc)
972     
973               do k = kstart3_all(iproc), kend3_all(iproc)
974                  do j = jstart3_all(iproc), jend3_all(iproc)
975                     do i = istart3_all(iproc), iend3_all(iproc)
976     
977                        ibuffer = funijk_proc(i,j,k,iproc) + ioffset
978                        isok_k = (kstartl <= k) .and. (k <=kendl)
979                        isok_j = (jstartl <= j) .and. (j <=jendl)
980                        isok_i = (istartl <= i) .and. (i <=iendl)
981     
982                        need_copy = isok_k .and. isok_j .and. isok_i
983     
984                        if (need_copy) then
985                           ijk_gl = funijk_gl(i,j,k)
986                           gbuf( ijk_gl ) = gbuf_pack(ibuffer)
987                        endif
988     
989                     enddo
990                  enddo
991               enddo
992               ioffset = ibuffer
993            enddo
994         endif
995     
996         deallocate(gbuf_pack)
997     #else
998         gbuf = lbuf
999     #endif
1000     
1001         return
1002       end subroutine gather_1i
1003     
1004     
1005       subroutine gather_2i( lbuf, gbuf, mroot, idebug )
1006         integer, intent(in), dimension(:,:) :: lbuf
1007         integer, intent(out), dimension(:,:) :: gbuf
1008         integer, optional, intent(in) :: mroot, idebug
1009     
1010     #ifdef MPI
1011         integer :: i,j,lroot, lidebug
1012     
1013         if (.not. present(mroot)) then
1014            lroot = 0
1015         else
1016            lroot = mroot
1017         endif
1018     
1019         if (.not. present(idebug)) then
1020            lidebug = 0
1021         else
1022            lidebug = idebug
1023         endif
1024     
1025         if(myPE.eq.lroot) then
1026            call assert( size(lbuf,2).eq.size(gbuf,2),  &
1027                 '** gather_2i: size(lbuf,2).ne.size(gbuf,2) ', &
1028                 size(lbuf,2), size(gbuf,2) )
1029         endif
1030     
1031         do j=lbound(lbuf,2),ubound(lbuf,2)
1032            call gather_1i( lbuf(:,j), gbuf(:,j), lroot, lidebug )
1033         enddo
1034     #else
1035         gbuf = lbuf
1036     #endif
1037     
1038         return
1039       end subroutine gather_2i
1040     
1041       subroutine gather_3i( lbuf, gbuf, mroot, idebug )
1042         integer, intent(in), dimension(:,:,:) :: lbuf
1043         integer, intent(out), dimension(:,:,:) :: gbuf
1044         integer, optional, intent(in) :: mroot, idebug
1045     
1046     #ifdef MPI
1047         integer :: j,k,lroot, lidebug
1048     
1049         if (.not. present(mroot)) then
1050            lroot = 0
1051         else
1052            lroot = mroot
1053         endif
1054     
1055         if (.not. present(idebug)) then
1056            lidebug = 0
1057         else
1058            lidebug = idebug
1059         endif
1060     
1061         if(myPE.eq.lroot) then
1062            call assert( size(lbuf,2).eq.size(gbuf,2),  &
1063                 '** gather_3i: size(lbuf,2).ne.size(gbuf,2) ', &
1064                 size(lbuf,2), size(gbuf,2) )
1065     
1066            call assert( size(lbuf,3).eq.size(gbuf,3),  &
1067                 '** gather_3i: size(lbuf,3).ne.size(gbuf,3) ', &
1068                 size(lbuf,3), size(gbuf,3) )
1069         endif
1070     
1071         do k=lbound(lbuf,3),ubound(lbuf,3)
1072            do j=lbound(lbuf,2),ubound(lbuf,2)
1073               call gather_1i( lbuf(:,j,k), gbuf(:,j,k), lroot, lidebug )
1074            enddo
1075         enddo
1076     #else
1077         gbuf = lbuf
1078     #endif
1079     
1080         return
1081       end subroutine gather_3i
1082     
1083       subroutine gather_1r( lbuf, gbuf, mroot, idebug )
1084     
1085         use functions
1086         implicit none
1087     
1088         real, intent(in), dimension(:) :: lbuf
1089         real, intent(out), dimension(:) :: gbuf
1090         integer, optional, intent(in) :: mroot, idebug
1091     
1092     #ifdef MPI
1093         real, allocatable, dimension(:) :: gbuf_pack
1094     
1095         integer :: recvtype, sendtype, ijk1,ijk2,sendcnt, ierr,lroot, lidebug
1096         integer :: i,j,k,ibuffer,iproc, ioffset
1097         integer :: ijk, ijk_gl
1098         integer :: istartl, iendl, jstartl, jendl, kstartl, kendl
1099         logical :: isok_k,isok_j,isok_i, isinterior
1100         logical :: isbc_k,isbc_j,isbc_i, isboundary, need_copy
1101     
1102         if (.not. present(mroot)) then
1103            lroot = 0
1104         else
1105            lroot = mroot
1106         endif
1107     
1108         if (.not. present(idebug)) then
1109            lidebug = 0
1110         else
1111            lidebug = idebug
1112         endif
1113     
1114         if(myPE.eq.lroot) then
1115            allocate(gbuf_pack(sum(ijksize3_all(:))))
1116         else
1117            allocate(gbuf_pack(10))
1118         endif
1119     
1120         recvtype = MPI_REAL
1121         sendtype = recvtype
1122     
1123         ijk1 = ijkstart3
1124         !        ijk2 = ijkend3
1125         ijk2 = max(ijkend3,BACKGROUND_IJKEND3)   !  For cell re-indexing
1126     
1127         sendcnt = ijk2-ijk1+1
1128     
1129         call MPI_Gatherv( lbuf, sendcnt, sendtype,  &
1130              gbuf_pack, ijksize3_all, displs, recvtype, &
1131              lroot, MPI_COMM_WORLD, ierr )
1132         call MPI_Check( 'gather_1r:MPI_Gatherv', ierr )
1133     
1134         if( myPE.eq.lroot) then
1135            ioffset = 0
1136            do iproc = 0,numPEs-1
1137               ibuffer = 0
1138               istartl = istart1_all(iproc)
1139               iendl = iend1_all(iproc)
1140               jstartl = jstart1_all(iproc)
1141               jendl = jend1_all(iproc)
1142               kstartl = kstart1_all(iproc)
1143               kendl = kend1_all(iproc)
1144     
1145               if(istart3_all(iproc).eq.imin3) istartl = istart3_all(iproc)
1146               if(iend3_all(iproc).eq.imax3) iendl = iend3_all(iproc)
1147               if(jstart3_all(iproc).eq.jmin3) jstartl = jstart3_all(iproc)
1148               if(jend3_all(iproc).eq.jmax3) jendl = jend3_all(iproc)
1149               if(kstart3_all(iproc).eq.kmin3) kstartl = kstart3_all(iproc)
1150               if(kend3_all(iproc).eq.kmax3) kendl = kend3_all(iproc)
1151     
1152               do k = kstart3_all(iproc), kend3_all(iproc)
1153                  do j = jstart3_all(iproc), jend3_all(iproc)
1154                     do i = istart3_all(iproc), iend3_all(iproc)
1155     
1156                        ibuffer = funijk_proc(i,j,k,iproc) + ioffset
1157                        isok_k = (kstartl <= k) .and. (k <=kendl)
1158                        isok_j = (jstartl <= j) .and. (j <=jendl)
1159                        isok_i = (istartl <= i) .and. (i <=iendl)
1160     
1161                        need_copy = isok_k .and. isok_j .and. isok_i
1162     
1163                        if (need_copy) then
1164                           ijk_gl = funijk_gl(i,j,k)
1165                           gbuf( ijk_gl ) = gbuf_pack(ibuffer)
1166                        endif
1167     
1168                     enddo
1169                  enddo
1170               enddo
1171               ioffset = ibuffer
1172            enddo
1173         endif
1174     
1175         deallocate(gbuf_pack)
1176     #else
1177         gbuf = lbuf
1178     #endif
1179     
1180         return
1181       end subroutine gather_1r
1182     
1183       subroutine gather_2r( lbuf, gbuf, mroot, idebug )
1184         real, intent(in), dimension(:,:) :: lbuf
1185         real, intent(out), dimension(:,:) :: gbuf
1186         integer, optional, intent(in) :: mroot, idebug
1187     
1188     #ifdef MPI
1189         integer :: i,j,lroot, lidebug
1190     
1191         if (.not. present(mroot)) then
1192            lroot = 0
1193         else
1194            lroot = mroot
1195         endif
1196     
1197         if (.not. present(idebug)) then
1198            lidebug = 0
1199         else
1200            lidebug = idebug
1201         endif
1202     
1203         if(myPE.eq.lroot) then
1204            call assert( size(lbuf,2).eq.size(gbuf,2),  &
1205                 '** gather_2r: size(lbuf,2).ne.size(gbuf,2) ', &
1206                 size(lbuf,2), size(gbuf,2) )
1207         endif
1208     
1209         do j=lbound(lbuf,2),ubound(lbuf,2)
1210            call gather_1r( lbuf(:,j), gbuf(:,j), lroot, lidebug )
1211         enddo
1212     #else
1213         gbuf = lbuf
1214     #endif
1215     
1216         return
1217       end subroutine gather_2r
1218     
1219       subroutine gather_3r( lbuf, gbuf, mroot, idebug )
1220         real, intent(in), dimension(:,:,:) :: lbuf
1221         real, intent(out), dimension(:,:,:) :: gbuf
1222         integer, optional, intent(in) :: mroot, idebug
1223     
1224     #ifdef MPI
1225         integer :: j,k,lroot, lidebug
1226     
1227         if (.not. present(mroot)) then
1228            lroot = 0
1229         else
1230            lroot = mroot
1231         endif
1232     
1233         if (.not. present(idebug)) then
1234            lidebug = 0
1235         else
1236            lidebug = idebug
1237         endif
1238     
1239         if(myPE.eq.lroot) then
1240            call assert( size(lbuf,2).eq.size(gbuf,2),  &
1241                 '** gather_3r: size(lbuf,2).ne.size(gbuf,2) ', &
1242                 size(lbuf,2), size(gbuf,2) )
1243     
1244            call assert( size(lbuf,3).eq.size(gbuf,3),  &
1245                 '** gather_3r: size(lbuf,3).ne.size(gbuf,3) ', &
1246                 size(lbuf,3), size(gbuf,3) )
1247         endif
1248     
1249         do k=lbound(lbuf,3),ubound(lbuf,3)
1250            do j=lbound(lbuf,2),ubound(lbuf,2)
1251               call gather_1r( lbuf(:,j,k), gbuf(:,j,k), lroot, lidebug )
1252            enddo
1253         enddo
1254     #else
1255         gbuf = lbuf
1256     #endif
1257     
1258         return
1259       end subroutine gather_3r
1260     
1261     
1262       subroutine gather_1d( lbuf, gbuf, mroot, idebug )
1263     
1264         use functions
1265         implicit none
1266     
1267         double precision, intent(in), dimension(:) :: lbuf
1268         double precision, intent(out), dimension(:) :: gbuf
1269         integer, optional, intent(in) :: mroot, idebug
1270     
1271     #ifdef MPI
1272         double precision, allocatable, dimension(:) :: gbuf_pack
1273     
1274         integer :: recvtype, sendtype, ijk1,ijk2,sendcnt, ierr,lroot, lidebug
1275         integer :: i,j,k,ibuffer,iproc, ioffset
1276         integer :: ijk, ijk_gl
1277         logical :: isok_k,isok_j,isok_i, isinterior
1278         logical :: isbc_k,isbc_j,isbc_i, isboundary, need_copy
1279         integer :: istartl, iendl, jstartl, jendl, kstartl, kendl
1280     
1281         if (.not. present(mroot)) then
1282            lroot = 0
1283         else
1284            lroot = mroot
1285         endif
1286     
1287         if (.not. present(idebug)) then
1288            lidebug = 0
1289         else
1290            lidebug = idebug
1291         endif
1292     
1293         if(myPE.eq.lroot) then
1294            allocate(gbuf_pack(sum(ijksize3_all(:))))
1295         else
1296            allocate(gbuf_pack(10))
1297         endif
1298     
1299         recvtype = MPI_DOUBLE_PRECISION
1300         sendtype = recvtype
1301     
1302         ijk1 = ijkstart3
1303         !        ijk2 = ijkend3
1304         ijk2 = max(ijkend3,BACKGROUND_IJKEND3)   !  For cell re-indexing
1305     
1306         sendcnt = ijk2-ijk1+1
1307     
1308         call MPI_Gatherv( lbuf, sendcnt, sendtype,  &
1309              gbuf_pack, ijksize3_all, displs, recvtype, &
1310              lroot, MPI_COMM_WORLD, ierr )
1311         call MPI_Check( 'gather_1d:MPI_Gatherv', ierr )
1312     
1313         if( myPE.eq.lroot) then
1314            ioffset = 0
1315            do iproc = 0,numPEs-1
1316               ibuffer = 0
1317               istartl = istart1_all(iproc)
1318               iendl = iend1_all(iproc)
1319               jstartl = jstart1_all(iproc)
1320               jendl = jend1_all(iproc)
1321               kstartl = kstart1_all(iproc)
1322               kendl = kend1_all(iproc)
1323     
1324               if(istart3_all(iproc).eq.imin3) istartl = istart3_all(iproc)
1325               if(iend3_all(iproc).eq.imax3) iendl = iend3_all(iproc)
1326               if(jstart3_all(iproc).eq.jmin3) jstartl = jstart3_all(iproc)
1327               if(jend3_all(iproc).eq.jmax3) jendl = jend3_all(iproc)
1328               if(kstart3_all(iproc).eq.kmin3) kstartl = kstart3_all(iproc)
1329               if(kend3_all(iproc).eq.kmax3) kendl = kend3_all(iproc)
1330     
1331               do k = kstart3_all(iproc), kend3_all(iproc)
1332                  do j = jstart3_all(iproc), jend3_all(iproc)
1333                     do i = istart3_all(iproc), iend3_all(iproc)
1334     
1335                        ibuffer = funijk_proc(i,j,k,iproc) + ioffset
1336                        isok_k = (kstartl <= k) .and. (k <=kendl)
1337                        isok_j = (jstartl <= j) .and. (j <=jendl)
1338                        isok_i = (istartl <= i) .and. (i <=iendl)
1339     
1340                        need_copy = isok_k .and. isok_j .and. isok_i
1341     
1342                        if (need_copy) then
1343                           ijk_gl = funijk_gl(i,j,k)
1344                           gbuf( ijk_gl ) = gbuf_pack(ibuffer)
1345                        endif
1346     
1347                     enddo
1348                  enddo
1349               enddo
1350               ioffset = ibuffer
1351            enddo
1352         endif
1353     
1354         deallocate(gbuf_pack)
1355     #else
1356         gbuf = lbuf
1357     #endif
1358     
1359         return
1360       end subroutine gather_1d
1361     
1362       subroutine gather_2d( lbuf, gbuf, mroot, idebug )
1363         double precision, intent(in), dimension(:,:) :: lbuf
1364         double precision, intent(out), dimension(:,:) :: gbuf
1365         integer, optional, intent(in) :: mroot, idebug
1366     
1367     #ifdef MPI
1368         integer :: i,j,lroot, lidebug
1369     
1370         if (.not. present(mroot)) then
1371            lroot = 0
1372         else
1373            lroot = mroot
1374         endif
1375     
1376         if (.not. present(idebug)) then
1377            lidebug = 0
1378         else
1379            lidebug = idebug
1380         endif
1381     
1382         if(myPE.eq.lroot) then
1383            call assert( size(lbuf,2).eq.size(gbuf,2),  &
1384                 '** gather_2d: size(lbuf,2).ne.size(gbuf,2) ', &
1385                 size(lbuf,2), size(gbuf,2) )
1386         endif
1387     
1388         do j=lbound(lbuf,2),ubound(lbuf,2)
1389            call gather_1d( lbuf(:,j), gbuf(:,j), lroot, lidebug )
1390         enddo
1391     #else
1392         gbuf = lbuf
1393     #endif
1394     
1395         return
1396       end subroutine gather_2d
1397     
1398       subroutine gather_3d( lbuf, gbuf, mroot, idebug )
1399         double precision, intent(in), dimension(:,:,:) :: lbuf
1400         double precision, intent(out), dimension(:,:,:) :: gbuf
1401         integer, optional, intent(in) :: mroot, idebug
1402     
1403     #ifdef MPI
1404         integer :: j,k,lroot, lidebug
1405     
1406         if (.not. present(mroot)) then
1407            lroot = 0
1408         else
1409            lroot = mroot
1410         endif
1411     
1412         if (.not. present(idebug)) then
1413            lidebug = 0
1414         else
1415            lidebug = idebug
1416         endif
1417     
1418         if(myPE.eq.lroot) then
1419            call assert( size(lbuf,2).eq.size(gbuf,2),  &
1420                 '** gather_3d: size(lbuf,2).ne.size(gbuf,2) ', &
1421                 size(lbuf,2), size(gbuf,2) )
1422     
1423            call assert( size(lbuf,3).eq.size(gbuf,3),  &
1424                 '** gather_3d: size(lbuf,3).ne.size(gbuf,3) ', &
1425                 size(lbuf,3), size(gbuf,3) )
1426         endif
1427     
1428         do k=lbound(lbuf,3),ubound(lbuf,3)
1429            do j=lbound(lbuf,2),ubound(lbuf,2)
1430               call gather_1d( lbuf(:,j,k), gbuf(:,j,k), lroot, lidebug )
1431            enddo
1432         enddo
1433     #else
1434         gbuf = lbuf
1435     #endif
1436     
1437         return
1438       end subroutine gather_3d
1439     
1440     
1441       subroutine gather_1c( lbuf, gbuf, mroot, idebug )
1442     
1443         use functions
1444         implicit none
1445     
1446         character(len=*), intent(in), dimension(:) :: lbuf
1447         character(len=*), intent(out), dimension(:) :: gbuf
1448         integer, optional, intent(in) :: mroot, idebug
1449     
1450     #ifdef MPI
1451         integer, allocatable, dimension(:,:) :: gbuf_pack,lbuf1
1452         character(len=80) :: string
1453     
1454         integer :: recvtype, sendtype, ijk1,ijk2,sendcnt, ierr,lroot, lidebug
1455         integer :: i,j,k,ibuffer,iproc, ioffset
1456         integer :: ijk, ijk_gl
1457         integer :: istartl, iendl, jstartl, jendl, kstartl, kendl
1458         integer :: lenchar, icount
1459         logical :: isok_k,isok_j,isok_i, isinterior
1460         logical :: isbc_k,isbc_j,isbc_i, isboundary, need_copy
1461     
1462         !       check to see whether there is root
1463     
1464         if (.not. present(mroot)) then
1465            lroot = 0
1466         else
1467            lroot = mroot
1468         endif
1469     
1470         if (.not. present(idebug)) then
1471            lidebug = 0
1472         else
1473            lidebug = idebug
1474         endif
1475     
1476     
1477         ijk1 = ijkstart3
1478         !        ijk2 = ijkend3
1479         ijk2 = max(ijkend3,BACKGROUND_IJKEND3)   !  For cell re-indexing
1480     
1481         lenchar = len(lbuf(1))
1482     
1483         if(myPE.eq.lroot) then
1484            allocate(gbuf_pack(ijkmax3,lenchar))
1485         else
1486            allocate(gbuf_pack(10,lenchar))
1487         endif
1488     
1489         allocate(lbuf1(ijk1:ijk2,lenchar))
1490     
1491         do i = ijk1,ijk2
1492            string = lbuf(i)(1:lenchar)
1493            do j = 1,lenchar
1494               lbuf1(i,j) = ichar(string(j:j))
1495            enddo
1496         enddo
1497     
1498         call gather_2i(lbuf1, gbuf_pack)
1499     
1500         if(myPE.eq.lroot) then
1501            do i = 1,ijkmax3
1502               do j = 1,lenchar
1503     
1504                  string(j:j) = char(gbuf_pack(i,j))
1505     
1506               enddo
1507               gbuf(i)(1:lenchar) = string(1:lenchar)
1508     
1509            enddo
1510         endif
1511     
1512         deallocate(gbuf_pack)
1513         deallocate(lbuf1)
1514     #else
1515         gbuf = lbuf
1516     #endif
1517     
1518         return
1519       end subroutine gather_1c
1520     
1521     
1522       subroutine gather_1l( lbuf, gbuf, mroot, idebug )
1523     
1524         use functions
1525         implicit none
1526     
1527         logical, intent(in), dimension(:) :: lbuf
1528         logical, intent(out), dimension(:) :: gbuf
1529         integer, optional, intent(in) :: mroot, idebug
1530     
1531     #ifdef MPI
1532         logical, allocatable, dimension(:) :: gbuf_pack
1533     
1534         integer :: recvtype, sendtype, ijk1,ijk2,sendcnt, ierr,lroot, lidebug
1535         integer :: i,j,k,ibuffer,iproc, ioffset
1536         integer :: ijk, ijk_gl
1537         integer :: istartl, iendl, jstartl, jendl, kstartl, kendl
1538         logical :: isok_k,isok_j,isok_i, isinterior
1539         logical :: isbc_k,isbc_j,isbc_i, isboundary, need_copy
1540     
1541         !       check to see whether there is root
1542     
1543         if (.not. present(mroot)) then
1544            lroot = 0
1545         else
1546            lroot = mroot
1547         endif
1548     
1549         if (.not. present(idebug)) then
1550            lidebug = 0
1551         else
1552            lidebug = idebug
1553         endif
1554     
1555         if(myPE.eq.lroot) then
1556            allocate(gbuf_pack(sum(ijksize3_all(:))))
1557         else
1558            allocate(gbuf_pack(10))
1559         endif
1560     
1561         recvtype = MPI_LOGICAL
1562         sendtype = recvtype
1563     
1564         ijk1 = ijkstart3
1565         !        ijk2 = ijkend3
1566         ijk2 = max(ijkend3,BACKGROUND_IJKEND3)   !  For cell re-indexing
1567     
1568         sendcnt = ijk2-ijk1+1
1569     
1570         call MPI_Gatherv( lbuf, sendcnt, sendtype,  &
1571              gbuf_pack, ijksize3_all, displs, recvtype, &
1572              lroot, MPI_COMM_WORLD, ierr )
1573         call MPI_Check( 'gather_1l:MPI_Gatherv', ierr )
1574     
1575         if( myPE.eq.lroot) then
1576            ioffset = 0
1577            do iproc = 0,numPEs-1
1578               ibuffer = 0
1579               istartl = istart1_all(iproc)
1580               iendl = iend1_all(iproc)
1581               jstartl = jstart1_all(iproc)
1582               jendl = jend1_all(iproc)
1583               kstartl = kstart1_all(iproc)
1584               kendl = kend1_all(iproc)
1585     
1586               if(istart3_all(iproc).eq.imin3) istartl = istart3_all(iproc)
1587               if(iend3_all(iproc).eq.imax3) iendl = iend3_all(iproc)
1588               if(jstart3_all(iproc).eq.jmin3) jstartl = jstart3_all(iproc)
1589               if(jend3_all(iproc).eq.jmax3) jendl = jend3_all(iproc)
1590               if(kstart3_all(iproc).eq.kmin3) kstartl = kstart3_all(iproc)
1591               if(kend3_all(iproc).eq.kmax3) kendl = kend3_all(iproc)
1592     
1593               do k = kstart3_all(iproc), kend3_all(iproc)
1594                  do j = jstart3_all(iproc), jend3_all(iproc)
1595                     do i = istart3_all(iproc), iend3_all(iproc)
1596     
1597                        ibuffer = funijk_proc(i,j,k,iproc) + ioffset
1598                        isok_k = (kstartl <= k) .and. (k <=kendl)
1599                        isok_j = (jstartl <= j) .and. (j <=jendl)
1600                        isok_i = (istartl <= i) .and. (i <=iendl)
1601     
1602                        need_copy = isok_k .and. isok_j .and. isok_i
1603     
1604                        if (need_copy) then
1605                           ijk_gl = funijk_gl(i,j,k)
1606                           gbuf( ijk_gl ) = gbuf_pack(ibuffer)
1607                        endif
1608     
1609                     enddo
1610                  enddo
1611               enddo
1612               ioffset = ibuffer
1613            enddo
1614         endif
1615     
1616         deallocate(gbuf_pack)
1617     #else
1618         gbuf = lbuf
1619     #endif
1620     
1621         return
1622       end subroutine gather_1l
1623     
1624     
1625       !       Routines to broadcast information from processor 0 in buffer to all
1626       !       the processors
1627     
1628       subroutine bcast_0i( buffer, mroot, idebug )
1629         integer, intent(inout) :: buffer
1630         integer, optional, intent(in) :: mroot, idebug
1631     
1632     #ifdef MPI
1633         integer :: datatype, count, ierr,lroot, lidebug
1634     
1635         if (.not. present(mroot)) then
1636            lroot = 0
1637         else
1638            lroot = mroot
1639         endif
1640     
1641         if (.not. present(idebug)) then
1642            lidebug = 0
1643         else
1644            lidebug = idebug
1645         endif
1646     
1647         datatype = MPI_INTEGER
1648     
1649         count = 1
1650     
1651         call MPI_Bcast( buffer, count, datatype, lroot, MPI_COMM_WORLD, ierr)
1652         call MPI_Check( 'bcast_0i:MPI_Bcast', ierr )
1653     #endif
1654     
1655         return
1656       end subroutine bcast_0i
1657     
1658     
1659       subroutine bcast_1i( buffer, mroot, idebug )
1660         integer, intent(inout), dimension(:) :: buffer
1661         integer, optional, intent(in) :: mroot, idebug
1662     
1663     #ifdef MPI
1664         integer :: datatype, count, ierr,lroot, lidebug
1665     
1666         if (.not. present(mroot)) then
1667            lroot = 0
1668         else
1669            lroot = mroot
1670         endif
1671     
1672         if (.not. present(idebug)) then
1673            lidebug = 0
1674         else
1675            lidebug = idebug
1676         endif
1677     
1678         datatype = MPI_INTEGER
1679     
1680         count = size(buffer,1)
1681     
1682         call MPI_Bcast( buffer, count, datatype, lroot, MPI_COMM_WORLD, ierr)
1683         call MPI_Check( 'bcast_1i:MPI_Bcast', ierr )
1684     #endif
1685     
1686         return
1687       end subroutine bcast_1i
1688     
1689     
1690       subroutine bcast_2i( buffer, mroot, idebug )
1691         integer, intent(inout), dimension(:,:) :: buffer
1692         integer, optional, intent(in) :: mroot, idebug
1693     
1694     #ifdef MPI
1695         integer :: i,j,lroot, lidebug
1696     
1697         if (.not. present(mroot)) then
1698            lroot = 0
1699         else
1700            lroot = mroot
1701         endif
1702     
1703         if (.not. present(idebug)) then
1704            lidebug = 0
1705         else
1706            lidebug = idebug
1707         endif
1708     
1709         do j=lbound(buffer,2),ubound(buffer,2)
1710            call bcast_1i( buffer(:,j), lroot, lidebug )
1711         enddo
1712     #endif
1713     
1714         return
1715       end subroutine bcast_2i
1716     
1717       subroutine bcast_3i( buffer, mroot, idebug )
1718         integer, intent(inout), dimension(:,:,:) :: buffer
1719         integer, optional, intent(in) :: mroot, idebug
1720     
1721     #ifdef MPI
1722         integer :: j,k,lroot, lidebug
1723     
1724         if (.not. present(mroot)) then
1725            lroot = 0
1726         else
1727            lroot = mroot
1728         endif
1729     
1730         if (.not. present(idebug)) then
1731            lidebug = 0
1732         else
1733            lidebug = idebug
1734         endif
1735     
1736         do k=lbound(buffer,3),ubound(buffer,3)
1737            do j=lbound(buffer,2),ubound(buffer,2)
1738               call bcast_1i( buffer(:,j,k), lroot, lidebug )
1739            enddo
1740         enddo
1741     #endif
1742     
1743         return
1744       end subroutine bcast_3i
1745     
1746       subroutine bcast_0r( buffer, mroot, idebug )
1747         real, intent(inout) :: buffer
1748         integer, optional, intent(in) :: mroot, idebug
1749     
1750     #ifdef MPI
1751         integer :: datatype, count, ierr,lroot, lidebug
1752     
1753         if (.not. present(mroot)) then
1754            lroot = 0
1755         else
1756            lroot = mroot
1757         endif
1758     
1759         if (.not. present(idebug)) then
1760            lidebug = 0
1761         else
1762            lidebug = idebug
1763         endif
1764     
1765         datatype = MPI_REAL
1766     
1767         count = 1
1768     
1769         call MPI_Bcast( buffer, count, datatype, lroot, MPI_COMM_WORLD, ierr)
1770         call MPI_Check( 'bcast_0r:MPI_Bcast', ierr )
1771     #endif
1772     
1773         return
1774       end subroutine bcast_0r
1775     
1776       subroutine bcast_1r( buffer, mroot, idebug )
1777         real, intent(inout), dimension(:) :: buffer
1778         integer, optional, intent(in) :: mroot, idebug
1779     
1780     #ifdef MPI
1781         integer :: datatype, count, ierr,lroot, lidebug
1782     
1783         if (.not. present(mroot)) then
1784            lroot = 0
1785         else
1786            lroot = mroot
1787         endif
1788     
1789         if (.not. present(idebug)) then
1790            lidebug = 0
1791         else
1792            lidebug = idebug
1793         endif
1794     
1795         datatype = MPI_REAL
1796     
1797         count = size(buffer,1)
1798     
1799         call MPI_Bcast( buffer, count, datatype, lroot, MPI_COMM_WORLD, ierr)
1800         call MPI_Check( 'bcast_1r:MPI_Bcast', ierr )
1801     #endif
1802     
1803         return
1804       end subroutine bcast_1r
1805     
1806       subroutine bcast_2r( buffer, mroot, idebug )
1807         real, intent(inout), dimension(:,:) :: buffer
1808         integer, optional, intent(in) :: mroot, idebug
1809     
1810     #ifdef MPI
1811         integer :: i,j,lroot, lidebug
1812     
1813         if (.not. present(mroot)) then
1814            lroot = 0
1815         else
1816            lroot = mroot
1817         endif
1818     
1819         if (.not. present(idebug)) then
1820            lidebug = 0
1821         else
1822            lidebug = idebug
1823         endif
1824     
1825         do j=lbound(buffer,2),ubound(buffer,2)
1826            call bcast_1r( buffer(:,j), lroot, lidebug )
1827         enddo
1828     #endif
1829     
1830         return
1831       end subroutine bcast_2r
1832     
1833       subroutine bcast_3r( buffer, mroot, idebug )
1834         real, intent(inout), dimension(:,:,:) :: buffer
1835         integer, optional, intent(in) :: mroot, idebug
1836     
1837     #ifdef MPI
1838         integer :: j,k,lroot, lidebug
1839     
1840         if (.not. present(mroot)) then
1841            lroot = 0
1842         else
1843            lroot = mroot
1844         endif
1845     
1846         if (.not. present(idebug)) then
1847            lidebug = 0
1848         else
1849            lidebug = idebug
1850         endif
1851     
1852         do k=lbound(buffer,3),ubound(buffer,3)
1853            do j=lbound(buffer,2),ubound(buffer,2)
1854               call bcast_1r( buffer(:,j,k), lroot, lidebug )
1855            enddo
1856         enddo
1857     #endif
1858     
1859         return
1860       end subroutine bcast_3r
1861     
1862       subroutine bcast_0d( buffer, mroot, idebug )
1863         double precision, intent(inout) :: buffer
1864         integer, optional, intent(in) :: mroot, idebug
1865     
1866     #ifdef MPI
1867         integer :: datatype, count, ierr,lroot, lidebug
1868     
1869         if (.not. present(mroot)) then
1870            lroot = 0
1871         else
1872            lroot = mroot
1873         endif
1874     
1875         if (.not. present(idebug)) then
1876            lidebug = 0
1877         else
1878            lidebug = idebug
1879         endif
1880     
1881         datatype = MPI_DOUBLE_PRECISION
1882     
1883         count = 1
1884     
1885         call MPI_Bcast( buffer, count, datatype, lroot, MPI_COMM_WORLD, ierr)
1886         call MPI_Check( 'bcast_0d:MPI_Bcast', ierr )
1887     #endif
1888     
1889         return
1890       end subroutine bcast_0d
1891     
1892     
1893       subroutine bcast_1d( buffer, mroot, idebug )
1894         double precision, intent(inout), dimension(:) :: buffer
1895         integer, optional, intent(in) :: mroot, idebug
1896     
1897     #ifdef MPI
1898         integer :: datatype, count, ierr,lroot, lidebug
1899     
1900         if (.not. present(mroot)) then
1901            lroot = 0
1902         else
1903            lroot = mroot
1904         endif
1905     
1906         if (.not. present(idebug)) then
1907            lidebug = 0
1908         else
1909            lidebug = idebug
1910         endif
1911     
1912         datatype = MPI_DOUBLE_PRECISION
1913     
1914         count = size(buffer,1)
1915     
1916         call MPI_Bcast( buffer, count, datatype, lroot, MPI_COMM_WORLD, ierr)
1917         call MPI_Check( 'bcast_1d:MPI_Bcast', ierr )
1918     #endif
1919     
1920         return
1921       end subroutine bcast_1d
1922     
1923       subroutine bcast_2d( buffer, mroot, idebug )
1924         double precision, intent(inout), dimension(:,:) :: buffer
1925         integer, optional, intent(in) :: mroot, idebug
1926     
1927     #ifdef MPI
1928         integer :: i,j,lroot, lidebug
1929     
1930         if (.not. present(mroot)) then
1931            lroot = 0
1932         else
1933            lroot = mroot
1934         endif
1935     
1936         if (.not. present(idebug)) then
1937            lidebug = 0
1938         else
1939            lidebug = idebug
1940         endif
1941     
1942         do j=lbound(buffer,2),ubound(buffer,2)
1943            call bcast_1d( buffer(:,j), lroot, lidebug )
1944         enddo
1945     #endif
1946     
1947         return
1948       end subroutine bcast_2d
1949     
1950       subroutine bcast_3d( buffer, mroot, idebug )
1951         double precision, intent(inout), dimension(:,:,:) :: buffer
1952         integer, optional, intent(in) :: mroot, idebug
1953     
1954     #ifdef MPI
1955         integer :: j,k,lroot, lidebug
1956     
1957         if (.not. present(mroot)) then
1958            lroot = 0
1959         else
1960            lroot = mroot
1961         endif
1962     
1963         if (.not. present(idebug)) then
1964            lidebug = 0
1965         else
1966            lidebug = idebug
1967         endif
1968     
1969         do k=lbound(buffer,3),ubound(buffer,3)
1970            do j=lbound(buffer,2),ubound(buffer,2)
1971               call bcast_1d( buffer(:,j,k), lroot, lidebug )
1972            enddo
1973         enddo
1974     #endif
1975     
1976         return
1977       end subroutine bcast_3d
1978     
1979       subroutine bcast_0c( buffer, mroot, idebug )
1980         character(len=*), intent(inout) :: buffer
1981         integer, optional, intent(in) :: mroot, idebug
1982         character, allocatable, dimension(:) :: buffer1
1983     
1984     #ifdef MPI
1985         integer :: datatype, count, ierr,lroot, lidebug
1986         integer :: lenchar,icount, i, j
1987     
1988         if (.not. present(mroot)) then
1989            lroot = 0
1990         else
1991            lroot = mroot
1992         endif
1993     
1994         if (.not. present(idebug)) then
1995            lidebug = 0
1996         else
1997            lidebug = idebug
1998         endif
1999     
2000         lenchar = len(buffer)
2001     
2002         allocate(buffer1(lenchar))
2003     
2004         icount = 0
2005         do j = 1,lenchar
2006     
2007            icount = icount+1
2008            buffer1(icount) = buffer(j:j)
2009     
2010         enddo
2011     
2012         datatype = MPI_CHARACTER
2013     
2014         count = 1
2015     
2016         call MPI_Bcast( buffer1, count*lenchar, datatype, lroot, MPI_COMM_WORLD, ierr)
2017         call MPI_Check( 'bcast_0c:MPI_Bcast', ierr )
2018     
2019         icount = 0
2020         do j = 1,lenchar
2021     
2022            icount = icount+1
2023            buffer(j:j) = buffer1(icount)
2024     
2025         enddo
2026     
2027         deallocate(buffer1)
2028     #endif
2029     
2030         return
2031       end subroutine bcast_0c
2032     
2033     
2034       subroutine bcast_1c( buffer, mroot, idebug )
2035         character(len=*), intent(inout), dimension(:) :: buffer
2036         integer, optional, intent(in) :: mroot, idebug
2037         character, allocatable, dimension(:) :: buffer1
2038     
2039     #ifdef MPI
2040         integer :: datatype, count, ierr,lroot, lidebug
2041         integer :: lenchar,icount, i, j
2042         character(len=len(buffer(1))) :: string
2043     
2044         if (.not. present(mroot)) then
2045            lroot = 0
2046         else
2047            lroot = mroot
2048         endif
2049     
2050         if (.not. present(idebug)) then
2051            lidebug = 0
2052         else
2053            lidebug = idebug
2054         endif
2055     
2056         lenchar = len(buffer(1))
2057     
2058         allocate(buffer1(size(buffer)*lenchar))
2059     
2060         icount = 0
2061         do i = 1,size(buffer)
2062            string = buffer(i)(1:lenchar)
2063            do j = 1,lenchar
2064     
2065               icount = icount+1
2066               buffer1(icount) = string(j:j)
2067     
2068            enddo
2069         enddo
2070     
2071         datatype = MPI_CHARACTER
2072     
2073         count = size(buffer,1)
2074     
2075         call MPI_Bcast( buffer1, count*lenchar, datatype, lroot, MPI_COMM_WORLD, ierr)
2076         call MPI_Check( 'bcast_1c:MPI_Bcast', ierr )
2077     
2078         icount = 0
2079         do i = 1,size(buffer)
2080            do j = 1,lenchar
2081     
2082               icount = icount+1
2083               string(j:j) = buffer1(icount)
2084     
2085            enddo
2086            buffer(i) = string
2087         enddo
2088     
2089         deallocate(buffer1)
2090     #endif
2091     
2092         return
2093       end subroutine bcast_1c
2094     
2095       subroutine bcast_0l( buffer, mroot, idebug )
2096         logical, intent(inout) :: buffer
2097         integer, optional, intent(in) :: mroot, idebug
2098     
2099     #ifdef MPI
2100         integer :: datatype, count, ierr,lroot, lidebug
2101     
2102         if (.not. present(mroot)) then
2103            lroot = 0
2104         else
2105            lroot = mroot
2106         endif
2107     
2108         if (.not. present(idebug)) then
2109            lidebug = 0
2110         else
2111            lidebug = idebug
2112         endif
2113     
2114         datatype = MPI_LOGICAL
2115     
2116         count = 1
2117     
2118         call MPI_Bcast( buffer, count, datatype, lroot, MPI_COMM_WORLD, ierr)
2119         call MPI_Check( 'bcast_0l:MPI_Bcast', ierr )
2120     #endif
2121     
2122         return
2123       end subroutine bcast_0l
2124     
2125     
2126       subroutine bcast_1l( buffer, mroot, idebug )
2127         logical, intent(inout), dimension(:) :: buffer
2128         integer, optional, intent(in) :: mroot, idebug
2129     
2130     #ifdef MPI
2131         integer :: datatype, count, ierr,lroot, lidebug
2132     
2133         if (.not. present(mroot)) then
2134            lroot = 0
2135         else
2136            lroot = mroot
2137         endif
2138     
2139         if (.not. present(idebug)) then
2140            lidebug = 0
2141         else
2142            lidebug = idebug
2143         endif
2144     
2145         datatype = MPI_LOGICAL
2146     
2147         count = size(buffer,1)
2148     
2149         call MPI_Bcast( buffer, count, datatype, lroot, MPI_COMM_WORLD, ierr)
2150         call MPI_Check( 'bcast_1l:MPI_Bcast', ierr )
2151     #endif
2152     
2153         return
2154       end subroutine bcast_1l
2155     
2156       !       Procedures to do global operations (Sum, Min, Max). _all_ routines
2157       !       send the information to all the processors otherwise they are
2158       !       kept on processor 0.
2159     
2160       subroutine global_sum_0i( lbuf, gbuf, mroot, idebug )
2161         integer, intent(in) :: lbuf
2162         integer, intent(out) :: gbuf
2163         integer, optional, intent(in) :: mroot, idebug
2164     
2165     #ifdef MPI
2166         integer :: recvtype, sendtype, sendcnt, ierr,lroot, lidebug
2167     
2168         if (.not. present(mroot)) then
2169            lroot = 0
2170         else
2171            lroot = mroot
2172         endif
2173     
2174         if (.not. present(idebug)) then
2175            lidebug = 0
2176         else
2177            lidebug = idebug
2178         endif
2179     
2180         recvtype = MPI_INTEGER
2181         sendtype = recvtype
2182     
2183         sendcnt = 1
2184     
2185         call MPI_Reduce( lbuf, gbuf, sendcnt, sendtype,  MPI_SUM, &
2186              lroot, MPI_COMM_WORLD, ierr )
2187         call MPI_Check( 'global_sum_0i:MPI_Reduce', ierr )
2188     #else
2189         gbuf = lbuf
2190     #endif
2191     
2192         return
2193       end subroutine global_sum_0i
2194     
2195     
2196       subroutine global_sum_1i( lbuf, gbuf, mroot, idebug )
2197         integer, intent(in), dimension(:) :: lbuf
2198         integer, intent(out), dimension(:) :: gbuf
2199         integer, optional, intent(in) :: mroot, idebug
2200     
2201     #ifdef MPI
2202         integer :: recvtype, sendtype, sendcnt, ierr,lroot, lidebug
2203     
2204         if (.not. present(mroot)) then
2205            lroot = 0
2206         else
2207            lroot = mroot
2208         endif
2209     
2210         if (.not. present(idebug)) then
2211            lidebug = 0
2212         else
2213            lidebug = idebug
2214         endif
2215     
2216         recvtype = MPI_INTEGER
2217         sendtype = recvtype
2218     
2219         sendcnt = size(lbuf)
2220     
2221         call MPI_Reduce( lbuf, gbuf, sendcnt, sendtype,  MPI_SUM, &
2222              lroot, MPI_COMM_WORLD, ierr )
2223         call MPI_Check( 'global_sum_1i:MPI_Reduce', ierr )
2224     #else
2225         gbuf = lbuf
2226     #endif
2227     
2228         return
2229       end subroutine global_sum_1i
2230     
2231       subroutine global_sum_2i( lbuf, gbuf, mroot, idebug )
2232         integer, intent(in), dimension(:,:) :: lbuf
2233         integer, intent(out), dimension(:,:) :: gbuf
2234         integer, optional, intent(in) :: mroot, idebug
2235     
2236     #ifdef MPI
2237         integer :: i,j,lroot, lidebug
2238     
2239         if (.not. present(mroot)) then
2240            lroot = 0
2241         else
2242            lroot = mroot
2243         endif
2244     
2245         if (.not. present(idebug)) then
2246            lidebug = 0
2247         else
2248            lidebug = idebug
2249         endif
2250     
2251         if(myPE.eq.lroot) then
2252            call assert( size(lbuf,2).eq.size(gbuf,2),  &
2253                 '** global_sum_2i: size(lbuf,2).ne.size(gbuf,2) ', &
2254                 size(lbuf,2), size(gbuf,2) )
2255         endif
2256     
2257         do j=lbound(lbuf,2),ubound(lbuf,2)
2258            call global_sum_1i( lbuf(:,j), gbuf(:,j), lroot, lidebug )
2259         enddo
2260     #else
2261         gbuf = lbuf
2262     #endif
2263     
2264         return
2265       end subroutine global_sum_2i
2266     
2267       subroutine global_sum_3i( lbuf, gbuf, mroot, idebug )
2268         integer, intent(in), dimension(:,:,:) :: lbuf
2269         integer, intent(out), dimension(:,:,:) :: gbuf
2270         integer, optional, intent(in) :: mroot, idebug
2271     
2272     #ifdef MPI
2273         integer :: j,k,lroot, lidebug
2274     
2275         if (.not. present(mroot)) then
2276            lroot = 0
2277         else
2278            lroot = mroot
2279         endif
2280     
2281         if (.not. present(idebug)) then
2282            lidebug = 0
2283         else
2284            lidebug = idebug
2285         endif
2286     
2287         if(myPE.eq.lroot) then
2288            call assert( size(lbuf,2).eq.size(gbuf,2),  &
2289                 '** global_sum_3i: size(lbuf,2).ne.size(gbuf,2) ', &
2290                 size(lbuf,2), size(gbuf,2) )
2291     
2292            call assert( size(lbuf,3).eq.size(gbuf,3),  &
2293                 '** global_sum_3i: size(lbuf,3).ne.size(gbuf,3) ', &
2294                 size(lbuf,3), size(gbuf,3) )
2295         endif
2296     
2297         do k=lbound(lbuf,3),ubound(lbuf,3)
2298            do j=lbound(lbuf,2),ubound(lbuf,2)
2299               call global_sum_1i( lbuf(:,j,k), gbuf(:,j,k), lroot, lidebug )
2300            enddo
2301         enddo
2302     #else
2303         gbuf = lbuf
2304     #endif
2305     
2306         return
2307       end subroutine global_sum_3i
2308     
2309       subroutine global_sum_0r( lbuf, gbuf, mroot, idebug )
2310         real, intent(in) :: lbuf
2311         real, intent(out) :: gbuf
2312         integer, optional, intent(in) :: mroot, idebug
2313     
2314     #ifdef MPI
2315         integer :: recvtype, sendtype, sendcnt, ierr,lroot, lidebug
2316     
2317         if (.not. present(mroot)) then
2318            lroot = 0
2319         else
2320            lroot = mroot
2321         endif
2322     
2323         if (.not. present(idebug)) then
2324            lidebug = 0
2325         else
2326            lidebug = idebug
2327         endif
2328     
2329         recvtype = MPI_REAL
2330         sendtype = recvtype
2331     
2332         sendcnt = 1
2333     
2334         call MPI_Reduce( lbuf, gbuf, sendcnt, sendtype,  MPI_SUM, &
2335              lroot, MPI_COMM_WORLD, ierr )
2336         call MPI_Check( 'global_sum_0r:MPI_Reduce', ierr )
2337     #else
2338         gbuf = lbuf
2339     #endif
2340     
2341         return
2342       end subroutine global_sum_0r
2343     
2344     
2345       subroutine global_sum_1r( lbuf, gbuf, mroot, idebug )
2346         real, intent(in), dimension(:) :: lbuf
2347         real, intent(out), dimension(:) :: gbuf
2348         integer, optional, intent(in) :: mroot, idebug
2349     
2350     #ifdef MPI
2351         integer :: recvtype, sendtype, sendcnt, ierr,lroot, lidebug
2352     
2353         if (.not. present(mroot)) then
2354            lroot = 0
2355         else
2356            lroot = mroot
2357         endif
2358     
2359         if (.not. present(idebug)) then
2360            lidebug = 0
2361         else
2362            lidebug = idebug
2363         endif
2364     
2365         recvtype = MPI_REAL
2366         sendtype = recvtype
2367     
2368         sendcnt = size(lbuf)
2369     
2370         call MPI_Reduce( lbuf, gbuf, sendcnt, sendtype,  MPI_SUM, &
2371              lroot, MPI_COMM_WORLD, ierr )
2372         call MPI_Check( 'global_sum_1r:MPI_Reduce', ierr )
2373     #else
2374         gbuf = lbuf
2375     #endif
2376     
2377         return
2378       end subroutine global_sum_1r
2379     
2380       subroutine global_sum_2r( lbuf, gbuf, mroot, idebug )
2381         real, intent(in), dimension(:,:) :: lbuf
2382         real, intent(out), dimension(:,:) :: gbuf
2383         integer, optional, intent(in) :: mroot, idebug
2384     
2385     #ifdef MPI
2386         integer :: i,j,lroot, lidebug
2387     
2388         if (.not. present(mroot)) then
2389            lroot = 0
2390         else
2391            lroot = mroot
2392         endif
2393     
2394         if (.not. present(idebug)) then
2395            lidebug = 0
2396         else
2397            lidebug = idebug
2398         endif
2399     
2400         if(myPE.eq.lroot) then
2401            call assert( size(lbuf,2).eq.size(gbuf,2),  &
2402                 '** global_sum_2r: size(lbuf,2).ne.size(gbuf,2) ', &
2403                 size(lbuf,2), size(gbuf,2) )
2404         endif
2405     
2406         do j=lbound(lbuf,2),ubound(lbuf,2)
2407            call global_sum_1r( lbuf(:,j), gbuf(:,j), lroot, lidebug )
2408         enddo
2409     #else
2410         gbuf = lbuf
2411     #endif
2412     
2413         return
2414       end subroutine global_sum_2r
2415     
2416       subroutine global_sum_3r( lbuf, gbuf, mroot, idebug )
2417         real, intent(in), dimension(:,:,:) :: lbuf
2418         real, intent(out), dimension(:,:,:) :: gbuf
2419         integer, optional, intent(in) :: mroot, idebug
2420     
2421     #ifdef MPI
2422         integer :: j,k,lroot, lidebug
2423     
2424         if (.not. present(mroot)) then
2425            lroot = 0
2426         else
2427            lroot = mroot
2428         endif
2429     
2430         if (.not. present(idebug)) then
2431            lidebug = 0
2432         else
2433            lidebug = idebug
2434         endif
2435     
2436         if(myPE.eq.lroot) then
2437            call assert( size(lbuf,2).eq.size(gbuf,2),  &
2438                 '** global_sum_3i: size(lbuf,2).ne.size(gbuf,2) ', &
2439                 size(lbuf,2), size(gbuf,2) )
2440     
2441            call assert( size(lbuf,3).eq.size(gbuf,3),  &
2442                 '** global_sum_3i: size(lbuf,3).ne.size(gbuf,3) ', &
2443                 size(lbuf,3), size(gbuf,3) )
2444         endif
2445     
2446         do k=lbound(lbuf,3),ubound(lbuf,3)
2447            do j=lbound(lbuf,2),ubound(lbuf,2)
2448               call global_sum_1r( lbuf(:,j,k), gbuf(:,j,k), lroot, lidebug )
2449            enddo
2450         enddo
2451     #else
2452         gbuf = lbuf
2453     #endif
2454     
2455         return
2456       end subroutine global_sum_3r
2457     
2458       subroutine global_sum_0d( lbuf, gbuf, mroot, idebug )
2459         double precision, intent(in) :: lbuf
2460         double precision, intent(out) :: gbuf
2461         integer, optional, intent(in) :: mroot, idebug
2462     
2463     #ifdef MPI
2464         integer :: recvtype, sendtype, sendcnt, ierr,lroot, lidebug
2465     
2466         if (.not. present(mroot)) then
2467            lroot = 0
2468         else
2469            lroot = mroot
2470         endif
2471     
2472         if (.not. present(idebug)) then
2473            lidebug = 0
2474         else
2475            lidebug = idebug
2476         endif
2477     
2478         recvtype = MPI_DOUBLE_PRECISION
2479         sendtype = recvtype
2480     
2481         sendcnt = 1
2482     
2483         call MPI_Reduce( lbuf, gbuf, sendcnt, sendtype,  MPI_SUM, &
2484              lroot, MPI_COMM_WORLD, ierr )
2485         call MPI_Check( 'global_sum_0d:MPI_Reduce', ierr )
2486     #else
2487         gbuf = lbuf
2488     #endif
2489     
2490         return
2491       end subroutine global_sum_0d
2492     
2493     
2494       subroutine global_sum_1d( lbuf, gbuf, mroot, idebug )
2495         double precision, intent(in), dimension(:) :: lbuf
2496         double precision, intent(out), dimension(:) :: gbuf
2497         integer, optional, intent(in) :: mroot, idebug
2498     
2499     #ifdef MPI
2500         integer :: recvtype, sendtype, sendcnt, ierr,lroot, lidebug
2501     
2502         if (.not. present(mroot)) then
2503            lroot = 0
2504         else
2505            lroot = mroot
2506         endif
2507     
2508         if (.not. present(idebug)) then
2509            lidebug = 0
2510         else
2511            lidebug = idebug
2512         endif
2513     
2514         recvtype = MPI_DOUBLE_PRECISION
2515         sendtype = recvtype
2516     
2517         sendcnt = size(lbuf)
2518     
2519         call MPI_Reduce( lbuf, gbuf, sendcnt, sendtype,  MPI_SUM, &
2520              lroot, MPI_COMM_WORLD, ierr )
2521         call MPI_Check( 'global_sum_1d:MPI_Reduce', ierr )
2522     #else
2523         gbuf = lbuf
2524     #endif
2525     
2526         return
2527       end subroutine global_sum_1d
2528     
2529       subroutine global_sum_2d( lbuf, gbuf, mroot, idebug )
2530         double precision, intent(in), dimension(:,:) :: lbuf
2531         double precision, intent(out), dimension(:,:) :: gbuf
2532         integer, optional, intent(in) :: mroot, idebug
2533     
2534     #ifdef MPI
2535         integer :: i,j,lroot, lidebug
2536     
2537         if (.not. present(mroot)) then
2538            lroot = 0
2539         else
2540            lroot = mroot
2541         endif
2542     
2543         if (.not. present(idebug)) then
2544            lidebug = 0
2545         else
2546            lidebug = idebug
2547         endif
2548     
2549         if(myPE.eq.lroot) then
2550            call assert( size(lbuf,2).eq.size(gbuf,2),  &
2551                 '** global_sum_2d: size(lbuf,2).ne.size(gbuf,2) ', &
2552                 size(lbuf,2), size(gbuf,2) )
2553         endif
2554     
2555         do j=lbound(lbuf,2),ubound(lbuf,2)
2556            call global_sum_1d( lbuf(:,j), gbuf(:,j), lroot, lidebug )
2557         enddo
2558     #else
2559         gbuf = lbuf
2560     #endif
2561     
2562         return
2563       end subroutine global_sum_2d
2564     
2565       subroutine global_sum_3d( lbuf, gbuf, mroot, idebug )
2566         double precision, intent(in), dimension(:,:,:) :: lbuf
2567         double precision, intent(out), dimension(:,:,:) :: gbuf
2568         integer, optional, intent(in) :: mroot, idebug
2569     
2570     #ifdef MPI
2571         integer :: j,k,lroot, lidebug
2572     
2573         if (.not. present(mroot)) then
2574            lroot = 0
2575         else
2576            lroot = mroot
2577         endif
2578     
2579         if (.not. present(idebug)) then
2580            lidebug = 0
2581         else
2582            lidebug = idebug
2583         endif
2584     
2585         if(myPE.eq.lroot) then
2586            call assert( size(lbuf,2).eq.size(gbuf,2),  &
2587                 '** global_sum_3i: size(lbuf,2).ne.size(gbuf,2) ', &
2588                 size(lbuf,2), size(gbuf,2) )
2589     
2590            call assert( size(lbuf,3).eq.size(gbuf,3),  &
2591                 '** global_sum_3i: size(lbuf,3).ne.size(gbuf,3) ', &
2592                 size(lbuf,3), size(gbuf,3) )
2593         endif
2594     
2595         do k=lbound(lbuf,3),ubound(lbuf,3)
2596            do j=lbound(lbuf,2),ubound(lbuf,2)
2597               call global_sum_1d( lbuf(:,j,k), gbuf(:,j,k), lroot, lidebug )
2598            enddo
2599         enddo
2600     #else
2601         gbuf = lbuf
2602     #endif
2603     
2604         return
2605       end subroutine global_sum_3d
2606     
2607       subroutine global_all_sum_onevar_0d( gbuf )
2608         doubleprecision, intent(inout) :: gbuf
2609         doubleprecision :: lbuf
2610     
2611     #ifdef MPI
2612         lbuf = gbuf
2613         call global_all_sum_0d( lbuf, gbuf )
2614     #endif
2615         return
2616       end subroutine global_all_sum_onevar_0d
2617     
2618     
2619       subroutine global_all_sum_onevar_1d( gbuf )
2620         doubleprecision, dimension(:), intent(inout) :: gbuf
2621         doubleprecision, dimension(size(gbuf)) :: lbuf
2622     
2623     #ifdef MPI
2624         lbuf = gbuf
2625         call global_all_sum_1d( lbuf, gbuf )
2626     #endif
2627         return
2628       end subroutine global_all_sum_onevar_1d
2629     
2630       subroutine global_all_sum_onevar_2d( gbuf )
2631         doubleprecision, dimension(:,:), intent(inout) :: gbuf
2632         doubleprecision, dimension(size(gbuf,1),size(gbuf,2)) :: lbuf
2633     
2634     #ifdef MPI
2635         lbuf = gbuf
2636         call global_all_sum_2d( lbuf, gbuf )
2637     #endif
2638         return
2639       end subroutine global_all_sum_onevar_2d
2640     
2641     
2642       subroutine global_all_sum_onevar_3d( gbuf )
2643         doubleprecision, dimension(:,:,:), intent(inout) :: gbuf
2644         doubleprecision, dimension(size(gbuf,1),size(gbuf,2),size(gbuf,3)) :: lbuf
2645     
2646     #ifdef MPI
2647         lbuf = gbuf
2648         call global_all_sum_3d( lbuf, gbuf )
2649     #endif
2650         return
2651       end subroutine global_all_sum_onevar_3d
2652     
2653     
2654     
2655       subroutine global_all_sum_onevar_0i( gbuf )
2656         integer, intent(inout) :: gbuf
2657         integer :: lbuf
2658     
2659     #ifdef MPI
2660         lbuf = gbuf
2661         call global_all_sum_0i( lbuf, gbuf )
2662     #endif
2663         return
2664       end subroutine global_all_sum_onevar_0i
2665     
2666       subroutine global_all_sum_onevar_1i( gbuf )
2667         integer, dimension(:), intent(inout) :: gbuf
2668         integer, dimension(size(gbuf)) :: lbuf
2669     
2670     #ifdef MPI
2671         lbuf = gbuf
2672         call global_all_sum_1i( lbuf, gbuf )
2673     #endif
2674         return
2675       end subroutine global_all_sum_onevar_1i
2676     
2677       subroutine global_all_sum_onevar_2i( gbuf )
2678         integer, dimension(:,:), intent(inout) :: gbuf
2679         integer, dimension(size(gbuf,1),size(gbuf,2)) :: lbuf
2680     
2681     #ifdef MPI
2682         lbuf = gbuf
2683         call global_all_sum_2i( lbuf, gbuf )
2684     #endif
2685         return
2686       end subroutine global_all_sum_onevar_2i
2687     
2688     
2689       subroutine global_all_sum_onevar_3i( gbuf )
2690         integer, dimension(:,:,:), intent(inout) :: gbuf
2691         integer, dimension(size(gbuf,1),size(gbuf,2),size(gbuf,3)) :: lbuf
2692     
2693     #ifdef MPI
2694         lbuf = gbuf
2695         call global_all_sum_3i( lbuf, gbuf )
2696     #endif
2697         return
2698       end subroutine global_all_sum_onevar_3i
2699     
2700       subroutine global_all_sum_onevar_0r( gbuf )
2701         real, intent(inout) :: gbuf
2702         real :: lbuf
2703     
2704     #ifdef MPI
2705         lbuf = gbuf
2706         call global_all_sum_0r( lbuf, gbuf )
2707     #endif
2708         return
2709       end subroutine global_all_sum_onevar_0r
2710     
2711     
2712       subroutine global_all_sum_onevar_1r( gbuf )
2713         real, dimension(:), intent(inout) :: gbuf
2714         real, dimension(size(gbuf)) :: lbuf
2715     
2716     #ifdef MPI
2717         lbuf = gbuf
2718         call global_all_sum_1r( lbuf, gbuf )
2719     #endif
2720         return
2721       end subroutine global_all_sum_onevar_1r
2722     
2723       subroutine global_all_sum_onevar_2r( gbuf )
2724         real, dimension(:,:), intent(inout) :: gbuf
2725         real, dimension(size(gbuf,1),size(gbuf,2)) :: lbuf
2726     
2727     #ifdef MPI
2728         lbuf = gbuf
2729         call global_all_sum_2r( lbuf, gbuf )
2730     #endif
2731         return
2732       end subroutine global_all_sum_onevar_2r
2733     
2734     
2735       subroutine global_all_sum_onevar_3r( gbuf )
2736         real, dimension(:,:,:), intent(inout) :: gbuf
2737         real, dimension(size(gbuf,1),size(gbuf,2),size(gbuf,3)) :: lbuf
2738     
2739     #ifdef MPI
2740         lbuf = gbuf
2741         call global_all_sum_3r( lbuf, gbuf )
2742     #endif
2743         return
2744       end subroutine global_all_sum_onevar_3r
2745     
2746     
2747       subroutine global_all_sum_0i( lbuf, gbuf, mroot, idebug )
2748         integer, intent(in) :: lbuf
2749         integer, intent(out) :: gbuf
2750         integer, optional, intent(in) :: mroot, idebug
2751     
2752     #ifdef MPI
2753         integer :: recvtype, sendtype, sendcnt, ierr,lroot, lidebug
2754     
2755         if (.not. present(mroot)) then
2756            lroot = 0
2757         else
2758            lroot = mroot
2759         endif
2760     
2761         if (.not. present(idebug)) then
2762            lidebug = 0
2763         else
2764            lidebug = idebug
2765         endif
2766     
2767         recvtype = MPI_INTEGER
2768         sendtype = recvtype
2769     
2770         sendcnt = 1
2771     
2772         call MPI_Allreduce( lbuf, gbuf, sendcnt, sendtype, MPI_SUM, &
2773              MPI_COMM_WORLD, ierr )
2774         call MPI_Check( 'global_all_sum_0i:MPI_Allreduce', ierr )
2775     #else
2776         gbuf = lbuf
2777     #endif
2778     
2779         return
2780       end subroutine global_all_sum_0i
2781     
2782     
2783       subroutine global_all_sum_1i( lbuf, gbuf, mroot, idebug )
2784         integer, intent(in), dimension(:) :: lbuf
2785         integer, intent(out), dimension(:) :: gbuf
2786         integer, optional, intent(in) :: mroot, idebug
2787     
2788     #ifdef MPI
2789         integer :: recvtype, sendtype, sendcnt, ierr,lroot, lidebug
2790     
2791         if (.not. present(mroot)) then
2792            lroot = 0
2793         else
2794            lroot = mroot
2795         endif
2796     
2797         if (.not. present(idebug)) then
2798            lidebug = 0
2799         else
2800            lidebug = idebug
2801         endif
2802     
2803         recvtype = MPI_INTEGER
2804         sendtype = recvtype
2805     
2806         sendcnt = size(lbuf)
2807     
2808         call MPI_Allreduce( lbuf, gbuf, sendcnt, sendtype, MPI_SUM, &
2809              MPI_COMM_WORLD, ierr )
2810         call MPI_Check( 'global_all_sum_1i:MPI_Allreduce', ierr )
2811     #else
2812         gbuf = lbuf
2813     #endif
2814     
2815         return
2816       end subroutine global_all_sum_1i
2817     
2818       subroutine global_all_sum_2i( lbuf, gbuf, mroot, idebug )
2819         integer, intent(in), dimension(:,:) :: lbuf
2820         integer, intent(out), dimension(:,:) :: gbuf
2821         integer, optional, intent(in) :: mroot, idebug
2822     
2823     #ifdef MPI
2824         integer :: i,j,lroot, lidebug
2825     
2826         if (.not. present(mroot)) then
2827            lroot = 0
2828         else
2829            lroot = mroot
2830         endif
2831     
2832         if (.not. present(idebug)) then
2833            lidebug = 0
2834         else
2835            lidebug = idebug
2836         endif
2837     
2838         call assert( size(lbuf,2).eq.size(gbuf,2),  &
2839              '** global_all_sum_2i: size(lbuf,2).ne.size(gbuf,2) ', &
2840              size(lbuf,2), size(gbuf,2) )
2841     
2842         do j=lbound(lbuf,2),ubound(lbuf,2)
2843            call global_all_sum_1i( lbuf(:,j), gbuf(:,j), lroot, lidebug )
2844         enddo
2845     #else
2846         gbuf = lbuf
2847     #endif
2848     
2849         return
2850       end subroutine global_all_sum_2i
2851     
2852       subroutine global_all_sum_3i( lbuf, gbuf, mroot, idebug )
2853         integer, intent(in), dimension(:,:,:) :: lbuf
2854         integer, intent(out), dimension(:,:,:) :: gbuf
2855         integer, optional, intent(in) :: mroot, idebug
2856     
2857     #ifdef MPI
2858         integer :: j,k,lroot, lidebug
2859     
2860         if (.not. present(mroot)) then
2861            lroot = 0
2862         else
2863            lroot = mroot
2864         endif
2865     
2866         if (.not. present(idebug)) then
2867            lidebug = 0
2868         else
2869            lidebug = idebug
2870         endif
2871     
2872         call assert( size(lbuf,2).eq.size(gbuf,2),  &
2873              '** global_all_sum_3i: size(lbuf,2).ne.size(gbuf,2) ', &
2874              size(lbuf,2), size(gbuf,2) )
2875     
2876         call assert( size(lbuf,3).eq.size(gbuf,3),  &
2877              '** global_all_sum_3i: size(lbuf,3).ne.size(gbuf,3) ', &
2878              size(lbuf,3), size(gbuf,3) )
2879     
2880         do k=lbound(lbuf,3),ubound(lbuf,3)
2881            do j=lbound(lbuf,2),ubound(lbuf,2)
2882               call global_all_sum_1i( lbuf(:,j,k), gbuf(:,j,k), lroot, lidebug )
2883            enddo
2884         enddo
2885     #else
2886         gbuf = lbuf
2887     #endif
2888     
2889         return
2890       end subroutine global_all_sum_3i
2891     
2892       subroutine global_all_sum_0r( lbuf, gbuf, mroot, idebug )
2893         real, intent(in) :: lbuf
2894         real, intent(out) :: gbuf
2895         integer, optional, intent(in) :: mroot, idebug
2896     
2897     #ifdef MPI
2898         integer :: recvtype, sendtype, sendcnt, ierr,lroot, lidebug
2899     
2900         if (.not. present(mroot)) then
2901            lroot = 0
2902         else
2903            lroot = mroot
2904         endif
2905     
2906         if (.not. present(idebug)) then
2907            lidebug = 0
2908         else
2909            lidebug = idebug
2910         endif
2911     
2912         recvtype = MPI_REAL
2913         sendtype = recvtype
2914     
2915         sendcnt = 1
2916     
2917         call MPI_Allreduce( lbuf, gbuf, sendcnt, sendtype, MPI_SUM, &
2918              MPI_COMM_WORLD, ierr )
2919         call MPI_Check( 'global_all_sum_0r:MPI_Allreduce', ierr )
2920     #else
2921         gbuf = lbuf
2922     #endif
2923     
2924         return
2925       end subroutine global_all_sum_0r
2926     
2927     
2928       subroutine global_all_sum_1r( lbuf, gbuf, mroot, idebug )
2929         real, intent(in), dimension(:) :: lbuf
2930         real, intent(out), dimension(:) :: gbuf
2931         integer, optional, intent(in) :: mroot, idebug
2932     
2933     #ifdef MPI
2934         integer :: recvtype, sendtype, sendcnt, ierr,lroot, lidebug
2935     
2936         if (.not. present(mroot)) then
2937            lroot = 0
2938         else
2939            lroot = mroot
2940         endif
2941     
2942         if (.not. present(idebug)) then
2943            lidebug = 0
2944         else
2945            lidebug = idebug
2946         endif
2947     
2948         recvtype = MPI_REAL
2949         sendtype = recvtype
2950     
2951         sendcnt = size(lbuf)
2952     
2953         call MPI_Allreduce( lbuf, gbuf, sendcnt, sendtype, MPI_SUM, &
2954              MPI_COMM_WORLD, ierr )
2955         call MPI_Check( 'global_all_sum_1r:MPI_Allreduce', ierr )
2956     #else
2957         gbuf = lbuf
2958     #endif
2959     
2960         return
2961       end subroutine global_all_sum_1r
2962     
2963       subroutine global_all_sum_2r( lbuf, gbuf, mroot, idebug )
2964         real, intent(in), dimension(:,:) :: lbuf
2965         real, intent(out), dimension(:,:) :: gbuf
2966         integer, optional, intent(in) :: mroot, idebug
2967     
2968     #ifdef MPI
2969         integer :: i,j,lroot, lidebug
2970     
2971         if (.not. present(mroot)) then
2972            lroot = 0
2973         else
2974            lroot = mroot
2975         endif
2976     
2977         if (.not. present(idebug)) then
2978            lidebug = 0
2979         else
2980            lidebug = idebug
2981         endif
2982     
2983         call assert( size(lbuf,2).eq.size(gbuf,2),  &
2984              '** global_all_sum_2r: size(lbuf,2).ne.size(gbuf,2) ', &
2985              size(lbuf,2), size(gbuf,2) )
2986     
2987         do j=lbound(lbuf,2),ubound(lbuf,2)
2988            call global_all_sum_1r( lbuf(:,j), gbuf(:,j), lroot, lidebug )
2989         enddo
2990     #else
2991         gbuf = lbuf
2992     #endif
2993     
2994         return
2995       end subroutine global_all_sum_2r
2996     
2997       subroutine global_all_sum_3r( lbuf, gbuf, mroot, idebug )
2998         real, intent(in), dimension(:,:,:) :: lbuf
2999         real, intent(out), dimension(:,:,:) :: gbuf
3000         integer, optional, intent(in) :: mroot, idebug
3001     
3002     #ifdef MPI
3003         integer :: j,k,lroot, lidebug
3004     
3005         if (.not. present(mroot)) then
3006            lroot = 0
3007         else
3008            lroot = mroot
3009         endif
3010     
3011         if (.not. present(idebug)) then
3012            lidebug = 0
3013         else
3014            lidebug = idebug
3015         endif
3016     
3017         call assert( size(lbuf,2).eq.size(gbuf,2),  &
3018              '** global_all_sum_3i: size(lbuf,2).ne.size(gbuf,2) ', &
3019              size(lbuf,2), size(gbuf,2) )
3020     
3021         call assert( size(lbuf,3).eq.size(gbuf,3),  &
3022              '** global_all_sum_3i: size(lbuf,3).ne.size(gbuf,3) ', &
3023              size(lbuf,3), size(gbuf,3) )
3024     
3025         do k=lbound(lbuf,3),ubound(lbuf,3)
3026            do j=lbound(lbuf,2),ubound(lbuf,2)
3027               call global_all_sum_1r( lbuf(:,j,k), gbuf(:,j,k), lroot, lidebug )
3028            enddo
3029         enddo
3030     #else
3031         gbuf = lbuf
3032     #endif
3033     
3034         return
3035       end subroutine global_all_sum_3r
3036     
3037       subroutine global_all_sum_0d( lbuf, gbuf, mroot, idebug )
3038         double precision, intent(in) :: lbuf
3039         double precision, intent(out) :: gbuf
3040         integer, optional, intent(in) :: mroot, idebug
3041     
3042     #ifdef MPI
3043         integer :: recvtype, sendtype, sendcnt, ierr,lroot, lidebug
3044     
3045         if (.not. present(mroot)) then
3046            lroot = 0
3047         else
3048            lroot = mroot
3049         endif
3050     
3051         if (.not. present(idebug)) then
3052            lidebug = 0
3053         else
3054            lidebug = idebug
3055         endif
3056     
3057         recvtype = MPI_DOUBLE_PRECISION
3058         sendtype = recvtype
3059     
3060         sendcnt = 1
3061     
3062         call MPI_Allreduce( lbuf, gbuf, sendcnt, sendtype, MPI_SUM, &
3063              MPI_COMM_WORLD, ierr )
3064         call MPI_Check( 'global_all_sum_0d:MPI_Allreduce', ierr )
3065     #else
3066         gbuf = lbuf
3067     #endif
3068     
3069         return
3070       end subroutine global_all_sum_0d
3071     
3072     
3073       subroutine global_all_sum_1d( lbuf, gbuf, mroot, idebug )
3074         double precision, intent(in), dimension(:) :: lbuf
3075         double precision, intent(out), dimension(:) :: gbuf
3076         integer, optional, intent(in) :: mroot, idebug
3077     
3078     #ifdef MPI
3079         integer :: recvtype, sendtype, sendcnt, ierr,lroot, lidebug
3080     
3081         if (.not. present(mroot)) then
3082            lroot = 0
3083         else
3084            lroot = mroot
3085         endif
3086     
3087         if (.not. present(idebug)) then
3088            lidebug = 0
3089         else
3090            lidebug = idebug
3091         endif
3092     
3093         recvtype = MPI_DOUBLE_PRECISION
3094         sendtype = recvtype
3095     
3096         sendcnt = size(lbuf)
3097     
3098         call MPI_Allreduce( lbuf, gbuf, sendcnt, sendtype,  MPI_SUM, &
3099              MPI_COMM_WORLD, ierr )
3100         call MPI_Check( 'global_all_sum_1d:MPI_Allreduce', ierr )
3101     #else
3102         gbuf = lbuf
3103     #endif
3104     
3105         return
3106       end subroutine global_all_sum_1d
3107     
3108       subroutine global_all_sum_2d( lbuf, gbuf, mroot, idebug )
3109         double precision, intent(in), dimension(:,:) :: lbuf
3110         double precision, intent(out), dimension(:,:) :: gbuf
3111         integer, optional, intent(in) :: mroot, idebug
3112     
3113     #ifdef MPI
3114         integer :: i,j,lroot, lidebug
3115     
3116         if (.not. present(mroot)) then
3117            lroot = 0
3118         else
3119            lroot = mroot
3120         endif
3121     
3122         if (.not. present(idebug)) then
3123            lidebug = 0
3124         else
3125            lidebug = idebug
3126         endif
3127     
3128         call assert( size(lbuf,2).eq.size(gbuf,2),  &
3129              '** global_all_sum_2d: size(lbuf,2).ne.size(gbuf,2) ', &
3130              size(lbuf,2), size(gbuf,2) )
3131     
3132         do j=lbound(lbuf,2),ubound(lbuf,2)
3133            call global_all_sum_1d( lbuf(:,j), gbuf(:,j), lroot, lidebug )
3134         enddo
3135     #else
3136         gbuf = lbuf
3137     #endif
3138     
3139         return
3140       end subroutine global_all_sum_2d
3141     
3142       subroutine global_all_sum_3d( lbuf, gbuf, mroot, idebug )
3143         double precision, intent(in), dimension(:,:,:) :: lbuf
3144         double precision, intent(out), dimension(:,:,:) :: gbuf
3145         integer, optional, intent(in) :: mroot, idebug
3146     
3147     #ifdef MPI
3148         integer :: j,k,lroot, lidebug
3149     
3150         if (.not. present(mroot)) then
3151            lroot = 0
3152         else
3153            lroot = mroot
3154         endif
3155     
3156         if (.not. present(idebug)) then
3157            lidebug = 0
3158         else
3159            lidebug = idebug
3160         endif
3161     
3162         call assert( size(lbuf,2).eq.size(gbuf,2),  &
3163              '** global_all_sum_3i: size(lbuf,2).ne.size(gbuf,2) ', &
3164              size(lbuf,2), size(gbuf,2) )
3165     
3166         call assert( size(lbuf,3).eq.size(gbuf,3),  &
3167              '** global_all_sum_3i: size(lbuf,3).ne.size(gbuf,3) ', &
3168              size(lbuf,3), size(gbuf,3) )
3169     
3170         do k=lbound(lbuf,3),ubound(lbuf,3)
3171            do j=lbound(lbuf,2),ubound(lbuf,2)
3172               call global_all_sum_1d( lbuf(:,j,k), gbuf(:,j,k), lroot, lidebug )
3173            enddo
3174         enddo
3175     #else
3176         gbuf = lbuf
3177     #endif
3178     
3179         return
3180       end subroutine global_all_sum_3d
3181     
3182       subroutine global_min_0i( lbuf, gbuf, mroot, idebug )
3183         integer, intent(in) :: lbuf
3184         integer, intent(out) :: gbuf
3185         integer, optional, intent(in) :: mroot, idebug
3186     
3187     #ifdef MPI
3188         integer :: recvtype, sendtype, sendcnt, ierr,lroot, lidebug
3189     
3190         if (.not. present(mroot)) then
3191            lroot = 0
3192         else
3193            lroot = mroot
3194         endif
3195     
3196         if (.not. present(idebug)) then
3197            lidebug = 0
3198         else
3199            lidebug = idebug
3200         endif
3201     
3202         recvtype = MPI_INTEGER
3203         sendtype = recvtype
3204     
3205         sendcnt = 1
3206     
3207         call MPI_Reduce( lbuf, gbuf, sendcnt, sendtype,  MPI_MIN, &
3208              lroot, MPI_COMM_WORLD, ierr )
3209         call MPI_Check( 'global_min_0i:MPI_Reduce', ierr )
3210     #else
3211         gbuf = lbuf
3212     #endif
3213     
3214         return
3215       end subroutine global_min_0i
3216     
3217     
3218       subroutine global_min_1i( lbuf, gbuf, mroot, idebug )
3219         integer, intent(in), dimension(:) :: lbuf
3220         integer, intent(out), dimension(:) :: gbuf
3221         integer, optional, intent(in) :: mroot, idebug
3222     
3223     #ifdef MPI
3224         integer :: recvtype, sendtype, sendcnt, ierr,lroot, lidebug
3225     
3226         if (.not. present(mroot)) then
3227            lroot = 0
3228         else
3229            lroot = mroot
3230         endif
3231     
3232         if (.not. present(idebug)) then
3233            lidebug = 0
3234         else
3235            lidebug = idebug
3236         endif
3237     
3238         recvtype = MPI_INTEGER
3239         sendtype = recvtype
3240     
3241         sendcnt = size(lbuf)
3242     
3243         call MPI_Reduce( lbuf, gbuf, sendcnt, sendtype,  MPI_MIN, &
3244              lroot, MPI_COMM_WORLD, ierr )
3245         call MPI_Check( 'global_min_1i:MPI_Reduce', ierr )
3246     #else
3247         gbuf = lbuf
3248     #endif
3249     
3250         return
3251       end subroutine global_min_1i
3252     
3253       subroutine global_min_2i( lbuf, gbuf, mroot, idebug )
3254         integer, intent(in), dimension(:,:) :: lbuf
3255         integer, intent(out), dimension(:,:) :: gbuf
3256         integer, optional, intent(in) :: mroot, idebug
3257     
3258     #ifdef MPI
3259         integer :: i,j,lroot, lidebug
3260     
3261         if (.not. present(mroot)) then
3262            lroot = 0
3263         else
3264            lroot = mroot
3265         endif
3266     
3267         if (.not. present(idebug)) then
3268            lidebug = 0
3269         else
3270            lidebug = idebug
3271         endif
3272     
3273         if(myPE.eq.lroot) then
3274            call assert( size(lbuf,2).eq.size(gbuf,2),  &
3275                 '** global_min_2i: size(lbuf,2).ne.size(gbuf,2) ', &
3276                 size(lbuf,2), size(gbuf,2) )
3277         endif
3278     
3279         do j=lbound(lbuf,2),ubound(lbuf,2)
3280            call global_min_1i( lbuf(:,j), gbuf(:,j), lroot, lidebug )
3281         enddo
3282     #else
3283         gbuf = lbuf
3284     #endif
3285     
3286         return
3287       end subroutine global_min_2i
3288     
3289       subroutine global_min_3i( lbuf, gbuf, mroot, idebug )
3290         integer, intent(in), dimension(:,:,:) :: lbuf
3291         integer, intent(out), dimension(:,:,:) :: gbuf
3292         integer, optional, intent(in) :: mroot, idebug
3293     
3294     #ifdef MPI
3295         integer :: j,k,lroot, lidebug
3296     
3297         if (.not. present(mroot)) then
3298            lroot = 0
3299         else
3300            lroot = mroot
3301         endif
3302     
3303         if (.not. present(idebug)) then
3304            lidebug = 0
3305         else
3306            lidebug = idebug
3307         endif
3308     
3309         if(myPE.eq.lroot) then
3310            call assert( size(lbuf,2).eq.size(gbuf,2),  &
3311                 '** global_min_3i: size(lbuf,2).ne.size(gbuf,2) ', &
3312                 size(lbuf,2), size(gbuf,2) )
3313     
3314            call assert( size(lbuf,3).eq.size(gbuf,3),  &
3315                 '** global_min_3i: size(lbuf,3).ne.size(gbuf,3) ', &
3316                 size(lbuf,3), size(gbuf,3) )
3317         endif
3318     
3319         do k=lbound(lbuf,3),ubound(lbuf,3)
3320            do j=lbound(lbuf,2),ubound(lbuf,2)
3321               call global_min_1i( lbuf(:,j,k), gbuf(:,j,k), lroot, lidebug )
3322            enddo
3323         enddo
3324     #else
3325         gbuf = lbuf
3326     #endif
3327     
3328         return
3329       end subroutine global_min_3i
3330     
3331       subroutine global_min_0r( lbuf, gbuf, mroot, idebug )
3332         real, intent(in) :: lbuf
3333         real, intent(out) :: gbuf
3334         integer, optional, intent(in) :: mroot, idebug
3335     
3336     #ifdef MPI
3337         integer :: recvtype, sendtype, sendcnt, ierr,lroot, lidebug
3338     
3339         if (.not. present(mroot)) then
3340            lroot = 0
3341         else
3342            lroot = mroot
3343         endif
3344     
3345         if (.not. present(idebug)) then
3346            lidebug = 0
3347         else
3348            lidebug = idebug
3349         endif
3350     
3351         recvtype = MPI_REAL
3352         sendtype = recvtype
3353     
3354         sendcnt = 1
3355     
3356         call MPI_Reduce( lbuf, gbuf, sendcnt, sendtype,  MPI_MIN, &
3357              lroot, MPI_COMM_WORLD, ierr )
3358         call MPI_Check( 'global_min_0r:MPI_Reduce', ierr )
3359     #else
3360         gbuf = lbuf
3361     #endif
3362     
3363         return
3364       end subroutine global_min_0r
3365     
3366     
3367       subroutine global_min_1r( lbuf, gbuf, mroot, idebug )
3368         real, intent(in), dimension(:) :: lbuf
3369         real, intent(out), dimension(:) :: gbuf
3370         integer, optional, intent(in) :: mroot, idebug
3371     
3372     #ifdef MPI
3373         integer :: recvtype, sendtype, sendcnt, ierr,lroot, lidebug
3374     
3375         if (.not. present(mroot)) then
3376            lroot = 0
3377         else
3378            lroot = mroot
3379         endif
3380     
3381         if (.not. present(idebug)) then
3382            lidebug = 0
3383         else
3384            lidebug = idebug
3385         endif
3386     
3387         recvtype = MPI_REAL
3388         sendtype = recvtype
3389     
3390         sendcnt = size(lbuf)
3391     
3392         call MPI_Reduce( lbuf, gbuf, sendcnt, sendtype,  MPI_MIN, &
3393              lroot, MPI_COMM_WORLD, ierr )
3394         call MPI_Check( 'global_min_1r:MPI_Reduce', ierr )
3395     #else
3396         gbuf = lbuf
3397     #endif
3398     
3399         return
3400       end subroutine global_min_1r
3401     
3402       subroutine global_min_2r( lbuf, gbuf, mroot, idebug )
3403         real, intent(in), dimension(:,:) :: lbuf
3404         real, intent(out), dimension(:,:) :: gbuf
3405         integer, optional, intent(in) :: mroot, idebug
3406     
3407     #ifdef MPI
3408         integer :: i,j,lroot, lidebug
3409     
3410         if (.not. present(mroot)) then
3411            lroot = 0
3412         else
3413            lroot = mroot
3414         endif
3415     
3416         if (.not. present(idebug)) then
3417            lidebug = 0
3418         else
3419            lidebug = idebug
3420         endif
3421     
3422         if(myPE.eq.lroot) then
3423            call assert( size(lbuf,2).eq.size(gbuf,2),  &
3424                 '** global_min_2r: size(lbuf,2).ne.size(gbuf,2) ', &
3425                 size(lbuf,2), size(gbuf,2) )
3426         endif
3427     
3428         do j=lbound(lbuf,2),ubound(lbuf,2)
3429            call global_min_1r( lbuf(:,j), gbuf(:,j), lroot, lidebug )
3430         enddo
3431     #else
3432         gbuf = lbuf
3433     #endif
3434     
3435         return
3436       end subroutine global_min_2r
3437     
3438       subroutine global_min_3r( lbuf, gbuf, mroot, idebug )
3439         real, intent(in), dimension(:,:,:) :: lbuf
3440         real, intent(out), dimension(:,:,:) :: gbuf
3441         integer, optional, intent(in) :: mroot, idebug
3442     
3443     #ifdef MPI
3444         integer :: j,k,lroot, lidebug
3445     
3446         if (.not. present(mroot)) then
3447            lroot = 0
3448         else
3449            lroot = mroot
3450         endif
3451     
3452         if (.not. present(idebug)) then
3453            lidebug = 0
3454         else
3455            lidebug = idebug
3456         endif
3457     
3458         if(myPE.eq.lroot) then
3459            call assert( size(lbuf,2).eq.size(gbuf,2),  &
3460                 '** global_min_3i: size(lbuf,2).ne.size(gbuf,2) ', &
3461                 size(lbuf,2), size(gbuf,2) )
3462     
3463            call assert( size(lbuf,3).eq.size(gbuf,3),  &
3464                 '** global_min_3i: size(lbuf,3).ne.size(gbuf,3) ', &
3465                 size(lbuf,3), size(gbuf,3) )
3466         endif
3467     
3468         do k=lbound(lbuf,3),ubound(lbuf,3)
3469            do j=lbound(lbuf,2),ubound(lbuf,2)
3470               call global_min_1r( lbuf(:,j,k), gbuf(:,j,k), lroot, lidebug )
3471            enddo
3472         enddo
3473     #else
3474         gbuf = lbuf
3475     #endif
3476     
3477         return
3478       end subroutine global_min_3r
3479     
3480       subroutine global_min_0d( lbuf, gbuf, mroot, idebug )
3481         double precision, intent(in) :: lbuf
3482         double precision, intent(out) :: gbuf
3483         integer, optional, intent(in) :: mroot, idebug
3484     
3485     #ifdef MPI
3486         integer :: recvtype, sendtype, sendcnt, ierr,lroot, lidebug
3487     
3488         if (.not. present(mroot)) then
3489            lroot = 0
3490         else
3491            lroot = mroot
3492         endif
3493     
3494         if (.not. present(idebug)) then
3495            lidebug = 0
3496         else
3497            lidebug = idebug
3498         endif
3499     
3500         recvtype = MPI_DOUBLE_PRECISION
3501         sendtype = recvtype
3502     
3503         sendcnt = 1
3504     
3505         call MPI_Reduce( lbuf, gbuf, sendcnt, sendtype,  MPI_MIN, &
3506              lroot, MPI_COMM_WORLD, ierr )
3507         call MPI_Check( 'global_min_0d:MPI_Reduce', ierr )
3508     #else
3509         gbuf = lbuf
3510     #endif
3511     
3512         return
3513       end subroutine global_min_0d
3514     
3515     
3516       subroutine global_min_1d( lbuf, gbuf, mroot, idebug )
3517         double precision, intent(in), dimension(:) :: lbuf
3518         double precision, intent(out), dimension(:) :: gbuf
3519         integer, optional, intent(in) :: mroot, idebug
3520     
3521     #ifdef MPI
3522         integer :: recvtype, sendtype, sendcnt, ierr,lroot, lidebug
3523     
3524         if (.not. present(mroot)) then
3525            lroot = 0
3526         else
3527            lroot = mroot
3528         endif
3529     
3530         if (.not. present(idebug)) then
3531            lidebug = 0
3532         else
3533            lidebug = idebug
3534         endif
3535     
3536         recvtype = MPI_DOUBLE_PRECISION
3537         sendtype = recvtype
3538     
3539         sendcnt = size(lbuf)
3540     
3541         call MPI_Reduce( lbuf, gbuf, sendcnt, sendtype,  MPI_MIN, &
3542              lroot, MPI_COMM_WORLD, ierr )
3543         call MPI_Check( 'global_min_1d:MPI_Reduce', ierr )
3544     #else
3545         gbuf = lbuf
3546     #endif
3547     
3548         return
3549       end subroutine global_min_1d
3550     
3551       subroutine global_min_2d( lbuf, gbuf, mroot, idebug )
3552         double precision, intent(in), dimension(:,:) :: lbuf
3553         double precision, intent(out), dimension(:,:) :: gbuf
3554         integer, optional, intent(in) :: mroot, idebug
3555     
3556     #ifdef MPI
3557         integer :: i,j,lroot, lidebug
3558     
3559         if (.not. present(mroot)) then
3560            lroot = 0
3561         else
3562            lroot = mroot
3563         endif
3564     
3565         if (.not. present(idebug)) then
3566            lidebug = 0
3567         else
3568            lidebug = idebug
3569         endif
3570     
3571         if(myPE.eq.lroot) then
3572            call assert( size(lbuf,2).eq.size(gbuf,2),  &
3573                 '** global_min_2d: size(lbuf,2).ne.size(gbuf,2) ', &
3574                 size(lbuf,2), size(gbuf,2) )
3575         endif
3576     
3577         do j=lbound(lbuf,2),ubound(lbuf,2)
3578            call global_min_1d( lbuf(:,j), gbuf(:,j), lroot, lidebug )
3579         enddo
3580     #else
3581         gbuf = lbuf
3582     #endif
3583     
3584         return
3585       end subroutine global_min_2d
3586     
3587       subroutine global_min_3d( lbuf, gbuf, mroot, idebug )
3588         double precision, intent(in), dimension(:,:,:) :: lbuf
3589         double precision, intent(out), dimension(:,:,:) :: gbuf
3590         integer, optional, intent(in) :: mroot, idebug
3591     
3592     #ifdef MPI
3593         integer :: j,k,lroot, lidebug
3594     
3595         if (.not. present(mroot)) then
3596            lroot = 0
3597         else
3598            lroot = mroot
3599         endif
3600     
3601         if (.not. present(idebug)) then
3602            lidebug = 0
3603         else
3604            lidebug = idebug
3605         endif
3606     
3607         if(myPE.eq.lroot) then
3608            call assert( size(lbuf,2).eq.size(gbuf,2),  &
3609                 '** global_min_3i: size(lbuf,2).ne.size(gbuf,2) ', &
3610                 size(lbuf,2), size(gbuf,2) )
3611     
3612            call assert( size(lbuf,3).eq.size(gbuf,3),  &
3613                 '** global_min_3i: size(lbuf,3).ne.size(gbuf,3) ', &
3614                 size(lbuf,3), size(gbuf,3) )
3615         endif
3616     
3617         do k=lbound(lbuf,3),ubound(lbuf,3)
3618            do j=lbound(lbuf,2),ubound(lbuf,2)
3619               call global_min_1d( lbuf(:,j,k), gbuf(:,j,k), lroot, lidebug )
3620            enddo
3621         enddo
3622     #else
3623         gbuf = lbuf
3624     #endif
3625     
3626         return
3627       end subroutine global_min_3d
3628     
3629       subroutine global_all_min_onevar_0d( gbuf )
3630         doubleprecision, intent(inout) :: gbuf
3631         doubleprecision :: lbuf
3632     
3633     #ifdef MPI
3634         lbuf = gbuf
3635         call global_all_min_0d( lbuf, gbuf )
3636     #endif
3637         return
3638       end subroutine global_all_min_onevar_0d
3639     
3640     
3641       subroutine global_all_min_onevar_1d( gbuf )
3642         doubleprecision, dimension(:), intent(inout) :: gbuf
3643         doubleprecision, dimension(size(gbuf)) :: lbuf
3644     
3645     #ifdef MPI
3646         lbuf = gbuf
3647         call global_all_min_1d( lbuf, gbuf )
3648     #endif
3649         return
3650       end subroutine global_all_min_onevar_1d
3651     
3652       subroutine global_all_min_onevar_2d( gbuf )
3653         doubleprecision, dimension(:,:), intent(inout) :: gbuf
3654         doubleprecision, dimension(size(gbuf,1),size(gbuf,2)) :: lbuf
3655     
3656     #ifdef MPI
3657         lbuf = gbuf
3658         call global_all_min_2d( lbuf, gbuf )
3659     #endif
3660         return
3661       end subroutine global_all_min_onevar_2d
3662     
3663     
3664       subroutine global_all_min_onevar_3d( gbuf )
3665         doubleprecision, dimension(:,:,:), intent(inout) :: gbuf
3666         doubleprecision, dimension(size(gbuf,1),size(gbuf,2),size(gbuf,3)) :: lbuf
3667     
3668     #ifdef MPI
3669         lbuf = gbuf
3670         call global_all_min_3d( lbuf, gbuf )
3671     #endif
3672         return
3673       end subroutine global_all_min_onevar_3d
3674     
3675     
3676     
3677     
3678       subroutine global_all_min_onevar_0i( gbuf )
3679         integer, intent(inout) :: gbuf
3680         integer :: lbuf
3681     
3682     #ifdef MPI
3683         lbuf = gbuf
3684         call global_all_min_0i( lbuf, gbuf )
3685     #endif
3686         return
3687       end subroutine global_all_min_onevar_0i
3688     
3689     
3690       subroutine global_all_min_onevar_1i( gbuf )
3691         integer, dimension(:), intent(inout) :: gbuf
3692         integer, dimension(size(gbuf)) :: lbuf
3693     
3694     #ifdef MPI
3695         lbuf = gbuf
3696         call global_all_min_1i( lbuf, gbuf )
3697     #endif
3698         return
3699       end subroutine global_all_min_onevar_1i
3700     
3701       subroutine global_all_min_onevar_2i( gbuf )
3702         integer, dimension(:,:), intent(inout) :: gbuf
3703         integer, dimension(size(gbuf,1),size(gbuf,2)) :: lbuf
3704     
3705     #ifdef MPI
3706         lbuf = gbuf
3707         call global_all_min_2i( lbuf, gbuf )
3708     #endif
3709         return
3710       end subroutine global_all_min_onevar_2i
3711     
3712     
3713       subroutine global_all_min_onevar_3i( gbuf )
3714         integer, dimension(:,:,:), intent(inout) :: gbuf
3715         integer, dimension(size(gbuf,1),size(gbuf,2),size(gbuf,3)) :: lbuf
3716     
3717     #ifdef MPI
3718         lbuf = gbuf
3719         call global_all_min_3i( lbuf, gbuf )
3720     #endif
3721         return
3722       end subroutine global_all_min_onevar_3i
3723     
3724       subroutine global_all_min_onevar_0r( gbuf )
3725         real, intent(inout) :: gbuf
3726         real :: lbuf
3727     
3728     #ifdef MPI
3729         lbuf = gbuf
3730         call global_all_min_0r( lbuf, gbuf )
3731     #endif
3732         return
3733       end subroutine global_all_min_onevar_0r
3734     
3735     
3736       subroutine global_all_min_onevar_1r( gbuf )
3737         real, dimension(:), intent(inout) :: gbuf
3738         real, dimension(size(gbuf)) :: lbuf
3739     
3740     #ifdef MPI
3741         lbuf = gbuf
3742         call global_all_min_1r( lbuf, gbuf )
3743     #endif
3744         return
3745       end subroutine global_all_min_onevar_1r
3746     
3747       subroutine global_all_min_onevar_2r( gbuf )
3748         real, dimension(:,:), intent(inout) :: gbuf
3749         real, dimension(size(gbuf,1),size(gbuf,2)) :: lbuf
3750     
3751     #ifdef MPI
3752         lbuf = gbuf
3753         call global_all_min_2r( lbuf, gbuf )
3754     #endif
3755         return
3756       end subroutine global_all_min_onevar_2r
3757     
3758     
3759     
3760       subroutine global_all_min_onevar_3r( gbuf )
3761         real, dimension(:,:,:), intent(inout) :: gbuf
3762         real, dimension(size(gbuf,1),size(gbuf,2),size(gbuf,3)) :: lbuf
3763     
3764     #ifdef MPI
3765         lbuf = gbuf
3766         call global_all_min_3r( lbuf, gbuf )
3767     #endif
3768         return
3769       end subroutine global_all_min_onevar_3r
3770     
3771     
3772       subroutine global_all_min_0i( lbuf, gbuf, mroot, idebug )
3773         integer, intent(in) :: lbuf
3774         integer, intent(out) :: gbuf
3775         integer, optional, intent(in) :: mroot, idebug
3776     
3777     #ifdef MPI
3778         integer :: recvtype, sendtype, sendcnt, ierr,lroot, lidebug
3779     
3780         if (.not. present(mroot)) then
3781            lroot = 0
3782         else
3783            lroot = mroot
3784         endif
3785     
3786         if (.not. present(idebug)) then
3787            lidebug = 0
3788         else
3789            lidebug = idebug
3790         endif
3791     
3792         recvtype = MPI_INTEGER
3793         sendtype = recvtype
3794     
3795         sendcnt = 1
3796     
3797         call MPI_Allreduce( lbuf, gbuf, sendcnt, sendtype, MPI_MIN, &
3798              MPI_COMM_WORLD, ierr )
3799         call MPI_Check( 'global_all_min_0i:MPI_Allreduce', ierr )
3800     #else
3801         gbuf = lbuf
3802     #endif
3803     
3804         return
3805       end subroutine global_all_min_0i
3806     
3807     
3808       subroutine global_all_min_1i( lbuf, gbuf, mroot, idebug )
3809         integer, intent(in), dimension(:) :: lbuf
3810         integer, intent(out), dimension(:) :: gbuf
3811         integer, optional, intent(in) :: mroot, idebug
3812     
3813     #ifdef MPI
3814         integer :: recvtype, sendtype, sendcnt, ierr,lroot, lidebug
3815     
3816         if (.not. present(mroot)) then
3817            lroot = 0
3818         else
3819            lroot = mroot
3820         endif
3821     
3822         if (.not. present(idebug)) then
3823            lidebug = 0
3824         else
3825            lidebug = idebug
3826         endif
3827     
3828         recvtype = MPI_INTEGER
3829         sendtype = recvtype
3830     
3831         sendcnt = size(lbuf)
3832     
3833         call MPI_Allreduce( lbuf, gbuf, sendcnt, sendtype, MPI_MIN, &
3834              MPI_COMM_WORLD, ierr )
3835         call MPI_Check( 'global_all_min_1i:MPI_Allreduce', ierr )
3836     #else
3837         gbuf = lbuf
3838     #endif
3839     
3840         return
3841       end subroutine global_all_min_1i
3842     
3843       subroutine global_all_min_2i( lbuf, gbuf, mroot, idebug )
3844         integer, intent(in), dimension(:,:) :: lbuf
3845         integer, intent(out), dimension(:,:) :: gbuf
3846         integer, optional, intent(in) :: mroot, idebug
3847     
3848     #ifdef MPI
3849         integer :: i,j,lroot, lidebug
3850     
3851         if (.not. present(mroot)) then
3852            lroot = 0
3853         else
3854            lroot = mroot
3855         endif
3856     
3857         if (.not. present(idebug)) then
3858            lidebug = 0
3859         else
3860            lidebug = idebug
3861         endif
3862     
3863         call assert( size(lbuf,2).eq.size(gbuf,2),  &
3864              '** global_all_min_2i: size(lbuf,2).ne.size(gbuf,2) ', &
3865              size(lbuf,2), size(gbuf,2) )
3866     
3867         do j=lbound(lbuf,2),ubound(lbuf,2)
3868            call global_all_min_1i( lbuf(:,j), gbuf(:,j), lroot, lidebug )
3869         enddo
3870     #else
3871         gbuf = lbuf
3872     #endif
3873     
3874         return
3875       end subroutine global_all_min_2i
3876     
3877       subroutine global_all_min_3i( lbuf, gbuf, mroot, idebug )
3878         integer, intent(in), dimension(:,:,:) :: lbuf
3879         integer, intent(out), dimension(:,:,:) :: gbuf
3880         integer, optional, intent(in) :: mroot, idebug
3881     
3882     #ifdef MPI
3883         integer :: j,k,lroot, lidebug
3884     
3885         if (.not. present(mroot)) then
3886            lroot = 0
3887         else
3888            lroot = mroot
3889         endif
3890     
3891         if (.not. present(idebug)) then
3892            lidebug = 0
3893         else
3894            lidebug = idebug
3895         endif
3896     
3897         call assert( size(lbuf,2).eq.size(gbuf,2),  &
3898              '** global_all_min_3i: size(lbuf,2).ne.size(gbuf,2) ', &
3899              size(lbuf,2), size(gbuf,2) )
3900     
3901         call assert( size(lbuf,3).eq.size(gbuf,3),  &
3902              '** global_all_min_3i: size(lbuf,3).ne.size(gbuf,3) ', &
3903              size(lbuf,3), size(gbuf,3) )
3904     
3905         do k=lbound(lbuf,3),ubound(lbuf,3)
3906            do j=lbound(lbuf,2),ubound(lbuf,2)
3907               call global_all_min_1i( lbuf(:,j,k), gbuf(:,j,k), lroot, lidebug )
3908            enddo
3909         enddo
3910     #else
3911         gbuf = lbuf
3912     #endif
3913     
3914         return
3915       end subroutine global_all_min_3i
3916     
3917       subroutine global_all_min_0r( lbuf, gbuf, mroot, idebug )
3918         real, intent(in) :: lbuf
3919         real, intent(out) :: gbuf
3920         integer, optional, intent(in) :: mroot, idebug
3921     
3922     #ifdef MPI
3923         integer :: recvtype, sendtype, sendcnt, ierr,lroot, lidebug
3924     
3925         if (.not. present(mroot)) then
3926            lroot = 0
3927         else
3928            lroot = mroot
3929         endif
3930     
3931         if (.not. present(idebug)) then
3932            lidebug = 0
3933         else
3934            lidebug = idebug
3935         endif
3936     
3937         recvtype = MPI_REAL
3938         sendtype = recvtype
3939     
3940         sendcnt = 1
3941     
3942         call MPI_Allreduce( lbuf, gbuf, sendcnt, sendtype, MPI_MIN, &
3943              MPI_COMM_WORLD, ierr )
3944         call MPI_Check( 'global_all_min_0r:MPI_Allreduce', ierr )
3945     #else
3946         gbuf = lbuf
3947     #endif
3948     
3949         return
3950       end subroutine global_all_min_0r
3951     
3952     
3953       subroutine global_all_min_1r( lbuf, gbuf, mroot, idebug )
3954         real, intent(in), dimension(:) :: lbuf
3955         real, intent(out), dimension(:) :: gbuf
3956         integer, optional, intent(in) :: mroot, idebug
3957     
3958     #ifdef MPI
3959         integer :: recvtype, sendtype, sendcnt, ierr,lroot, lidebug
3960     
3961         if (.not. present(mroot)) then
3962            lroot = 0
3963         else
3964            lroot = mroot
3965         endif
3966     
3967         if (.not. present(idebug)) then
3968            lidebug = 0
3969         else
3970            lidebug = idebug
3971         endif
3972     
3973         recvtype = MPI_REAL
3974         sendtype = recvtype
3975     
3976         sendcnt = size(lbuf)
3977     
3978         call MPI_Allreduce( lbuf, gbuf, sendcnt, sendtype, MPI_MIN, &
3979              MPI_COMM_WORLD, ierr )
3980         call MPI_Check( 'global_all_min_1r:MPI_Allreduce', ierr )
3981     #else
3982         gbuf = lbuf
3983     #endif
3984     
3985         return
3986       end subroutine global_all_min_1r
3987     
3988       subroutine global_all_min_2r( lbuf, gbuf, mroot, idebug )
3989         real, intent(in), dimension(:,:) :: lbuf
3990         real, intent(out), dimension(:,:) :: gbuf
3991         integer, optional, intent(in) :: mroot, idebug
3992     
3993     #ifdef MPI
3994         integer :: i,j,lroot, lidebug
3995     
3996         if (.not. present(mroot)) then
3997            lroot = 0
3998         else
3999            lroot = mroot
4000         endif
4001     
4002         if (.not. present(idebug)) then
4003            lidebug = 0
4004         else
4005            lidebug = idebug
4006         endif
4007     
4008         call assert( size(lbuf,2).eq.size(gbuf,2),  &
4009              '** global_all_min_2r: size(lbuf,2).ne.size(gbuf,2) ', &
4010              size(lbuf,2), size(gbuf,2) )
4011     
4012         do j=lbound(lbuf,2),ubound(lbuf,2)
4013            call global_all_min_1r( lbuf(:,j), gbuf(:,j), lroot, lidebug )
4014         enddo
4015     #else
4016         gbuf = lbuf
4017     #endif
4018     
4019         return
4020       end subroutine global_all_min_2r
4021     
4022       subroutine global_all_min_3r( lbuf, gbuf, mroot, idebug )
4023         real, intent(in), dimension(:,:,:) :: lbuf
4024         real, intent(out), dimension(:,:,:) :: gbuf
4025         integer, optional, intent(in) :: mroot, idebug
4026     
4027     #ifdef MPI
4028         integer :: j,k,lroot, lidebug
4029     
4030         if (.not. present(mroot)) then
4031            lroot = 0
4032         else
4033            lroot = mroot
4034         endif
4035     
4036         if (.not. present(idebug)) then
4037            lidebug = 0
4038         else
4039            lidebug = idebug
4040         endif
4041     
4042         call assert( size(lbuf,2).eq.size(gbuf,2),  &
4043              '** global_all_min_3i: size(lbuf,2).ne.size(gbuf,2) ', &
4044              size(lbuf,2), size(gbuf,2) )
4045     
4046         call assert( size(lbuf,3).eq.size(gbuf,3),  &
4047              '** global_all_min_3i: size(lbuf,3).ne.size(gbuf,3) ', &
4048              size(lbuf,3), size(gbuf,3) )
4049     
4050         do k=lbound(lbuf,3),ubound(lbuf,3)
4051            do j=lbound(lbuf,2),ubound(lbuf,2)
4052               call global_all_min_1r( lbuf(:,j,k), gbuf(:,j,k), lroot, lidebug )
4053            enddo
4054         enddo
4055     #else
4056         gbuf = lbuf
4057     #endif
4058     
4059         return
4060       end subroutine global_all_min_3r
4061     
4062       subroutine global_all_min_0d( lbuf, gbuf, mroot, idebug )
4063         double precision, intent(in) :: lbuf
4064         double precision, intent(out) :: gbuf
4065         integer, optional, intent(in) :: mroot, idebug
4066     
4067     #ifdef MPI
4068         integer :: recvtype, sendtype, sendcnt, ierr,lroot, lidebug
4069     
4070         if (.not. present(mroot)) then
4071            lroot = 0
4072         else
4073            lroot = mroot
4074         endif
4075     
4076         if (.not. present(idebug)) then
4077            lidebug = 0
4078         else
4079            lidebug = idebug
4080         endif
4081     
4082         recvtype = MPI_DOUBLE_PRECISION
4083         sendtype = recvtype
4084     
4085         sendcnt = 1
4086     
4087         call MPI_Allreduce( lbuf, gbuf, sendcnt, sendtype, MPI_MIN, &
4088              MPI_COMM_WORLD, ierr )
4089         call MPI_Check( 'global_all_min_0d:MPI_Allreduce', ierr )
4090     #else
4091         gbuf = lbuf
4092     #endif
4093     
4094         return
4095       end subroutine global_all_min_0d
4096     
4097     
4098       subroutine global_all_min_1d( lbuf, gbuf, mroot, idebug )
4099         double precision, intent(in), dimension(:) :: lbuf
4100         double precision, intent(out), dimension(:) :: gbuf
4101         integer, optional, intent(in) :: mroot, idebug
4102     
4103     #ifdef MPI
4104         integer :: recvtype, sendtype, sendcnt, ierr,lroot, lidebug
4105     
4106         if (.not. present(mroot)) then
4107            lroot = 0
4108         else
4109            lroot = mroot
4110         endif
4111     
4112         if (.not. present(idebug)) then
4113            lidebug = 0
4114         else
4115            lidebug = idebug
4116         endif
4117     
4118         recvtype = MPI_DOUBLE_PRECISION
4119         sendtype = recvtype
4120     
4121         sendcnt = size(lbuf)
4122     
4123         call MPI_Allreduce( lbuf, gbuf, sendcnt, sendtype,  MPI_MIN, &
4124              MPI_COMM_WORLD, ierr )
4125         call MPI_Check( 'global_all_min_1d:MPI_Allreduce', ierr )
4126     #else
4127         gbuf = lbuf
4128     #endif
4129     
4130         return
4131       end subroutine global_all_min_1d
4132     
4133       subroutine global_all_min_2d( lbuf, gbuf, mroot, idebug )
4134         double precision, intent(in), dimension(:,:) :: lbuf
4135         double precision, intent(out), dimension(:,:) :: gbuf
4136         integer, optional, intent(in) :: mroot, idebug
4137     
4138     #ifdef MPI
4139         integer :: i,j,lroot, lidebug
4140     
4141         if (.not. present(mroot)) then
4142            lroot = 0
4143         else
4144            lroot = mroot
4145         endif
4146     
4147         if (.not. present(idebug)) then
4148            lidebug = 0
4149         else
4150            lidebug = idebug
4151         endif
4152     
4153         call assert( size(lbuf,2).eq.size(gbuf,2),  &
4154              '** global_all_min_2d: size(lbuf,2).ne.size(gbuf,2) ', &
4155              size(lbuf,2), size(gbuf,2) )
4156     
4157         do j=lbound(lbuf,2),ubound(lbuf,2)
4158            call global_all_min_1d( lbuf(:,j), gbuf(:,j), lroot, lidebug )
4159         enddo
4160     #else
4161         gbuf = lbuf
4162     #endif
4163     
4164         return
4165       end subroutine global_all_min_2d
4166     
4167       subroutine global_all_min_3d( lbuf, gbuf, mroot, idebug )
4168         double precision, intent(in), dimension(:,:,:) :: lbuf
4169         double precision, intent(out), dimension(:,:,:) :: gbuf
4170         integer, optional, intent(in) :: mroot, idebug
4171     
4172     #ifdef MPI
4173         integer :: j,k,lroot, lidebug
4174     
4175         if (.not. present(mroot)) then
4176            lroot = 0
4177         else
4178            lroot = mroot
4179         endif
4180     
4181         if (.not. present(idebug)) then
4182            lidebug = 0
4183         else
4184            lidebug = idebug
4185         endif
4186     
4187         call assert( size(lbuf,2).eq.size(gbuf,2),  &
4188              '** global_all_min_3i: size(lbuf,2).ne.size(gbuf,2) ', &
4189              size(lbuf,2), size(gbuf,2) )
4190     
4191         call assert( size(lbuf,3).eq.size(gbuf,3),  &
4192              '** global_all_min_3i: size(lbuf,3).ne.size(gbuf,3) ', &
4193              size(lbuf,3), size(gbuf,3) )
4194     
4195         do k=lbound(lbuf,3),ubound(lbuf,3)
4196            do j=lbound(lbuf,2),ubound(lbuf,2)
4197               call global_all_min_1d( lbuf(:,j,k), gbuf(:,j,k), lroot, lidebug )
4198            enddo
4199         enddo
4200     #else
4201         gbuf = lbuf
4202     #endif
4203     
4204         return
4205       end subroutine global_all_min_3d
4206     
4207       subroutine global_max_0i( lbuf, gbuf, mroot, idebug )
4208         integer, intent(in) :: lbuf
4209         integer, intent(out) :: gbuf
4210         integer, optional, intent(in) :: mroot, idebug
4211     
4212     #ifdef MPI
4213         integer :: recvtype, sendtype, sendcnt, ierr,lroot, lidebug
4214     
4215         if (.not. present(mroot)) then
4216            lroot = 0
4217         else
4218            lroot = mroot
4219         endif
4220     
4221         if (.not. present(idebug)) then
4222            lidebug = 0
4223         else
4224            lidebug = idebug
4225         endif
4226     
4227         recvtype = MPI_INTEGER
4228         sendtype = recvtype
4229     
4230         sendcnt = 1
4231     
4232         call MPI_Reduce( lbuf, gbuf, sendcnt, sendtype,  MPI_MAX, &
4233              lroot, MPI_COMM_WORLD, ierr )
4234         call MPI_Check( 'global_max_0i:MPI_Reduce', ierr )
4235     #else
4236         gbuf = lbuf
4237     #endif
4238     
4239         return
4240       end subroutine global_max_0i
4241     
4242     
4243       subroutine global_max_1i( lbuf, gbuf, mroot, idebug )
4244         integer, intent(in), dimension(:) :: lbuf
4245         integer, intent(out), dimension(:) :: gbuf
4246         integer, optional, intent(in) :: mroot, idebug
4247     
4248     #ifdef MPI
4249         integer :: recvtype, sendtype, sendcnt, ierr,lroot, lidebug
4250     
4251         if (.not. present(mroot)) then
4252            lroot = 0
4253         else
4254            lroot = mroot
4255         endif
4256     
4257         if (.not. present(idebug)) then
4258            lidebug = 0
4259         else
4260            lidebug = idebug
4261         endif
4262     
4263         recvtype = MPI_INTEGER
4264         sendtype = recvtype
4265     
4266         sendcnt = size(lbuf)
4267     
4268         call MPI_Reduce( lbuf, gbuf, sendcnt, sendtype,  MPI_MAX, &
4269              lroot, MPI_COMM_WORLD, ierr )
4270         call MPI_Check( 'global_max_1i:MPI_Reduce', ierr )
4271     #else
4272         gbuf = lbuf
4273     #endif
4274     
4275         return
4276       end subroutine global_max_1i
4277     
4278       subroutine global_max_2i( lbuf, gbuf, mroot, idebug )
4279         integer, intent(in), dimension(:,:) :: lbuf
4280         integer, intent(out), dimension(:,:) :: gbuf
4281         integer, optional, intent(in) :: mroot, idebug
4282     
4283     #ifdef MPI
4284         integer :: i,j,lroot, lidebug
4285     
4286         if (.not. present(mroot)) then
4287            lroot = 0
4288         else
4289            lroot = mroot
4290         endif
4291     
4292         if (.not. present(idebug)) then
4293            lidebug = 0
4294         else
4295            lidebug = idebug
4296         endif
4297     
4298         if(myPE.eq.lroot) then
4299            call assert( size(lbuf,2).eq.size(gbuf,2),  &
4300                 '** global_max_2i: size(lbuf,2).ne.size(gbuf,2) ', &
4301                 size(lbuf,2), size(gbuf,2) )
4302         endif
4303     
4304         do j=lbound(lbuf,2),ubound(lbuf,2)
4305            call global_max_1i( lbuf(:,j), gbuf(:,j), lroot, lidebug )
4306         enddo
4307     #else
4308         gbuf = lbuf
4309     #endif
4310     
4311         return
4312       end subroutine global_max_2i
4313     
4314       subroutine global_max_3i( lbuf, gbuf, mroot, idebug )
4315         integer, intent(in), dimension(:,:,:) :: lbuf
4316         integer, intent(out), dimension(:,:,:) :: gbuf
4317         integer, optional, intent(in) :: mroot, idebug
4318     
4319     #ifdef MPI
4320         integer :: j,k,lroot, lidebug
4321     
4322         if (.not. present(mroot)) then
4323            lroot = 0
4324         else
4325            lroot = mroot
4326         endif
4327     
4328         if (.not. present(idebug)) then
4329            lidebug = 0
4330         else
4331            lidebug = idebug
4332         endif
4333     
4334         if(myPE.eq.lroot) then
4335            call assert( size(lbuf,2).eq.size(gbuf,2),  &
4336                 '** global_max_3i: size(lbuf,2).ne.size(gbuf,2) ', &
4337                 size(lbuf,2), size(gbuf,2) )
4338     
4339            call assert( size(lbuf,3).eq.size(gbuf,3),  &
4340                 '** global_max_3i: size(lbuf,3).ne.size(gbuf,3) ', &
4341                 size(lbuf,3), size(gbuf,3) )
4342         endif
4343     
4344         do k=lbound(lbuf,3),ubound(lbuf,3)
4345            do j=lbound(lbuf,2),ubound(lbuf,2)
4346               call global_max_1i( lbuf(:,j,k), gbuf(:,j,k), lroot, lidebug )
4347            enddo
4348         enddo
4349     #else
4350         gbuf = lbuf
4351     #endif
4352     
4353         return
4354       end subroutine global_max_3i
4355     
4356       subroutine global_max_0r( lbuf, gbuf, mroot, idebug )
4357         real, intent(in) :: lbuf
4358         real, intent(out) :: gbuf
4359         integer, optional, intent(in) :: mroot, idebug
4360     
4361     #ifdef MPI
4362         integer :: recvtype, sendtype, sendcnt, ierr,lroot, lidebug
4363     
4364         if (.not. present(mroot)) then
4365            lroot = 0
4366         else
4367            lroot = mroot
4368         endif
4369     
4370         if (.not. present(idebug)) then
4371            lidebug = 0
4372         else
4373            lidebug = idebug
4374         endif
4375     
4376         recvtype = MPI_REAL
4377         sendtype = recvtype
4378     
4379         sendcnt = 1
4380     
4381         call MPI_Reduce( lbuf, gbuf, sendcnt, sendtype,  MPI_MAX, &
4382              lroot, MPI_COMM_WORLD, ierr )
4383         call MPI_Check( 'global_max_0r:MPI_Reduce', ierr )
4384     #else
4385         gbuf = lbuf
4386     #endif
4387     
4388         return
4389       end subroutine global_max_0r
4390     
4391     
4392       subroutine global_max_1r( lbuf, gbuf, mroot, idebug )
4393         real, intent(in), dimension(:) :: lbuf
4394         real, intent(out), dimension(:) :: gbuf
4395         integer, optional, intent(in) :: mroot, idebug
4396     
4397     #ifdef MPI
4398         integer :: recvtype, sendtype, sendcnt, ierr,lroot, lidebug
4399     
4400         if (.not. present(mroot)) then
4401            lroot = 0
4402         else
4403            lroot = mroot
4404         endif
4405     
4406         if (.not. present(idebug)) then
4407            lidebug = 0
4408         else
4409            lidebug = idebug
4410         endif
4411     
4412         recvtype = MPI_REAL
4413         sendtype = recvtype
4414     
4415         sendcnt = size(lbuf)
4416     
4417         call MPI_Reduce( lbuf, gbuf, sendcnt, sendtype,  MPI_MAX, &
4418              lroot, MPI_COMM_WORLD, ierr )
4419         call MPI_Check( 'global_max_1r:MPI_Reduce', ierr )
4420     #else
4421         gbuf = lbuf
4422     #endif
4423     
4424         return
4425       end subroutine global_max_1r
4426     
4427       subroutine global_max_2r( lbuf, gbuf, mroot, idebug )
4428         real, intent(in), dimension(:,:) :: lbuf
4429         real, intent(out), dimension(:,:) :: gbuf
4430         integer, optional, intent(in) :: mroot, idebug
4431     
4432     #ifdef MPI
4433         integer :: i,j,lroot, lidebug
4434     
4435         if (.not. present(mroot)) then
4436            lroot = 0
4437         else
4438            lroot = mroot
4439         endif
4440     
4441         if (.not. present(idebug)) then
4442            lidebug = 0
4443         else
4444            lidebug = idebug
4445         endif
4446     
4447         if(myPE.eq.lroot) then
4448            call assert( size(lbuf,2).eq.size(gbuf,2),  &
4449                 '** global_max_2r: size(lbuf,2).ne.size(gbuf,2) ', &
4450                 size(lbuf,2), size(gbuf,2) )
4451         endif
4452     
4453         do j=lbound(lbuf,2),ubound(lbuf,2)
4454            call global_max_1r( lbuf(:,j), gbuf(:,j), lroot, lidebug )
4455         enddo
4456     #else
4457         gbuf = lbuf
4458     #endif
4459     
4460         return
4461       end subroutine global_max_2r
4462     
4463       subroutine global_max_3r( lbuf, gbuf, mroot, idebug )
4464         real, intent(in), dimension(:,:,:) :: lbuf
4465         real, intent(out), dimension(:,:,:) :: gbuf
4466         integer, optional, intent(in) :: mroot, idebug
4467     
4468     #ifdef MPI
4469         integer :: j,k,lroot, lidebug
4470     
4471         if (.not. present(mroot)) then
4472            lroot = 0
4473         else
4474            lroot = mroot
4475         endif
4476     
4477         if (.not. present(idebug)) then
4478            lidebug = 0
4479         else
4480            lidebug = idebug
4481         endif
4482     
4483         if(myPE.eq.lroot) then
4484            call assert( size(lbuf,2).eq.size(gbuf,2),  &
4485                 '** global_max_3i: size(lbuf,2).ne.size(gbuf,2) ', &
4486                 size(lbuf,2), size(gbuf,2) )
4487     
4488            call assert( size(lbuf,3).eq.size(gbuf,3),  &
4489                 '** global_max_3i: size(lbuf,3).ne.size(gbuf,3) ', &
4490                 size(lbuf,3), size(gbuf,3) )
4491         endif
4492     
4493         do k=lbound(lbuf,3),ubound(lbuf,3)
4494            do j=lbound(lbuf,2),ubound(lbuf,2)
4495               call global_max_1r( lbuf(:,j,k), gbuf(:,j,k), lroot, lidebug )
4496            enddo
4497         enddo
4498     #else
4499         gbuf = lbuf
4500     #endif
4501     
4502         return
4503       end subroutine global_max_3r
4504     
4505       subroutine global_max_0d( lbuf, gbuf, mroot, idebug )
4506         double precision, intent(in) :: lbuf
4507         double precision, intent(out) :: gbuf
4508         integer, optional, intent(in) :: mroot, idebug
4509     
4510     #ifdef MPI
4511         integer :: recvtype, sendtype, sendcnt, ierr,lroot, lidebug
4512     
4513         if (.not. present(mroot)) then
4514            lroot = 0
4515         else
4516            lroot = mroot
4517         endif
4518     
4519         if (.not. present(idebug)) then
4520            lidebug = 0
4521         else
4522            lidebug = idebug
4523         endif
4524     
4525         recvtype = MPI_DOUBLE_PRECISION
4526         sendtype = recvtype
4527     
4528         sendcnt = 1
4529     
4530         call MPI_Reduce( lbuf, gbuf, sendcnt, sendtype,  MPI_MAX, &
4531              lroot, MPI_COMM_WORLD, ierr )
4532         call MPI_Check( 'global_max_0d:MPI_Reduce', ierr )
4533     #else
4534         gbuf = lbuf
4535     #endif
4536     
4537         return
4538       end subroutine global_max_0d
4539     
4540     
4541       subroutine global_max_1d( lbuf, gbuf, mroot, idebug )
4542         double precision, intent(in), dimension(:) :: lbuf
4543         double precision, intent(out), dimension(:) :: gbuf
4544         integer, optional, intent(in) :: mroot, idebug
4545     
4546     #ifdef MPI
4547         integer :: recvtype, sendtype, sendcnt, ierr,lroot, lidebug
4548     
4549         if (.not. present(mroot)) then
4550            lroot = 0
4551         else
4552            lroot = mroot
4553         endif
4554     
4555         if (.not. present(idebug)) then
4556            lidebug = 0
4557         else
4558            lidebug = idebug
4559         endif
4560     
4561         recvtype = MPI_DOUBLE_PRECISION
4562         sendtype = recvtype
4563     
4564         sendcnt = size(lbuf)
4565     
4566         call MPI_Reduce( lbuf, gbuf, sendcnt, sendtype,  MPI_MAX, &
4567              lroot, MPI_COMM_WORLD, ierr )
4568         call MPI_Check( 'global_max_1d:MPI_Reduce', ierr )
4569     #else
4570         gbuf = lbuf
4571     #endif
4572     
4573         return
4574       end subroutine global_max_1d
4575     
4576       subroutine global_max_2d( lbuf, gbuf, mroot, idebug )
4577         double precision, intent(in), dimension(:,:) :: lbuf
4578         double precision, intent(out), dimension(:,:) :: gbuf
4579         integer, optional, intent(in) :: mroot, idebug
4580     
4581     #ifdef MPI
4582         integer :: i,j,lroot, lidebug
4583     
4584         if (.not. present(mroot)) then
4585            lroot = 0
4586         else
4587            lroot = mroot
4588         endif
4589     
4590         if (.not. present(idebug)) then
4591            lidebug = 0
4592         else
4593            lidebug = idebug
4594         endif
4595     
4596         if(myPE.eq.lroot) then
4597            call assert( size(lbuf,2).eq.size(gbuf,2),  &
4598                 '** global_max_2d: size(lbuf,2).ne.size(gbuf,2) ', &
4599                 size(lbuf,2), size(gbuf,2) )
4600         endif
4601     
4602         do j=lbound(lbuf,2),ubound(lbuf,2)
4603            call global_max_1d( lbuf(:,j), gbuf(:,j), lroot, lidebug )
4604         enddo
4605     #else
4606         gbuf = lbuf
4607     #endif
4608     
4609         return
4610       end subroutine global_max_2d
4611     
4612       subroutine global_max_3d( lbuf, gbuf, mroot, idebug )
4613         double precision, intent(in), dimension(:,:,:) :: lbuf
4614         double precision, intent(out), dimension(:,:,:) :: gbuf
4615         integer, optional, intent(in) :: mroot, idebug
4616     
4617     #ifdef MPI
4618         integer :: j,k,lroot, lidebug
4619     
4620         if (.not. present(mroot)) then
4621            lroot = 0
4622         else
4623            lroot = mroot
4624         endif
4625     
4626         if (.not. present(idebug)) then
4627            lidebug = 0
4628         else
4629            lidebug = idebug
4630         endif
4631     
4632         if(myPE.eq.lroot) then
4633            call assert( size(lbuf,2).eq.size(gbuf,2),  &
4634                 '** global_max_3i: size(lbuf,2).ne.size(gbuf,2) ', &
4635                 size(lbuf,2), size(gbuf,2) )
4636     
4637            call assert( size(lbuf,3).eq.size(gbuf,3),  &
4638                 '** global_max_3i: size(lbuf,3).ne.size(gbuf,3) ', &
4639                 size(lbuf,3), size(gbuf,3) )
4640         endif
4641     
4642         do k=lbound(lbuf,3),ubound(lbuf,3)
4643            do j=lbound(lbuf,2),ubound(lbuf,2)
4644               call global_max_1d( lbuf(:,j,k), gbuf(:,j,k), lroot, lidebug )
4645            enddo
4646         enddo
4647     #else
4648         gbuf = lbuf
4649     #endif
4650     
4651         return
4652       end subroutine global_max_3d
4653     
4654       subroutine global_all_max_onevar_0d( gbuf )
4655         doubleprecision, intent(inout) :: gbuf
4656         doubleprecision :: lbuf
4657     
4658     #ifdef MPI
4659         lbuf = gbuf
4660         call global_all_max_0d( lbuf, gbuf )
4661     #endif
4662         return
4663       end subroutine global_all_max_onevar_0d
4664     
4665     
4666       subroutine global_all_max_onevar_1d( gbuf )
4667         doubleprecision, dimension(:), intent(inout) :: gbuf
4668         doubleprecision, dimension(size(gbuf)) :: lbuf
4669     
4670     #ifdef MPI
4671         lbuf = gbuf
4672         call global_all_max_1d( lbuf, gbuf )
4673     #endif
4674         return
4675       end subroutine global_all_max_onevar_1d
4676     
4677       subroutine global_all_max_onevar_2d( gbuf )
4678         doubleprecision, dimension(:,:), intent(inout) :: gbuf
4679         doubleprecision, dimension(size(gbuf,1),size(gbuf,2)) :: lbuf
4680     
4681     #ifdef MPI
4682         lbuf = gbuf
4683         call global_all_max_2d( lbuf, gbuf )
4684     #endif
4685         return
4686       end subroutine global_all_max_onevar_2d
4687     
4688     
4689       subroutine global_all_max_onevar_3d( gbuf )
4690         doubleprecision, dimension(:,:,:), intent(inout) :: gbuf
4691         doubleprecision, dimension(size(gbuf,1),size(gbuf,2),size(gbuf,3)) :: lbuf
4692     
4693     #ifdef MPI
4694         lbuf = gbuf
4695         call global_all_max_3d( lbuf, gbuf )
4696     #endif
4697         return
4698       end subroutine global_all_max_onevar_3d
4699     
4700     
4701     
4702     
4703       subroutine global_all_max_onevar_0i( gbuf )
4704         integer, intent(inout) :: gbuf
4705         integer :: lbuf
4706     
4707     #ifdef MPI
4708         lbuf = gbuf
4709         call global_all_max_0i( lbuf, gbuf )
4710     #endif
4711         return
4712       end subroutine global_all_max_onevar_0i
4713     
4714     
4715       subroutine global_all_max_onevar_1i( gbuf )
4716         integer, dimension(:), intent(inout) :: gbuf
4717         integer, dimension(size(gbuf)) :: lbuf
4718     
4719     #ifdef MPI
4720         lbuf = gbuf
4721         call global_all_max_1i( lbuf, gbuf )
4722     #endif
4723         return
4724       end subroutine global_all_max_onevar_1i
4725     
4726       subroutine global_all_max_onevar_2i( gbuf )
4727         integer, dimension(:,:), intent(inout) :: gbuf
4728         integer, dimension(size(gbuf,1),size(gbuf,2)) :: lbuf
4729     
4730     #ifdef MPI
4731         lbuf = gbuf
4732         call global_all_max_2i( lbuf, gbuf )
4733     #endif
4734         return
4735       end subroutine global_all_max_onevar_2i
4736     
4737     
4738       subroutine global_all_max_onevar_3i( gbuf )
4739         integer, dimension(:,:,:), intent(inout) :: gbuf
4740         integer, dimension(size(gbuf,1),size(gbuf,2),size(gbuf,3)) :: lbuf
4741     
4742     #ifdef MPI
4743         lbuf = gbuf
4744         call global_all_max_3i( lbuf, gbuf )
4745     #endif
4746         return
4747       end subroutine global_all_max_onevar_3i
4748     
4749       subroutine global_all_max_onevar_0r( gbuf )
4750         real, intent(inout) :: gbuf
4751         real :: lbuf
4752     
4753     #ifdef MPI
4754         lbuf = gbuf
4755         call global_all_max_0r( lbuf, gbuf )
4756     #endif
4757         return
4758       end subroutine global_all_max_onevar_0r
4759     
4760       subroutine global_all_max_onevar_1r( gbuf )
4761         real, dimension(:), intent(inout) :: gbuf
4762         real, dimension(size(gbuf)) :: lbuf
4763     
4764     #ifdef MPI
4765         lbuf = gbuf
4766         call global_all_max_1r( lbuf, gbuf )
4767     #endif
4768         return
4769       end subroutine global_all_max_onevar_1r
4770     
4771       subroutine global_all_max_onevar_2r( gbuf )
4772         real, dimension(:,:), intent(inout) :: gbuf
4773         real, dimension(size(gbuf,1),size(gbuf,2)) :: lbuf
4774     
4775     #ifdef MPI
4776         lbuf = gbuf
4777         call global_all_max_2r( lbuf, gbuf )
4778     #endif
4779         return
4780       end subroutine global_all_max_onevar_2r
4781     
4782     
4783       subroutine global_all_max_onevar_3r( gbuf )
4784         real, dimension(:,:,:), intent(inout) :: gbuf
4785         real, dimension(size(gbuf,1),size(gbuf,2),size(gbuf,3)) :: lbuf
4786     
4787     #ifdef MPI
4788         lbuf = gbuf
4789         call global_all_max_3r( lbuf, gbuf )
4790     #endif
4791         return
4792       end subroutine global_all_max_onevar_3r
4793     
4794     
4795     
4796       subroutine global_all_max_0i( lbuf, gbuf, mroot, idebug )
4797         integer, intent(in) :: lbuf
4798         integer, intent(out) :: gbuf
4799         integer, optional, intent(in) :: mroot, idebug
4800     
4801     #ifdef MPI
4802         integer :: recvtype, sendtype, sendcnt, ierr,lroot, lidebug
4803     
4804         if (.not. present(mroot)) then
4805            lroot = 0
4806         else
4807            lroot = mroot
4808         endif
4809     
4810         if (.not. present(idebug)) then
4811            lidebug = 0
4812         else
4813            lidebug = idebug
4814         endif
4815     
4816         recvtype = MPI_INTEGER
4817         sendtype = recvtype
4818     
4819         sendcnt = 1
4820     
4821         call MPI_Allreduce( lbuf, gbuf, sendcnt, sendtype, MPI_MAX, &
4822              MPI_COMM_WORLD, ierr )
4823         call MPI_Check( 'global_all_max_0i:MPI_Allreduce', ierr )
4824     #else
4825         gbuf = lbuf
4826     #endif
4827     
4828         return
4829       end subroutine global_all_max_0i
4830     
4831     
4832       subroutine global_all_max_1i( lbuf, gbuf, mroot, idebug )
4833         integer, intent(in), dimension(:) :: lbuf
4834         integer, intent(out), dimension(:) :: gbuf
4835         integer, optional, intent(in) :: mroot, idebug
4836     
4837     #ifdef MPI
4838         integer :: recvtype, sendtype, sendcnt, ierr,lroot, lidebug
4839     
4840         if (.not. present(mroot)) then
4841            lroot = 0
4842         else
4843            lroot = mroot
4844         endif
4845     
4846         if (.not. present(idebug)) then
4847            lidebug = 0
4848         else
4849            lidebug = idebug
4850         endif
4851     
4852         recvtype = MPI_INTEGER
4853         sendtype = recvtype
4854     
4855         sendcnt = size(lbuf)
4856     
4857         call MPI_Allreduce( lbuf, gbuf, sendcnt, sendtype, MPI_MAX, &
4858              MPI_COMM_WORLD, ierr )
4859         call MPI_Check( 'global_all_max_1i:MPI_Allreduce', ierr )
4860     #else
4861         gbuf = lbuf
4862     #endif
4863     
4864         return
4865       end subroutine global_all_max_1i
4866     
4867       subroutine global_all_max_2i( lbuf, gbuf, mroot, idebug )
4868         integer, intent(in), dimension(:,:) :: lbuf
4869         integer, intent(out), dimension(:,:) :: gbuf
4870         integer, optional, intent(in) :: mroot, idebug
4871     
4872     #ifdef MPI
4873         integer :: i,j,lroot, lidebug
4874     
4875         if (.not. present(mroot)) then
4876            lroot = 0
4877         else
4878            lroot = mroot
4879         endif
4880     
4881         if (.not. present(idebug)) then
4882            lidebug = 0
4883         else
4884            lidebug = idebug
4885         endif
4886     
4887         call assert( size(lbuf,2).eq.size(gbuf,2),  &
4888              '** global_all_max_2i: size(lbuf,2).ne.size(gbuf,2) ', &
4889              size(lbuf,2), size(gbuf,2) )
4890     
4891         do j=lbound(lbuf,2),ubound(lbuf,2)
4892            call global_all_max_1i( lbuf(:,j), gbuf(:,j), lroot, lidebug )
4893         enddo
4894     #else
4895         gbuf = lbuf
4896     #endif
4897     
4898         return
4899       end subroutine global_all_max_2i
4900     
4901       subroutine global_all_max_3i( lbuf, gbuf, mroot, idebug )
4902         integer, intent(in), dimension(:,:,:) :: lbuf
4903         integer, intent(out), dimension(:,:,:) :: gbuf
4904         integer, optional, intent(in) :: mroot, idebug
4905     
4906     #ifdef MPI
4907         integer :: j,k,lroot, lidebug
4908     
4909         if (.not. present(mroot)) then
4910            lroot = 0
4911         else
4912            lroot = mroot
4913         endif
4914     
4915         if (.not. present(idebug)) then
4916            lidebug = 0
4917         else
4918            lidebug = idebug
4919         endif
4920     
4921         call assert( size(lbuf,2).eq.size(gbuf,2),  &
4922              '** global_all_max_3i: size(lbuf,2).ne.size(gbuf,2) ', &
4923              size(lbuf,2), size(gbuf,2) )
4924     
4925         call assert( size(lbuf,3).eq.size(gbuf,3),  &
4926              '** global_all_max_3i: size(lbuf,3).ne.size(gbuf,3) ', &
4927              size(lbuf,3), size(gbuf,3) )
4928     
4929         do k=lbound(lbuf,3),ubound(lbuf,3)
4930            do j=lbound(lbuf,2),ubound(lbuf,2)
4931               call global_all_max_1i( lbuf(:,j,k), gbuf(:,j,k), lroot, lidebug )
4932            enddo
4933         enddo
4934     #else
4935         gbuf = lbuf
4936     #endif
4937     
4938         return
4939       end subroutine global_all_max_3i
4940     
4941       subroutine global_all_max_0r( lbuf, gbuf, mroot, idebug )
4942         real, intent(in) :: lbuf
4943         real, intent(out) :: gbuf
4944         integer, optional, intent(in) :: mroot, idebug
4945     
4946     #ifdef MPI
4947         integer :: recvtype, sendtype, sendcnt, ierr,lroot, lidebug
4948     
4949         if (.not. present(mroot)) then
4950            lroot = 0
4951         else
4952            lroot = mroot
4953         endif
4954     
4955         if (.not. present(idebug)) then
4956            lidebug = 0
4957         else
4958            lidebug = idebug
4959         endif
4960     
4961         recvtype = MPI_REAL
4962         sendtype = recvtype
4963     
4964         sendcnt = 1
4965     
4966         call MPI_Allreduce( lbuf, gbuf, sendcnt, sendtype, MPI_MAX, &
4967              MPI_COMM_WORLD, ierr )
4968         call MPI_Check( 'global_all_max_0r:MPI_Allreduce', ierr )
4969     #else
4970         gbuf = lbuf
4971     #endif
4972     
4973         return
4974       end subroutine global_all_max_0r
4975     
4976     
4977       subroutine global_all_max_1r( lbuf, gbuf, mroot, idebug )
4978         real, intent(in), dimension(:) :: lbuf
4979         real, intent(out), dimension(:) :: gbuf
4980         integer, optional, intent(in) :: mroot, idebug
4981     
4982     #ifdef MPI
4983         integer :: recvtype, sendtype, sendcnt, ierr,lroot, lidebug
4984     
4985         if (.not. present(mroot)) then
4986            lroot = 0
4987         else
4988            lroot = mroot
4989         endif
4990     
4991         if (.not. present(idebug)) then
4992            lidebug = 0
4993         else
4994            lidebug = idebug
4995         endif
4996     
4997         recvtype = MPI_REAL
4998         sendtype = recvtype
4999     
5000         sendcnt = size(lbuf)
5001     
5002         call MPI_Allreduce( lbuf, gbuf, sendcnt, sendtype, MPI_MAX, &
5003              MPI_COMM_WORLD, ierr )
5004         call MPI_Check( 'global_all_max_1r:MPI_Allreduce', ierr )
5005     #else
5006         gbuf = lbuf
5007     #endif
5008     
5009         return
5010       end subroutine global_all_max_1r
5011     
5012       subroutine global_all_max_2r( lbuf, gbuf, mroot, idebug )
5013         real, intent(in), dimension(:,:) :: lbuf
5014         real, intent(out), dimension(:,:) :: gbuf
5015         integer, optional, intent(in) :: mroot, idebug
5016     
5017     #ifdef MPI
5018         integer :: i,j,lroot, lidebug
5019     
5020         if (.not. present(mroot)) then
5021            lroot = 0
5022         else
5023            lroot = mroot
5024         endif
5025     
5026         if (.not. present(idebug)) then
5027            lidebug = 0
5028         else
5029            lidebug = idebug
5030         endif
5031     
5032         call assert( size(lbuf,2).eq.size(gbuf,2),  &
5033              '** global_all_max_2r: size(lbuf,2).ne.size(gbuf,2) ', &
5034              size(lbuf,2), size(gbuf,2) )
5035     
5036         do j=lbound(lbuf,2),ubound(lbuf,2)
5037            call global_all_max_1r( lbuf(:,j), gbuf(:,j), lroot, lidebug )
5038         enddo
5039     #else
5040         gbuf = lbuf
5041     #endif
5042     
5043         return
5044       end subroutine global_all_max_2r
5045     
5046       subroutine global_all_max_3r( lbuf, gbuf, mroot, idebug )
5047         real, intent(in), dimension(:,:,:) :: lbuf
5048         real, intent(out), dimension(:,:,:) :: gbuf
5049         integer, optional, intent(in) :: mroot, idebug
5050     
5051     #ifdef MPI
5052         integer :: j,k,lroot, lidebug
5053     
5054         if (.not. present(mroot)) then
5055            lroot = 0
5056         else
5057            lroot = mroot
5058         endif
5059     
5060         if (.not. present(idebug)) then
5061            lidebug = 0
5062         else
5063            lidebug = idebug
5064         endif
5065     
5066         call assert( size(lbuf,2).eq.size(gbuf,2),  &
5067              '** global_all_max_3i: size(lbuf,2).ne.size(gbuf,2) ', &
5068              size(lbuf,2), size(gbuf,2) )
5069     
5070         call assert( size(lbuf,3).eq.size(gbuf,3),  &
5071              '** global_all_max_3i: size(lbuf,3).ne.size(gbuf,3) ', &
5072              size(lbuf,3), size(gbuf,3) )
5073     
5074         do k=lbound(lbuf,3),ubound(lbuf,3)
5075            do j=lbound(lbuf,2),ubound(lbuf,2)
5076               call global_all_max_1r( lbuf(:,j,k), gbuf(:,j,k), lroot, lidebug )
5077            enddo
5078         enddo
5079     #else
5080         gbuf = lbuf
5081     #endif
5082     
5083         return
5084       end subroutine global_all_max_3r
5085     
5086       subroutine global_all_max_0d( lbuf, gbuf, mroot, idebug )
5087         double precision, intent(in) :: lbuf
5088         double precision, intent(out) :: gbuf
5089         integer, optional, intent(in) :: mroot, idebug
5090     
5091     #ifdef MPI
5092         integer :: recvtype, sendtype, sendcnt, ierr,lroot, lidebug
5093     
5094         if (.not. present(mroot)) then
5095            lroot = 0
5096         else
5097            lroot = mroot
5098         endif
5099     
5100         if (.not. present(idebug)) then
5101            lidebug = 0
5102         else
5103            lidebug = idebug
5104         endif
5105     
5106         recvtype = MPI_DOUBLE_PRECISION
5107         sendtype = recvtype
5108     
5109         sendcnt = 1
5110     
5111         call MPI_Allreduce( lbuf, gbuf, sendcnt, sendtype, MPI_MAX, &
5112              MPI_COMM_WORLD, ierr )
5113         call MPI_Check( 'global_all_max_0d:MPI_Allreduce', ierr )
5114     #else
5115         gbuf = lbuf
5116     #endif
5117     
5118         return
5119       end subroutine global_all_max_0d
5120     
5121     
5122       subroutine global_all_max_1d( lbuf, gbuf, mroot, idebug )
5123         double precision, intent(in), dimension(:) :: lbuf
5124         double precision, intent(out), dimension(:) :: gbuf
5125         integer, optional, intent(in) :: mroot, idebug
5126     
5127     #ifdef MPI
5128         integer :: recvtype, sendtype, sendcnt, ierr,lroot, lidebug
5129     
5130         if (.not. present(mroot)) then
5131            lroot = 0
5132         else
5133            lroot = mroot
5134         endif
5135     
5136         if (.not. present(idebug)) then
5137            lidebug = 0
5138         else
5139            lidebug = idebug
5140         endif
5141     
5142         recvtype = MPI_DOUBLE_PRECISION
5143         sendtype = recvtype
5144     
5145         sendcnt = size(lbuf)
5146     
5147         call MPI_Allreduce( lbuf, gbuf, sendcnt, sendtype,  MPI_MAX, &
5148              MPI_COMM_WORLD, ierr )
5149         call MPI_Check( 'global_all_max_1d:MPI_Allreduce', ierr )
5150     #else
5151         gbuf = lbuf
5152     #endif
5153     
5154         return
5155       end subroutine global_all_max_1d
5156     
5157       subroutine global_all_max_2d( lbuf, gbuf, mroot, idebug )
5158         double precision, intent(in), dimension(:,:) :: lbuf
5159         double precision, intent(out), dimension(:,:) :: gbuf
5160         integer, optional, intent(in) :: mroot, idebug
5161     
5162     #ifdef MPI
5163         integer :: i,j,lroot, lidebug
5164     
5165         if (.not. present(mroot)) then
5166            lroot = 0
5167         else
5168            lroot = mroot
5169         endif
5170     
5171         if (.not. present(idebug)) then
5172            lidebug = 0
5173         else
5174            lidebug = idebug
5175         endif
5176     
5177         call assert( size(lbuf,2).eq.size(gbuf,2),  &
5178              '** global_all_max_2d: size(lbuf,2).ne.size(gbuf,2) ', &
5179              size(lbuf,2), size(gbuf,2) )
5180     
5181         do j=lbound(lbuf,2),ubound(lbuf,2)
5182            call global_all_max_1d( lbuf(:,j), gbuf(:,j), lroot, lidebug )
5183         enddo
5184     #else
5185         gbuf = lbuf
5186     #endif
5187     
5188         return
5189       end subroutine global_all_max_2d
5190     
5191       subroutine global_all_max_3d( lbuf, gbuf, mroot, idebug )
5192         double precision, intent(in), dimension(:,:,:) :: lbuf
5193         double precision, intent(out), dimension(:,:,:) :: gbuf
5194         integer, optional, intent(in) :: mroot, idebug
5195     
5196     #ifdef MPI
5197         integer :: j,k,lroot, lidebug
5198     
5199         if (.not. present(mroot)) then
5200            lroot = 0
5201         else
5202            lroot = mroot
5203         endif
5204     
5205         if (.not. present(idebug)) then
5206            lidebug = 0
5207         else
5208            lidebug = idebug
5209         endif
5210     
5211         call assert( size(lbuf,2).eq.size(gbuf,2),  &
5212              '** global_all_max_3i: size(lbuf,2).ne.size(gbuf,2) ', &
5213              size(lbuf,2), size(gbuf,2) )
5214     
5215         call assert( size(lbuf,3).eq.size(gbuf,3),  &
5216              '** global_all_max_3i: size(lbuf,3).ne.size(gbuf,3) ', &
5217              size(lbuf,3), size(gbuf,3) )
5218     
5219         do k=lbound(lbuf,3),ubound(lbuf,3)
5220            do j=lbound(lbuf,2),ubound(lbuf,2)
5221               call global_all_max_1d( lbuf(:,j,k), gbuf(:,j,k), lroot, lidebug )
5222            enddo
5223         enddo
5224     #else
5225         gbuf = lbuf
5226     #endif
5227     
5228         return
5229       end subroutine global_all_max_3d
5230     
5231     
5232     
5233       subroutine global_all_or_onevar_0d( gbuf )
5234         logical, intent(inout) :: gbuf
5235         logical :: lbuf
5236     
5237     #ifdef MPI
5238         lbuf = gbuf
5239         call global_all_or_0d( lbuf, gbuf )
5240     #endif
5241         return
5242       end subroutine global_all_or_onevar_0d
5243     
5244       subroutine global_all_or_onevar_1d( gbuf )
5245         logical, dimension(:), intent(inout) :: gbuf
5246         logical, dimension(size(gbuf)) :: lbuf
5247     
5248     #ifdef MPI
5249         lbuf = gbuf
5250         call global_all_or_1d( lbuf, gbuf )
5251     #endif
5252         return
5253       end subroutine global_all_or_onevar_1d
5254     
5255       subroutine global_all_and_onevar_0d( gbuf )
5256         logical, intent(inout) :: gbuf
5257         logical :: lbuf
5258     
5259     #ifdef MPI
5260         lbuf = gbuf
5261         call global_all_and_0d( lbuf, gbuf )
5262     #endif
5263         return
5264       end subroutine global_all_and_onevar_0d
5265     
5266       subroutine global_all_and_onevar_1d( gbuf )
5267         logical, dimension(:), intent(inout) :: gbuf
5268         logical, dimension(size(gbuf)) :: lbuf
5269     
5270     #ifdef MPI
5271         lbuf = gbuf
5272         call global_all_and_1d( lbuf, gbuf )
5273     #endif
5274         return
5275       end subroutine global_all_and_onevar_1d
5276     
5277     
5278       subroutine global_all_and_0d( lvalue, gvalue, mroot, idebug )
5279         logical, intent(in) :: lvalue
5280         logical, intent(out) :: gvalue
5281         integer, optional, intent(in) :: mroot, idebug
5282     
5283     #ifdef MPI
5284         !       ---------------
5285         !       local variables
5286         !       ---------------
5287         integer :: ierror, icount
5288         integer :: lroot, lidebug
5289     
5290         if (.not. present(mroot)) then
5291            lroot = 0
5292         else
5293            lroot = mroot
5294         endif
5295     
5296         if (.not. present(idebug)) then
5297            lidebug = 0
5298         else
5299            lidebug = idebug
5300         endif
5301     
5302         icount = 1
5303     
5304         call MPI_Allreduce( lvalue, gvalue, icount, MPI_LOGICAL, &
5305              MPI_LAND, MPI_COMM_WORLD, ierror )
5306     
5307         call MPI_Check( 'global_all_and_0d ', ierror )
5308     #else
5309         gvalue = lvalue
5310     #endif
5311         return
5312       end subroutine  global_all_and_0d
5313     
5314     
5315       subroutine global_all_and_1d( lvalue, gvalue, mroot, idebug )
5316         logical, intent(in), dimension(:) :: lvalue
5317         logical, intent(out), dimension(:) :: gvalue
5318         integer, optional, intent(in) :: mroot, idebug
5319     
5320     #ifdef MPI
5321         !       ---------------
5322         !       local variables
5323         !       ---------------
5324         integer :: ierror, icount
5325         integer :: lroot, lidebug
5326     
5327         if (.not. present(mroot)) then
5328            lroot = 0
5329         else
5330            lroot = mroot
5331         endif
5332     
5333         if (.not. present(idebug)) then
5334            lidebug = 0
5335         else
5336            lidebug = idebug
5337         endif
5338     
5339     
5340         icount = size( lvalue )
5341     
5342         call MPI_Allreduce( lvalue, gvalue, icount, MPI_LOGICAL, &
5343              MPI_LAND, MPI_COMM_WORLD, ierror )
5344     
5345         call MPI_Check( 'global_all_and_1d ', ierror )
5346     #else
5347         gvalue = lvalue
5348     #endif
5349         return
5350       end subroutine global_all_and_1d
5351     
5352     
5353     
5354       subroutine global_all_or_0d( lvalue, gvalue, mroot, idebug )
5355         logical, intent(in) :: lvalue
5356         logical, intent(out) :: gvalue
5357         integer, optional, intent(in) :: mroot, idebug
5358     
5359     #ifdef MPI
5360         !       ---------------
5361         !       local variables
5362         !       ---------------
5363         integer :: ierror, icount
5364         integer :: lroot, lidebug
5365     
5366         if (.not. present(mroot)) then
5367            lroot = 0
5368         else
5369            lroot = mroot
5370         endif
5371     
5372         if (.not. present(idebug)) then
5373            lidebug = 0
5374         else
5375            lidebug = idebug
5376         endif
5377     
5378     
5379         icount = 1
5380     
5381         call MPI_Allreduce( lvalue, gvalue, icount, MPI_LOGICAL, &
5382              MPI_LOR, MPI_COMM_WORLD, ierror )
5383     
5384         call MPI_Check( 'global_all_or_0d ', ierror )
5385     #else
5386         gvalue = lvalue
5387     #endif
5388         return
5389       end subroutine global_all_or_0d
5390     
5391     
5392       subroutine global_all_or_1d( lvalue, gvalue, mroot, idebug )
5393         logical, intent(in), dimension(:) :: lvalue
5394         logical, intent(out), dimension(:) :: gvalue
5395         integer, optional, intent(in) :: mroot, idebug
5396     
5397     #ifdef MPI
5398         !       ---------------
5399         !       local variables
5400         !       ---------------
5401         integer :: ierror, icount
5402         integer :: lroot, lidebug
5403     
5404         if (.not. present(mroot)) then
5405            lroot = 0
5406         else
5407            lroot = mroot
5408         endif
5409     
5410         if (.not. present(idebug)) then
5411            lidebug = 0
5412         else
5413            lidebug = idebug
5414         endif
5415     
5416     
5417         icount = size( lvalue )
5418     
5419         call MPI_Allreduce( lvalue, gvalue, icount, MPI_LOGICAL, &
5420              MPI_LOR, MPI_COMM_WORLD, ierror )
5421     
5422         call MPI_Check( 'global_all_or_1d ', ierror )
5423     #else
5424         gvalue = lvalue
5425     #endif
5426         return
5427       end subroutine global_all_or_1d
5428     
5429     
5430       !``````````````````````````````````````````````````````````````````````!
5431       ! Subroutine: ExitMPI                                                  !
5432       !                                                                      !
5433       ! Purpose: Clean abort from a parallel run. This routine is invoked by !
5434       ! by calling MFIX_EXIT                                                 !
5435       !......................................................................!
5436       SUBROUTINE ExitMPI(myid)
5437     
5438         USE funits, only: UNIT_LOG
5439         USE funits, only: DMP_LOG
5440     
5441         implicit none
5442     
5443         INTEGER, optional, intent(in) :: MyID
5444     
5445     #ifdef MPI
5446         INTEGER :: MyID_l
5447         INTEGER :: ERRORCODE
5448     
5449         ! Flag to call MPI_ABORT and bypass the call to MPI_Finalize.
5450         ! This is only needed if debugging a 'deadlocked' run.
5451         LOGICAL, PARAMETER :: FORCED_ABORT = .FALSE.
5452     
5453         ! Process ID (myPE) converted to a character string.
5454         CHARACTER(len=64) :: myID_c
5455     
5456     
5457         ! Set the ID of the caller.
5458         myID_l= merge(MyID, myPE, PRESENT(MyID))
5459         myID_c=''; WRITE(myID_c,*) myID_l
5460     
5461         ! Hard abort. If you need this functionality, then you need to figure
5462         ! out why the code has deadlocked. Most likely, a call to MFIX_EXIT
5463         ! was put inside of a logical branch that only a few ranks execute.
5464         ! DON'T JUST USE A FORCED ABORT --> FIX THE CODE CAUSE DEADLOCK <--
5465         IF(FORCED_ABORT) THEN
5466            ERRORCODE = 100 + myPE
5467            CALL MPI_ABORT(MPI_COMM_WORLD, ERRORCODE, MPIERR)
5468            WRITE(*,2000) myID_c, MPIERR
5469     
5470            ! Calls to ExitMPI (via MFIX_EXIT) should be made by all processes
5471            ! and therefore calling MPI_Finalize should be sufficient to exit
5472            ! a failed run. However, a FORCED_ABORT can be issued if deadlock
5473            ! is an issue.
5474         ELSE
5475            CALL MPI_BARRIER(MPI_COMM_WORLD, MPIERR)
5476            CALL MPI_Finalize(MPIERR)
5477         ENDIF
5478     
5479         ! Notify that MPI was cleanly terminated. This point will not be
5480         ! reached if MPI is aborted.
5481         IF(myPE == PE_IO) WRITE(*,1000)
5482     #endif
5483     
5484         RETURN
5485     
5486     1000 FORMAT(2/,1X,'MPI Terminated.')
5487     2000 FORMAT(2/,1X,'Rank ',A,' :: MPI_ABORT CODE = ',I4)
5488     
5489       END SUBROUTINE exitMPI
5490     
5491     end module mpi_utility
5492