source: palm/trunk/SOURCE/global_min_max.f90 @ 866

Last change on this file since 866 was 866, checked in by raasch, 12 years ago

bugfix for timestep calculation in case of Galilei transformation

  • Property svn:keywords set to Id
File size: 9.1 KB
RevLine 
[866]1 SUBROUTINE global_min_max( i1, i2, j1, j2, k1, k2, ar, mode, offset, value, &
[1]2                            value_ijk, value1, value1_ijk )
3
[623]4!------------------------------------------------------------------------------!
[484]5! Current revisions:
[866]6! ------------------
7! new mode "absoff" accounts for an offset in the respective array
[1]8!
9! Former revisions:
10! -----------------
[3]11! $Id: global_min_max.f90 866 2012-03-28 06:44:41Z raasch $
[623]12!
[668]13! 667 2010-12-23 12:06:00Z suehring/gryschka
14! Adapting of the index arrays, because MINLOC assumes lowerbound at 1 and not
15! at nbgp.
16!
[623]17! 622 2010-12-10 08:08:13Z raasch
18! optional barriers included in order to speed up collective operations
19!
20! Feb. 2007
[3]21! RCS Log replace by Id keyword, revision history cleaned up
22!
[1]23! Revision 1.11  2003/04/16 12:56:58  raasch
24! Index values of the extrema are limited to the range 0..nx, 0..ny
25!
26! Revision 1.1  1997/07/24 11:14:03  raasch
27! Initial revision
28!
29!
30! Description:
31! ------------
32! Determine the array minimum/maximum and the corresponding indices.
[623]33!------------------------------------------------------------------------------!
[1]34
35    USE indices
36    USE pegrid
37
38    IMPLICIT NONE
39
40    CHARACTER (LEN=*) ::  mode
41
42    INTEGER           ::  i, i1, i2, id_fmax, id_fmin, j, j1, j2, k, k1, k2, &
43                          fmax_ijk(3), fmax_ijk_l(3), fmin_ijk(3), &
44                          fmin_ijk_l(3), value_ijk(3)
45    INTEGER, OPTIONAL ::  value1_ijk(3)
[866]46    REAL              ::  offset, value, &
[1]47                          ar(i1:i2,j1:j2,k1:k2)
48#if defined( __ibm )
49    REAL (KIND=4)     ::  fmax(2), fmax_l(2), fmin(2), fmin_l(2)  ! on 32bit-
50                          ! machines MPI_2REAL must not be replaced by
51                          ! MPI_2DOUBLE_PRECISION
52#else
53    REAL              ::  fmax(2), fmax_l(2), fmin(2), fmin_l(2)
54#endif
55    REAL, OPTIONAL    ::  value1
56
57
58!
59!-- Determine array minimum
60    IF ( mode == 'min'  .OR.  mode == 'minmax' )  THEN
61
62!
63!--    Determine the local minimum
64       fmin_ijk_l = MINLOC( ar )
[866]65       fmin_ijk_l(1) = i1 + fmin_ijk_l(1) - nbgp ! MINLOC assumes lowerbound = 1
[667]66       fmin_ijk_l(2) = j1 + fmin_ijk_l(2) - nbgp
[1]67       fmin_ijk_l(3) = k1 + fmin_ijk_l(3) - 1
68       fmin_l(1)  = ar(fmin_ijk_l(1),fmin_ijk_l(2),fmin_ijk_l(3))
69
70#if defined( __parallel )
71       fmin_l(2)  = myid
[622]72       IF ( collective_wait )  CALL MPI_BARRIER( comm2d, ierr )
[623]73       CALL MPI_ALLREDUCE( fmin_l, fmin, 1, MPI_2REAL, MPI_MINLOC, comm2d, &
74                           ierr )
[1]75
76!
77!--    Determine the global minimum. Result stored on PE0.
78       id_fmin = fmin(2)
79       IF ( id_fmin /= 0 )  THEN
80          IF ( myid == 0 )  THEN
81             CALL MPI_RECV( fmin_ijk, 3, MPI_INTEGER, id_fmin, 0, comm2d, &
82                            status, ierr )
83          ELSEIF ( myid == id_fmin )  THEN
84             CALL MPI_SEND( fmin_ijk_l, 3, MPI_INTEGER, 0, 0, comm2d, ierr )
85          ENDIF
86       ELSE
87          fmin_ijk = fmin_ijk_l
88       ENDIF
89!
90!--    Send the indices of the just determined array minimum to other PEs
91       CALL MPI_BCAST( fmin_ijk, 3, MPI_INTEGER, 0, comm2d, ierr )
92#else
93       fmin(1)  = fmin_l(1)
94       fmin_ijk = fmin_ijk_l
95#endif
96
97    ENDIF
98
99!
100!-- Determine array maximum
101    IF ( mode == 'max'  .OR.  mode == 'minmax' )  THEN
102
103!
104!--    Determine the local maximum
105       fmax_ijk_l = MAXLOC( ar )
[866]106       fmax_ijk_l(1) = i1 + fmax_ijk_l(1) - nbgp ! MAXLOC assumes lowerbound = 1
[667]107       fmax_ijk_l(2) = j1 + fmax_ijk_l(2) - nbgp
[1]108       fmax_ijk_l(3) = k1 + fmax_ijk_l(3) - 1
109       fmax_l(1) = ar(fmax_ijk_l(1),fmax_ijk_l(2),fmax_ijk_l(3))
110
111#if defined( __parallel )
112       fmax_l(2)  = myid
[622]113       IF ( collective_wait )  CALL MPI_BARRIER( comm2d, ierr )
[623]114       CALL MPI_ALLREDUCE( fmax_l, fmax, 1, MPI_2REAL, MPI_MAXLOC, comm2d, &
115                           ierr )
[1]116
117!
118!--    Determine the global maximum. Result stored on PE0.
119       id_fmax = fmax(2)
120       IF ( id_fmax /= 0 )  THEN
121          IF ( myid == 0 )  THEN
122             CALL MPI_RECV( fmax_ijk, 3, MPI_INTEGER, id_fmax, 0, comm2d, &
123                            status, ierr )
124          ELSEIF ( myid == id_fmax )  THEN
125             CALL MPI_SEND( fmax_ijk_l, 3, MPI_INTEGER, 0, 0, comm2d, ierr )
126          ENDIF
127       ELSE
128          fmax_ijk = fmax_ijk_l
129       ENDIF
130!
131!--    send the indices of the just determined array maximum to other PEs
132       CALL MPI_BCAST( fmax_ijk, 3, MPI_INTEGER, 0, comm2d, ierr )
133#else
134       fmax(1)  = fmax_l(1)
135       fmax_ijk = fmax_ijk_l
136#endif
137
138    ENDIF
139
140!
141!-- Determine absolute array maximum
142    IF ( mode == 'abs' )  THEN
143
144!
145!--    Determine the local absolut maximum
146       fmax_l(1)     = 0.0
147       fmax_ijk_l(1) =  i1
148       fmax_ijk_l(2) =  j1
149       fmax_ijk_l(3) =  k1
150       DO  k = k1, k2
151          DO  j = j1, j2
152             DO  i = i1, i2
153                IF ( ABS( ar(i,j,k) ) > fmax_l(1) )  THEN
154                   fmax_l(1) = ABS( ar(i,j,k) )
155                   fmax_ijk_l(1) = i
156                   fmax_ijk_l(2) = j
157                   fmax_ijk_l(3) = k
158                ENDIF
159             ENDDO
160          ENDDO
161       ENDDO
162
163!
164!--    Set a flag in case that the determined value is negative.
165!--    A constant offset has to be subtracted in order to handle the special
166!--    case i=0 correctly
167       IF ( ar(fmax_ijk_l(1),fmax_ijk_l(2),fmax_ijk_l(3)) < 0.0 )  THEN
168          fmax_ijk_l(1) = -fmax_ijk_l(1) - 10
169       ENDIF
170
171#if defined( __parallel )
172       fmax_l(2)  = myid
[622]173       IF ( collective_wait )  CALL MPI_BARRIER( comm2d, ierr )
[1]174       CALL MPI_ALLREDUCE( fmax_l, fmax, 1, MPI_2REAL, MPI_MAXLOC, comm2d, &
175                           ierr )
176
177!
178!--    Determine the global absolut maximum. Result stored on PE0.
179       id_fmax = fmax(2)
180       IF ( id_fmax /= 0 )  THEN
181          IF ( myid == 0 )  THEN
182             CALL MPI_RECV( fmax_ijk, 3, MPI_INTEGER, id_fmax, 0, comm2d, &
183                            status, ierr )
184          ELSEIF ( myid == id_fmax )  THEN
185             CALL MPI_SEND( fmax_ijk_l, 3, MPI_INTEGER, 0, 0, comm2d, ierr )
186          ENDIF
187       ELSE
188          fmax_ijk = fmax_ijk_l
189       ENDIF
190!
191!--    Send the indices of the just determined absolut maximum to other PEs
192       CALL MPI_BCAST( fmax_ijk, 3, MPI_INTEGER, 0, comm2d, ierr )
193#else
194       fmax(1)  = fmax_l(1)
195       fmax_ijk = fmax_ijk_l
196#endif
197
198    ENDIF
199
200!
[866]201!-- Determine absolute maximum of ( array - offset )
202    IF ( mode == 'absoff' )  THEN
203
204!
205!--    Determine the local absolut maximum
206       fmax_l(1)     = 0.0
207       fmax_ijk_l(1) =  i1
208       fmax_ijk_l(2) =  j1
209       fmax_ijk_l(3) =  k1
210       DO  k = k1, k2
211          DO  j = j1, j2
212!
213!--          Attention: the lowest gridpoint is excluded here, because there
214!--          ---------  is no advection at nzb=0 and mode 'absoff' is only
215!--                     used for calculating u,v extrema for CFL-criteria
216             DO  i = i1+1, i2
217                IF ( ABS( ar(i,j,k) - offset ) > fmax_l(1) )  THEN
218                   fmax_l(1) = ABS( ar(i,j,k) - offset )
219                   fmax_ijk_l(1) = i
220                   fmax_ijk_l(2) = j
221                   fmax_ijk_l(3) = k
222                ENDIF
223             ENDDO
224          ENDDO
225       ENDDO
226
227!
228!--    Set a flag in case that the determined value is negative.
229!--    A constant offset has to be subtracted in order to handle the special
230!--    case i=0 correctly
231       IF ( ar(fmax_ijk_l(1),fmax_ijk_l(2),fmax_ijk_l(3)) < 0.0 )  THEN
232          fmax_ijk_l(1) = -fmax_ijk_l(1) - 10
233       ENDIF
234
235#if defined( __parallel )
236       fmax_l(2)  = myid
237       IF ( collective_wait )  CALL MPI_BARRIER( comm2d, ierr )
238       CALL MPI_ALLREDUCE( fmax_l, fmax, 1, MPI_2REAL, MPI_MAXLOC, comm2d, &
239                           ierr )
240
241!
242!--    Determine the global absolut maximum. Result stored on PE0.
243       id_fmax = fmax(2)
244       IF ( id_fmax /= 0 )  THEN
245          IF ( myid == 0 )  THEN
246             CALL MPI_RECV( fmax_ijk, 3, MPI_INTEGER, id_fmax, 0, comm2d, &
247                            status, ierr )
248          ELSEIF ( myid == id_fmax )  THEN
249             CALL MPI_SEND( fmax_ijk_l, 3, MPI_INTEGER, 0, 0, comm2d, ierr )
250          ENDIF
251       ELSE
252          fmax_ijk = fmax_ijk_l
253       ENDIF
254!
255!--    Send the indices of the just determined absolut maximum to other PEs
256       CALL MPI_BCAST( fmax_ijk, 3, MPI_INTEGER, 0, comm2d, ierr )
257#else
258       fmax(1)  = fmax_l(1)
259       fmax_ijk = fmax_ijk_l
260#endif
261
262    ENDIF
263
264!
[1]265!-- Determine output parameters
266    SELECT CASE( mode )
267
268       CASE( 'min' )
269
270          value     = fmin(1)
271          value_ijk = fmin_ijk
272
273       CASE( 'max' )
274
275          value     = fmax(1)
276          value_ijk = fmax_ijk
277
278       CASE( 'minmax' )
279
280          value      = fmin(1)
281          value_ijk  = fmin_ijk
282          value1     = fmax(1)
283          value1_ijk = fmax_ijk
284
[866]285       CASE( 'abs', 'absoff' )
[1]286
287          value     = fmax(1)
288          value_ijk = fmax_ijk
289          IF ( fmax_ijk(1) < 0 )  THEN
290             value        = -value
[667]291             value_ijk(1) = -value_ijk(1) - 10         !???
[1]292          ENDIF
293
294    END SELECT
295
296!
297!-- Limit index values to the range 0..nx, 0..ny
[667]298    IF ( value_ijk(3) < 0  ) value_ijk(3) = nx +1 + value_ijk(3)
299    IF ( value_ijk(3) > nx ) value_ijk(3) = value_ijk(3) - (nx+1)
300    IF ( value_ijk(2) < 0  ) value_ijk(2) = ny +1 + value_ijk(2)
301    IF ( value_ijk(2) > ny ) value_ijk(2) = value_ijk(2) - (ny+1)
[1]302
303
304 END SUBROUTINE global_min_max
Note: See TracBrowser for help on using the repository browser.