module_tiles.F

References to this file elsewhere.
1 !WRF:DRIVER_LAYER:TILING
2 !
3 
4 MODULE module_tiles
5 
6   USE module_configure
7 
8   INTERFACE set_tiles
9     MODULE PROCEDURE set_tiles1 , set_tiles2, set_tiles3
10   END INTERFACE
11 
12 CONTAINS
13 
14 ! CPP macro for error checking
15 #define ERROR_TEST(A,O,B) IF( A O B )THEN;WRITE(mess,'(3A4)')'A','O','B';CALL WRF_ERROR_FATAL(mess);ENDIF
16 
17 ! this version is used to compute only on a boundary of some width
18 ! The ids, ide, jds, and jde arguments specify the edge of the boundary (a way of
19 ! accounting for staggering, and the bdyw gives the number of cells
20 ! (idea: if bdyw is negative, have it do the reverse and specify the 
21 ! interior, less the boundary.
22 
23   SUBROUTINE set_tiles1 ( grid , ids , ide , jds , jde , bdyw )
24 
25      USE module_domain
26      USE module_driver_constants
27      USE module_machine
28      USE module_wrf_error
29 
30      IMPLICIT NONE
31   
32      !  Input data.
33   
34      TYPE(domain)                   , INTENT(INOUT)  :: grid
35      INTEGER                        , INTENT(IN)     :: ids , ide , jds , jde , bdyw
36 
37      !  Local data
38 
39      INTEGER                                :: spx, epx, spy, epy, t, tt, ts, te
40      INTEGER                                :: smx, emx, smy, emy
41      INTEGER                                :: ntiles , num_tiles
42 
43      CHARACTER*80              :: mess
44 
45      data_ordering : SELECT CASE ( model_data_order )
46        CASE  ( DATA_ORDER_XYZ )
47          spx = grid%sp31 ; epx = grid%ep31 ; spy = grid%sp32 ; epy = grid%ep32
48        CASE  ( DATA_ORDER_YXZ )
49          spx = grid%sp32 ; epx = grid%ep32 ; spy = grid%sp31 ; epy = grid%ep31
50        CASE  ( DATA_ORDER_ZXY )
51          spx = grid%sp32 ; epx = grid%ep32 ; spy = grid%sp33 ; epy = grid%ep33
52        CASE  ( DATA_ORDER_ZYX )
53          spx = grid%sp33 ; epx = grid%ep33 ; spy = grid%sp32 ; epy = grid%ep32
54        CASE  ( DATA_ORDER_XZY )
55          spx = grid%sp31 ; epx = grid%ep31 ; spy = grid%sp33 ; epy = grid%ep33
56        CASE  ( DATA_ORDER_YZX )
57          spx = grid%sp33 ; epx = grid%ep33 ; spy = grid%sp31 ; epy = grid%ep31
58      END SELECT data_ordering
59 
60      num_tiles = 4
61 
62      IF ( num_tiles > grid%max_tiles ) THEN
63        IF ( ASSOCIATED(grid%i_start) ) THEN ; DEALLOCATE( grid%i_start ) ; NULLIFY( grid%i_start ) ; ENDIF
64        IF ( ASSOCIATED(grid%i_end) )   THEN ; DEALLOCATE( grid%i_end   ) ; NULLIFY( grid%i_end   ) ; ENDIF
65        IF ( ASSOCIATED(grid%j_start) ) THEN ; DEALLOCATE( grid%j_start ) ; NULLIFY( grid%j_start ) ; ENDIF
66        IF ( ASSOCIATED(grid%j_end) )   THEN ; DEALLOCATE( grid%j_end   ) ; NULLIFY( grid%j_end   ) ; ENDIF
67        ALLOCATE(grid%i_start(num_tiles))
68        ALLOCATE(grid%i_end(num_tiles))
69        ALLOCATE(grid%j_start(num_tiles))
70        ALLOCATE(grid%j_end(num_tiles))
71        grid%max_tiles = num_tiles
72      ENDIF
73 
74 ! XS boundary
75      IF      ( ids .ge. spx .and. ids .le. epx ) THEN
76         grid%i_start(1) = ids
77         grid%i_end(1)   = min( ids+bdyw-1 , epx )
78         grid%j_start(1) = max( spy , jds )
79         grid%j_end(1)   = min( epy , jde )
80      ELSEIF  ( (ids+bdyw-1) .ge. spx .and. (ids+bdyw-1) .le. epx ) THEN
81         grid%i_start(1) = max( ids , spx )
82         grid%i_end(1)   = ids+bdyw-1
83         grid%j_start(1) = max( spy , jds )
84         grid%j_end(1)   = min( epy , jde )
85      ELSE
86         grid%i_start(1) = 1
87         grid%i_end(1)   = -1
88         grid%j_start(1) = 1
89         grid%j_end(1)   = -1
90      ENDIF
91 
92 ! XE boundary
93      IF      ( ide .ge. spx .and. ide .le. epx ) THEN
94         grid%i_start(2) = max( ide-bdyw+1 , spx )
95         grid%i_end(2)   = ide
96         grid%j_start(2) = max( spy , jds )
97         grid%j_end(2)   = min( epy , jde )
98      ELSEIF  ( (ide-bdyw+1) .ge. spx .and. (ide-bdyw+1) .le. epx ) THEN
99         grid%i_start(2) = ide-bdyw+1
100         grid%i_end(2)   = min( ide , epx )
101         grid%j_start(2) = max( spy , jds )
102         grid%j_end(2)   = min( epy , jde )
103      ELSE
104         grid%i_start(2) = 1
105         grid%i_end(2)   = -1
106         grid%j_start(2) = 1
107         grid%j_end(2)   = -1
108      ENDIF
109 
110 ! YS boundary (note that the corners may already be done by XS and XE)
111      IF      ( jds .ge. spy .and. jds .le. epy ) THEN
112         grid%j_start(3) = jds
113         grid%j_end(3)   = min( jds+bdyw-1 , epy )
114         grid%i_start(3) = max( spx , ids+bdyw )
115         grid%i_end(3)   = min( epx , ide-bdyw )
116      ELSEIF  ( (jds+bdyw-1) .ge. spy .and. (jds+bdyw-1) .le. epy ) THEN
117         grid%j_start(3) = max( jds , spy )
118         grid%j_end(3)   = jds+bdyw-1
119         grid%i_start(3) = max( spx , ids+bdyw )
120         grid%i_end(3)   = min( epx , ide-bdyw )
121      ELSE
122         grid%j_start(3) = 1
123         grid%j_end(3)   = -1
124         grid%i_start(3) = 1
125         grid%i_end(3)   = -1
126      ENDIF
127 
128 ! YE boundary (note that the corners may already be done by XS and XE)
129      IF      ( jde .ge. spy .and. jde .le. epy ) THEN
130         grid%j_start(4) = max( jde-bdyw+1 , spy )
131         grid%j_end(4)   = jde
132         grid%i_start(4) = max( spx , ids+bdyw )
133         grid%i_end(4)   = min( epx , ide-bdyw )
134      ELSEIF  ( (jde-bdyw+1) .ge. spy .and. (jde-bdyw+1) .le. epy ) THEN
135         grid%j_start(4) = jde-bdyw+1
136         grid%j_end(4)   = min( jde , epy )
137         grid%i_start(4) = max( spx , ids+bdyw )
138         grid%i_end(4)   = min( epx , ide-bdyw )
139      ELSE
140         grid%j_start(4) = 1
141         grid%j_end(4)   = -1
142         grid%i_start(4) = 1
143         grid%i_end(4)   = -1
144      ENDIF
145 
146      grid%num_tiles = num_tiles
147 
148      RETURN
149   END SUBROUTINE set_tiles1
150 
151 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
152 ! this version is used to limit the domain or compute onto halos
153   SUBROUTINE set_tiles2 ( grid , ids , ide , jds , jde , ips , ipe , jps , jpe )
154      USE module_domain
155      USE module_driver_constants
156      USE module_machine
157      USE module_wrf_error
158 
159      IMPLICIT NONE
160   
161      !  Input data.
162   
163      TYPE(domain)                   , INTENT(INOUT)  :: grid
164      INTEGER                        , INTENT(IN)     :: ids , ide , jds , jde
165      INTEGER                        , INTENT(IN)     :: ips , ipe , jps , jpe
166 
167      !  Output data.
168 
169      !  Local data.
170   
171      INTEGER                                :: num_tiles_x, num_tiles_y, num_tiles
172      INTEGER                                :: tile_sz_x, tile_sz_y
173      INTEGER                                :: spx, epx, spy, epy, t, tt, ts, te
174      INTEGER                                :: smx, emx, smy, emy
175      INTEGER                                :: ntiles
176      INTEGER                                :: one
177 #ifdef _OPENMP
178      INTEGER , EXTERNAL        :: omp_get_max_threads
179 #endif
180      CHARACTER*80              :: mess
181 
182      data_ordering : SELECT CASE ( model_data_order )
183        CASE  ( DATA_ORDER_XYZ )
184          spx = grid%sp31 ; epx = grid%ep31 ; spy = grid%sp32 ; epy = grid%ep32
185          smx = grid%sm31 ; emx = grid%em31 ; smy = grid%sm32 ; emy = grid%em32
186        CASE  ( DATA_ORDER_YXZ )
187          spx = grid%sp32 ; epx = grid%ep32 ; spy = grid%sp31 ; epy = grid%ep31
188          smx = grid%sm32 ; emx = grid%em32 ; smy = grid%sm31 ; emy = grid%em31
189        CASE  ( DATA_ORDER_ZXY )
190          spx = grid%sp32 ; epx = grid%ep32 ; spy = grid%sp33 ; epy = grid%ep33
191          smx = grid%sm32 ; emx = grid%em32 ; smy = grid%sm33 ; emy = grid%em33
192        CASE  ( DATA_ORDER_ZYX )
193          spx = grid%sp33 ; epx = grid%ep33 ; spy = grid%sp32 ; epy = grid%ep32
194          smx = grid%sm33 ; emx = grid%em33 ; smy = grid%sm32 ; emy = grid%em32
195        CASE  ( DATA_ORDER_XZY )
196          spx = grid%sp31 ; epx = grid%ep31 ; spy = grid%sp33 ; epy = grid%ep33
197          smx = grid%sm31 ; emx = grid%em31 ; smy = grid%sm33 ; emy = grid%em33
198        CASE  ( DATA_ORDER_YZX )
199          spx = grid%sp33 ; epx = grid%ep33 ; spy = grid%sp31 ; epy = grid%ep31
200          smx = grid%sm33 ; emx = grid%em33 ; smy = grid%sm31 ; emy = grid%em31
201      END SELECT data_ordering
202 
203      ERROR_TEST(ips,<,smx)
204      ERROR_TEST(ipe,>,emx)
205      ERROR_TEST(jps,<,smy)
206      ERROR_TEST(jpe,>,emy)
207 
208      ! Here's how the number of tiles is arrived at:
209      !
210      !          if tile sizes are specified use those otherwise
211      !          if num_tiles is specified use that otherwise
212      !          if omp provides a value use that otherwise
213      !          use 1.
214      !
215 
216      IF ( grid%num_tiles_spec .EQ. 0 ) THEN
217        CALL nl_get_numtiles( 1, num_tiles )
218        IF ( num_tiles .EQ. 1 ) THEN
219 #ifdef _OPENMP
220          num_tiles = omp_get_max_threads()
221          WRITE(mess,'("WRF NUMBER OF TILES FROM OMP_GET_MAX_THREADS = ",I3)')num_tiles
222          CALL WRF_MESSAGE ( mess )
223 #else
224          num_tiles = 1
225 #endif
226        ENDIF
227 ! override num_tiles setting (however gotten) if tile sizes are specified
228        CALL nl_get_tile_sz_x( 1, tile_sz_x )
229        CALL nl_get_tile_sz_y( 1, tile_sz_y )
230        IF ( tile_sz_x >= 1 .and. tile_sz_y >= 1 ) THEN
231         ! figure number of whole tiles and add 1 for any partials in each dim
232           num_tiles_x = (epx-spx+1) / tile_sz_x
233           if ( tile_sz_x*num_tiles_x < epx-spx+1 ) num_tiles_x = num_tiles_x + 1
234           num_tiles_y = (epy-spy+1) / tile_sz_y
235           if ( tile_sz_y*num_tiles_y < epy-spy+1 ) num_tiles_y = num_tiles_y + 1
236           num_tiles = num_tiles_x * num_tiles_y
237        ELSE
238          IF      ( machine_info%tile_strategy == TILE_X ) THEN
239            num_tiles_x = num_tiles
240            num_tiles_y = 1
241          ELSE IF ( machine_info%tile_strategy == TILE_Y ) THEN
242            num_tiles_x = 1
243            num_tiles_y = num_tiles
244          ELSE ! Default ( machine_info%tile_strategy == TILE_XY ) THEN
245            one = 1
246            call least_aspect( num_tiles, one, one, num_tiles_y, num_tiles_x )
247          ENDIF
248        ENDIF
249        grid%num_tiles_spec = num_tiles
250        grid%num_tiles_x = num_tiles_x
251        grid%num_tiles_y = num_tiles_y
252        WRITE(mess,'("WRF NUMBER OF TILES = ",I3)')num_tiles
253        CALL WRF_MESSAGE ( mess )
254      ENDIF
255 
256      num_tiles   = grid%num_tiles_spec
257      num_tiles_x = grid%num_tiles_x
258      num_tiles_y = grid%num_tiles_y
259 
260      IF ( num_tiles > grid%max_tiles ) THEN
261        IF ( ASSOCIATED(grid%i_start) ) THEN ; DEALLOCATE( grid%i_start ) ; NULLIFY( grid%i_start ) ; ENDIF
262        IF ( ASSOCIATED(grid%i_end) )   THEN ; DEALLOCATE( grid%i_end   ) ; NULLIFY( grid%i_end   ) ; ENDIF
263        IF ( ASSOCIATED(grid%j_start) ) THEN ; DEALLOCATE( grid%j_start ) ; NULLIFY( grid%j_start ) ; ENDIF
264        IF ( ASSOCIATED(grid%j_end) )   THEN ; DEALLOCATE( grid%j_end   ) ; NULLIFY( grid%j_end   ) ; ENDIF
265        ALLOCATE(grid%i_start(num_tiles))
266        ALLOCATE(grid%i_end(num_tiles))
267        ALLOCATE(grid%j_start(num_tiles))
268        ALLOCATE(grid%j_end(num_tiles))
269        grid%max_tiles = num_tiles
270      ENDIF
271 
272      DO t = 0, num_tiles-1
273         ntiles = mod(t,num_tiles_x)
274         CALL region_bounds( spx, epx,                                  &
275                             num_tiles_x, ntiles,                       &
276                             ts, te )
277 !!!
278 ! This bit allows the user to specify execution out onto the halo region
279 ! in the call to set_tiles. If the low patch boundary specified by the arguments
280 ! is less than what the model already knows to be the patch boundary and if
281 ! the user hasn't erred by specifying something that would fall off memory
282 ! (safety tests are higher up in this routine, outside the IF) then adjust
283 ! the tile boundary of the low edge tiles accordingly. Likewise for high edges.
284         IF ( ips .lt. spx .and. ts .eq. spx ) ts = ips ;
285         IF ( ipe .gt. epx .and. te .eq. epx ) te = ipe ;
286 !!!
287         grid%i_start(t+1) = max ( ts , ids )
288         grid%i_end(t+1)   = min ( te , ide )
289         ntiles = t / num_tiles_x
290         CALL region_bounds( spy, epy,                                  &
291                             num_tiles_y, ntiles,                       &
292                             ts, te )
293 !
294         IF ( jps .lt. spy .and. ts .eq. spy ) ts = jps ;
295         IF ( jpe .gt. epy .and. te .eq. epy ) te = jpe ;
296 !
297         grid%j_start(t+1) = max ( ts , jds )
298         grid%j_end(t+1)   = min ( te , jde )
299      END DO
300      grid%num_tiles = num_tiles
301 
302      RETURN
303   END SUBROUTINE set_tiles2
304   
305 
306 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
307 ! this version sets the tiles based on a passed in integer mask
308 ! the assumption here is that the mask is relatively straigthforward
309 ! and coverable with 2 or three rectangles. No weird stuff...
310 
311   SUBROUTINE set_tiles3 ( grid , imask, ims, ime, jms, jme, ips, ipe, jps, jpe )
312      USE module_domain
313      USE module_driver_constants
314      USE module_machine
315      USE module_wrf_error
316 
317      IMPLICIT NONE
318   
319      !  Input data.
320   
321      TYPE(domain)                   , INTENT(INOUT)  :: grid
322      INTEGER                        , INTENT(IN)     :: ims , ime , jms , jme
323      INTEGER                        , INTENT(IN)     :: ips , ipe , jps , jpe
324      INTEGER, DIMENSION(ims:ime,jms:jme), INTENT(IN) :: imask
325      INTEGER                :: num_tiles
326      INTEGER, DIMENSION(50) :: i_start, i_end, j_start, j_end
327 
328      !  Output data.
329 
330      !  Local data.
331   
332      CHARACTER*80              :: mess
333 
334      CALL set_tiles_masked ( imask, ims, ime, jms, jme, ips, ipe, jps, jpe, &
335                              num_tiles, i_start, i_end, j_start, j_end )
336 
337      IF ( num_tiles > grid%max_tiles ) THEN
338        IF ( ASSOCIATED(grid%i_start) ) THEN ; DEALLOCATE( grid%i_start ) ; NULLIFY( grid%i_start ) ; ENDIF
339        IF ( ASSOCIATED(grid%i_end) )   THEN ; DEALLOCATE( grid%i_end   ) ; NULLIFY( grid%i_end   ) ; ENDIF
340        IF ( ASSOCIATED(grid%j_start) ) THEN ; DEALLOCATE( grid%j_start ) ; NULLIFY( grid%j_start ) ; ENDIF
341        IF ( ASSOCIATED(grid%j_end) )   THEN ; DEALLOCATE( grid%j_end   ) ; NULLIFY( grid%j_end   ) ; ENDIF
342        ALLOCATE(grid%i_start(num_tiles))
343        ALLOCATE(grid%i_end(num_tiles))
344        ALLOCATE(grid%j_start(num_tiles))
345        ALLOCATE(grid%j_end(num_tiles))
346        grid%max_tiles = num_tiles
347      ENDIF
348      grid%num_tiles = num_tiles
349      grid%i_start(1:num_tiles) = i_start(1:num_tiles)
350      grid%i_end(1:num_tiles)   = i_end(1:num_tiles)
351      grid%j_start(1:num_tiles) = j_start(1:num_tiles)
352      grid%j_end(1:num_tiles)   = j_end(1:num_tiles)
353 
354      RETURN
355   END SUBROUTINE set_tiles3
356 
357   SUBROUTINE set_tiles_masked ( imask, ims, ime, jms, jme, ips, ipe, jps, jpe, &
358                                 num_tiles, istarts, iends, jstarts, jends )
359 
360       IMPLICIT NONE
361 
362       !  Arguments
363 
364       INTEGER                        , INTENT(IN)     :: ims , ime , jms , jme
365       INTEGER, DIMENSION(ims:ime,jms:jme), INTENT(IN) :: imask
366       INTEGER                        , INTENT(IN)     :: ips , ipe , jps , jpe
367       INTEGER                        , INTENT(OUT)    :: num_tiles
368       INTEGER, DIMENSION(*)          , INTENT(OUT)    :: istarts, iends
369       INTEGER, DIMENSION(*)          , INTENT(OUT)    :: jstarts, jends
370 
371       !  Output data.
372 
373       !  Local data.
374       CHARACTER*80              :: mess
375       INTEGER :: i, j, ir, jr
376       INTEGER :: imaskcopy(ips:ipe,jps:jpe)    ! copy of imask to write on
377 
378       imaskcopy = imask(ips:ipe,jps:jpe)
379       num_tiles = 0
380       ! simple multi-pass scheme, optimize later...
381       DO WHILE (ANY(imaskcopy == 1))
382         DO j = jps,jpe
383           DO i = ips,ipe
384             ! find first "1" and build a rectangle from it
385             IF ( imaskcopy(i,j) == 1 ) THEN
386               num_tiles = num_tiles + 1
387               istarts(num_tiles) = i
388               iends(num_tiles)   = i
389               jstarts(num_tiles) = j
390               jends(num_tiles)   = j
391               ! don't check this point again
392               imaskcopy(i,j) = 0
393               ! find length of first row
394               DO ir = istarts(num_tiles)+1,ipe
395                 IF ( imaskcopy(ir,j) == 1 ) THEN
396                   iends(num_tiles) = ir
397                   ! don't check this point again
398                   imaskcopy(ir,j) = 0
399                 ELSE
400                   EXIT
401                 ENDIF
402               ENDDO
403               ! find number of rows
404               DO jr = jstarts(num_tiles)+1,jpe
405                 IF (ALL(imaskcopy(istarts(num_tiles):iends(num_tiles),jr) == 1)) THEN
406                   jends(num_tiles) = jr
407                   ! don't check these points again
408                   imaskcopy(istarts(num_tiles):iends(num_tiles),jr) = 0
409                 ELSE
410                   EXIT
411                 ENDIF
412               ENDDO
413             ENDIF   ! if ( imaskcopy(i,j) == 1 )
414           ENDDO
415         ENDDO
416       ENDDO
417       RETURN
418   END SUBROUTINE set_tiles_masked
419 
420   
421   SUBROUTINE init_module_tiles
422   END SUBROUTINE init_module_tiles
423 
424 END MODULE module_tiles
425