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