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

Last change on this file since 1320 was 1320, checked in by raasch, 7 years ago

ONLY-attribute added to USE-statements,
kind-parameters added to all INTEGER and REAL declaration statements,
kinds are defined in new module kinds,
old module precision_kind is removed,
revision history before 2012 removed,
comment fields (!:) to be used for variable explanations added to all variable declaration statements

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