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

Last change on this file since 1806 was 1683, checked in by knoop, 9 years ago

last commit documented

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