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

Last change on this file since 2699 was 2696, checked in by kanani, 7 years ago

Merge of branch palm4u into trunk

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