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

Last change on this file since 1189 was 1189, checked in by heinze, 8 years ago

last commit documented

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