SUBROUTINE global_min_max( i1, i2, j1, j2, k1, k2, ar, mode, value, & value_ijk, value1, value1_ijk ) !-------------------------------------------------------------------------------! ! Current revisions: ! ----------------- ! ! ! Former revisions: ! ----------------- ! $Id: global_min_max.f90 484 2010-02-05 07:36:54Z maronga $ ! RCS Log replace by Id keyword, revision history cleaned up ! ! Revision 1.11 2003/04/16 12:56:58 raasch ! Index values of the extrema are limited to the range 0..nx, 0..ny ! ! Revision 1.1 1997/07/24 11:14:03 raasch ! Initial revision ! ! ! Description: ! ------------ ! Determine the array minimum/maximum and the corresponding indices. !-------------------------------------------------------------------------------! USE indices USE pegrid IMPLICIT NONE CHARACTER (LEN=*) :: mode INTEGER :: i, i1, i2, id_fmax, id_fmin, j, j1, j2, k, k1, k2, & fmax_ijk(3), fmax_ijk_l(3), fmin_ijk(3), & fmin_ijk_l(3), value_ijk(3) INTEGER, OPTIONAL :: value1_ijk(3) REAL :: value, & ar(i1:i2,j1:j2,k1:k2) #if defined( __ibm ) REAL (KIND=4) :: fmax(2), fmax_l(2), fmin(2), fmin_l(2) ! on 32bit- ! machines MPI_2REAL must not be replaced by ! MPI_2DOUBLE_PRECISION #else REAL :: fmax(2), fmax_l(2), fmin(2), fmin_l(2) #endif REAL, OPTIONAL :: value1 ! !-- Determine array minimum IF ( mode == 'min' .OR. mode == 'minmax' ) THEN ! !-- Determine the local minimum fmin_ijk_l = MINLOC( ar ) fmin_ijk_l(1) = i1 + fmin_ijk_l(1) - 1 ! MINLOC assumes lowerbound = 1 fmin_ijk_l(2) = j1 + fmin_ijk_l(2) - 1 fmin_ijk_l(3) = k1 + fmin_ijk_l(3) - 1 fmin_l(1) = ar(fmin_ijk_l(1),fmin_ijk_l(2),fmin_ijk_l(3)) #if defined( __parallel ) fmin_l(2) = myid CALL MPI_ALLREDUCE( fmin_l, fmin, 1, MPI_2REAL, MPI_MINLOC, comm2d, ierr ) ! !-- Determine the global minimum. Result stored on PE0. id_fmin = fmin(2) IF ( id_fmin /= 0 ) THEN IF ( myid == 0 ) THEN CALL MPI_RECV( fmin_ijk, 3, MPI_INTEGER, id_fmin, 0, comm2d, & status, ierr ) ELSEIF ( myid == id_fmin ) THEN CALL MPI_SEND( fmin_ijk_l, 3, MPI_INTEGER, 0, 0, comm2d, ierr ) ENDIF ELSE fmin_ijk = fmin_ijk_l ENDIF ! !-- Send the indices of the just determined array minimum to other PEs CALL MPI_BCAST( fmin_ijk, 3, MPI_INTEGER, 0, comm2d, ierr ) #else fmin(1) = fmin_l(1) fmin_ijk = fmin_ijk_l #endif ENDIF ! !-- Determine array maximum IF ( mode == 'max' .OR. mode == 'minmax' ) THEN ! !-- Determine the local maximum fmax_ijk_l = MAXLOC( ar ) fmax_ijk_l(1) = i1 + fmax_ijk_l(1) - 1 ! MAXLOC assumes lowerbound = 1 fmax_ijk_l(2) = j1 + fmax_ijk_l(2) - 1 fmax_ijk_l(3) = k1 + fmax_ijk_l(3) - 1 fmax_l(1) = ar(fmax_ijk_l(1),fmax_ijk_l(2),fmax_ijk_l(3)) #if defined( __parallel ) fmax_l(2) = myid CALL MPI_ALLREDUCE( fmax_l, fmax, 1, MPI_2REAL, MPI_MAXLOC, comm2d, ierr ) ! !-- Determine the global maximum. Result stored on PE0. id_fmax = fmax(2) IF ( id_fmax /= 0 ) THEN IF ( myid == 0 ) THEN CALL MPI_RECV( fmax_ijk, 3, MPI_INTEGER, id_fmax, 0, comm2d, & status, ierr ) ELSEIF ( myid == id_fmax ) THEN CALL MPI_SEND( fmax_ijk_l, 3, MPI_INTEGER, 0, 0, comm2d, ierr ) ENDIF ELSE fmax_ijk = fmax_ijk_l ENDIF ! !-- send the indices of the just determined array maximum to other PEs CALL MPI_BCAST( fmax_ijk, 3, MPI_INTEGER, 0, comm2d, ierr ) #else fmax(1) = fmax_l(1) fmax_ijk = fmax_ijk_l #endif ENDIF ! !-- Determine absolute array maximum IF ( mode == 'abs' ) THEN ! !-- Determine the local absolut maximum fmax_l(1) = 0.0 fmax_ijk_l(1) = i1 fmax_ijk_l(2) = j1 fmax_ijk_l(3) = k1 DO k = k1, k2 DO j = j1, j2 DO i = i1, i2 IF ( ABS( ar(i,j,k) ) > fmax_l(1) ) THEN fmax_l(1) = ABS( ar(i,j,k) ) fmax_ijk_l(1) = i fmax_ijk_l(2) = j fmax_ijk_l(3) = k ENDIF ENDDO ENDDO ENDDO ! !-- Set a flag in case that the determined value is negative. !-- A constant offset has to be subtracted in order to handle the special !-- case i=0 correctly IF ( ar(fmax_ijk_l(1),fmax_ijk_l(2),fmax_ijk_l(3)) < 0.0 ) THEN fmax_ijk_l(1) = -fmax_ijk_l(1) - 10 ENDIF #if defined( __parallel ) fmax_l(2) = myid CALL MPI_ALLREDUCE( fmax_l, fmax, 1, MPI_2REAL, MPI_MAXLOC, comm2d, & ierr ) ! !-- Determine the global absolut maximum. Result stored on PE0. id_fmax = fmax(2) IF ( id_fmax /= 0 ) THEN IF ( myid == 0 ) THEN CALL MPI_RECV( fmax_ijk, 3, MPI_INTEGER, id_fmax, 0, comm2d, & status, ierr ) ELSEIF ( myid == id_fmax ) THEN CALL MPI_SEND( fmax_ijk_l, 3, MPI_INTEGER, 0, 0, comm2d, ierr ) ENDIF ELSE fmax_ijk = fmax_ijk_l ENDIF ! !-- Send the indices of the just determined absolut maximum to other PEs CALL MPI_BCAST( fmax_ijk, 3, MPI_INTEGER, 0, comm2d, ierr ) #else fmax(1) = fmax_l(1) fmax_ijk = fmax_ijk_l #endif ENDIF ! !-- Determine output parameters SELECT CASE( mode ) CASE( 'min' ) value = fmin(1) value_ijk = fmin_ijk CASE( 'max' ) value = fmax(1) value_ijk = fmax_ijk CASE( 'minmax' ) value = fmin(1) value_ijk = fmin_ijk value1 = fmax(1) value1_ijk = fmax_ijk CASE( 'abs' ) value = fmax(1) value_ijk = fmax_ijk IF ( fmax_ijk(1) < 0 ) THEN value = -value value_ijk(1) = -value_ijk(1) - 10 ENDIF END SELECT ! !-- Limit index values to the range 0..nx, 0..ny IF ( value_ijk(3) == -1 ) value_ijk(3) = nx IF ( value_ijk(3) == nx+1 ) value_ijk(3) = 0 IF ( value_ijk(2) == -1 ) value_ijk(2) = ny IF ( value_ijk(2) == ny+1 ) value_ijk(2) = 0 END SUBROUTINE global_min_max