<p><b>duda</b> 2010-02-24 13:40:38 -0700 (Wed, 24 Feb 2010)</p><p>Add additional global reductions from ocean branch (sum, min, <br>
and max for real and integer arrays).<br>
<br>
M    module_dmpar.F<br>
</p><hr noshade><pre><font color="gray">Modified: branches/hyd_model/src/module_dmpar.F
===================================================================
--- branches/hyd_model/src/module_dmpar.F        2010-02-10 23:58:45 UTC (rev 117)
+++ branches/hyd_model/src/module_dmpar.F        2010-02-24 20:40:38 UTC (rev 118)
@@ -242,6 +242,221 @@
    end subroutine dmpar_sum_int
 
 
+   subroutine dmpar_sum_real(dminfo, r, rsum)
+
+      implicit none
+
+      type (dm_info), intent(in) :: dminfo
+      real(kind=RKIND), intent(in) :: r
+      real(kind=RKIND), intent(out) :: rsum
+
+      integer :: mpi_ierr
+
+#ifdef _MPI
+      call MPI_Allreduce(r, rsum, 1, MPI_REALKIND, MPI_SUM, dminfo % comm, mpi_ierr)
+#else
+      rsum = r
+#endif
+
+   end subroutine dmpar_sum_real
+
+
+   subroutine dmpar_min_int(dminfo, i, imin)
+
+      implicit none
+
+      type (dm_info), intent(in) :: dminfo
+      integer, intent(in) :: i
+      integer, intent(out) :: imin
+
+      integer :: mpi_ierr
+
+#ifdef _MPI
+      call MPI_Allreduce(i, imin, 1, MPI_INTEGER, MPI_MIN, dminfo % comm, mpi_ierr)
+#else
+      imin = i
+#endif
+
+   end subroutine dmpar_min_int
+
+
+   subroutine dmpar_min_real(dminfo, r, rmin)
+
+      implicit none
+
+      type (dm_info), intent(in) :: dminfo
+      real(kind=RKIND), intent(in) :: r
+      real(kind=RKIND), intent(out) :: rmin
+
+      integer :: mpi_ierr
+
+#ifdef _MPI
+      call MPI_Allreduce(r, rmin, 1, MPI_REALKIND, MPI_MIN, dminfo % comm, mpi_ierr)
+#else
+      rmin = r
+#endif
+
+   end subroutine dmpar_min_real
+
+
+   subroutine dmpar_max_int(dminfo, i, imax)
+
+      implicit none
+
+      type (dm_info), intent(in) :: dminfo
+      integer, intent(in) :: i
+      integer, intent(out) :: imax
+      
+      integer :: mpi_ierr 
+      
+#ifdef _MPI
+      call MPI_Allreduce(i, imax, 1, MPI_INTEGER, MPI_MAX, dminfo % comm, mpi_ierr)
+#else
+      imax = i
+#endif
+
+   end subroutine dmpar_max_int
+
+
+   subroutine dmpar_max_real(dminfo, r, rmax)
+
+      implicit none
+
+      type (dm_info), intent(in) :: dminfo
+      real(kind=RKIND), intent(in) :: r
+      real(kind=RKIND), intent(out) :: rmax
+
+      integer :: mpi_ierr
+
+#ifdef _MPI
+      call MPI_Allreduce(r, rmax, 1, MPI_REALKIND, MPI_MAX, dminfo % comm, mpi_ierr)
+#else
+      rmax = r
+#endif
+
+   end subroutine dmpar_max_real
+
+
+   subroutine dmpar_sum_int_array(dminfo, nElements, inArray, outArray)
+
+      implicit none
+   
+      type (dm_info), intent(in) :: dminfo
+      integer, intent(in) :: nElements
+      integer, dimension(nElements), intent(in) :: inArray
+      integer, dimension(nElements), intent(out) :: outArray
+      
+      integer :: mpi_ierr
+
+#ifdef _MPI
+      call MPI_Allreduce(inArray, outArray, nElements, MPI_INTEGER, MPI_SUM, dminfo % comm, mpi_ierr)
+#else
+      outArray = inArray
+#endif
+
+   end subroutine dmpar_sum_int_array
+
+
+   subroutine dmpar_min_int_array(dminfo, nElements, inArray, outArray)
+   
+      implicit none
+      
+      type (dm_info), intent(in) :: dminfo
+      integer, intent(in) :: nElements
+      integer, dimension(nElements), intent(in) :: inArray
+      integer, dimension(nElements), intent(out) :: outArray
+
+      integer :: mpi_ierr
+
+#ifdef _MPI
+      call MPI_Allreduce(inArray, outArray, nElements, MPI_INTEGER, MPI_MIN, dminfo % comm, mpi_ierr)
+#else
+      outArray = inArray
+#endif
+
+   end subroutine dmpar_min_int_array
+
+
+   subroutine dmpar_max_int_array(dminfo, nElements, inArray, outArray)
+
+      implicit none
+
+      type (dm_info), intent(in) :: dminfo
+      integer, intent(in) :: nElements
+      integer, dimension(nElements), intent(in) :: inArray
+      integer, dimension(nElements), intent(out) :: outArray
+
+      integer :: mpi_ierr
+
+#ifdef _MPI
+      call MPI_Allreduce(inArray, outArray, nElements, MPI_INTEGER, MPI_MAX, dminfo % comm, mpi_ierr)
+#else
+      outArray = inArray
+#endif
+
+   end subroutine dmpar_max_int_array
+
+
+   subroutine dmpar_sum_real_array(dminfo, nElements, inArray, outArray)
+
+      implicit none
+
+      type (dm_info), intent(in) :: dminfo
+      integer, intent(in) :: nElements
+      real(kind=RKIND), dimension(nElements), intent(in) :: inArray
+      real(kind=RKIND), dimension(nElements), intent(out) :: outArray
+
+      integer :: mpi_ierr
+
+#ifdef _MPI
+      call MPI_Allreduce(inArray, outArray, nElements, MPI_REALKIND, MPI_SUM, dminfo % comm, mpi_ierr)
+#else
+      outArray = inArray
+#endif
+
+   end subroutine dmpar_sum_real_array
+
+
+   subroutine dmpar_min_real_array(dminfo, nElements, inArray, outArray)
+
+      implicit none
+
+      type (dm_info), intent(in) :: dminfo
+      integer, intent(in) :: nElements
+      real(kind=RKIND), dimension(nElements), intent(in) :: inArray
+      real(kind=RKIND), dimension(nElements), intent(out) :: outArray
+
+      integer :: mpi_ierr
+
+#ifdef _MPI
+      call MPI_Allreduce(inArray, outArray, nElements, MPI_REALKIND, MPI_MIN, dminfo % comm, mpi_ierr)
+#else
+      outArray = inArray
+#endif
+
+   end subroutine dmpar_min_real_array
+
+
+   subroutine dmpar_max_real_array(dminfo, nElements, inArray, outArray)
+
+      implicit none
+
+      type (dm_info), intent(in) :: dminfo
+      integer, intent(in) :: nElements
+      real(kind=RKIND), dimension(nElements), intent(in) :: inArray
+      real(kind=RKIND), dimension(nElements), intent(out) :: outArray
+
+      integer :: mpi_ierr
+
+#ifdef _MPI
+      call MPI_Allreduce(inArray, outArray, nElements, MPI_REALKIND, MPI_MAX, dminfo % comm, mpi_ierr)
+#else
+      outArray = inArray
+#endif
+
+   end subroutine dmpar_max_real_array
+
+
    subroutine dmpar_scatter_ints(dminfo, nprocs, noutlist, displs, counts, inlist, outlist)
 
       implicit none

</font>
</pre>