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