source: trunk/WRF.COMMON/WRFV2/frame/module_tiles.F @ 3094

Last change on this file since 3094 was 11, checked in by aslmd, 14 years ago

spiga@svn-planeto:ajoute le modele meso-echelle martien

File size: 16.1 KB
Line 
1!WRF:DRIVER_LAYER:TILING
2!
3
4MODULE 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
12CONTAINS
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
424END MODULE module_tiles
425
Note: See TracBrowser for help on using the repository browser.