source: trunk/WRF.COMMON/WRFV3/share/mediation_nest_move.F @ 3567

Last change on this file since 3567 was 2759, checked in by aslmd, 2 years ago

adding unmodified code from WRFV3.0.1.1, expurged from useless data +1M size

File size: 34.1 KB
Line 
1
2SUBROUTINE med_nest_move ( parent, nest )
3  ! Driver layer
4   USE module_domain
5   USE module_timing
6   USE module_configure
7   USE module_io_domain
8   USE module_dm
9   TYPE(domain) , POINTER                     :: parent, nest, grid
10   INTEGER dx, dy       ! number of parent domain points to move
11#ifdef MOVE_NESTS
12  ! Local
13   CHARACTER*256 mess
14   INTEGER i, j, p, parent_grid_ratio
15   INTEGER px, py       ! number and direction of nd points to move
16   INTEGER                         :: ids , ide , jds , jde , kds , kde , &
17                                      ims , ime , jms , jme , kms , kme , &
18                                      ips , ipe , jps , jpe , kps , kpe
19   INTEGER ierr, fid
20   LOGICAL input_from_hires
21   LOGICAL saved_restart_value
22   TYPE (grid_config_rec_type)   :: config_flags
23   LOGICAL, EXTERNAL :: wrf_dm_on_monitor
24   LOGICAL, EXTERNAL :: should_not_move
25
26   INTERFACE
27     SUBROUTINE med_interp_domain ( parent , nest )
28        USE module_domain
29        TYPE(domain) , POINTER                 :: parent , nest
30     END SUBROUTINE med_interp_domain
31     SUBROUTINE start_domain ( grid , allowed_to_move )
32        USE module_domain
33        TYPE(domain) :: grid
34        LOGICAL, INTENT(IN) :: allowed_to_move
35     END SUBROUTINE start_domain
36#if ( EM_CORE == 1 )
37     SUBROUTINE shift_domain_em ( grid, disp_x, disp_y  &
38!
39# include <dummy_new_args.inc>
40!
41                           )
42        USE module_domain
43        USE module_configure
44        USE module_timing
45        IMPLICIT NONE
46        INTEGER disp_x, disp_y
47        TYPE(domain) , POINTER                 :: grid
48# include <dummy_new_decl.inc>
49     END SUBROUTINE shift_domain_em
50#endif
51#if ( NMM_CORE == 1 )
52#endif
53     LOGICAL FUNCTION time_for_move ( parent , nest , dx , dy )
54        USE module_domain
55        USE module_utility
56        TYPE(domain) , POINTER    :: parent , nest
57        INTEGER, INTENT(OUT)      :: dx , dy
58     END FUNCTION time_for_move
59     SUBROUTINE  input_terrain_rsmas ( grid ,                  &
60                           ids , ide , jds , jde , kds , kde , &
61                           ims , ime , jms , jme , kms , kme , &
62                           ips , ipe , jps , jpe , kps , kpe )
63       USE module_domain
64       TYPE ( domain ) :: grid
65       INTEGER                           :: ids , ide , jds , jde , kds , kde , &
66                                            ims , ime , jms , jme , kms , kme , &
67                                            ips , ipe , jps , jpe , kps , kpe
68     END SUBROUTINE input_terrain_rsmas
69     SUBROUTINE med_nest_feedback ( parent , nest , config_flags )
70       USE module_domain
71       USE module_configure
72       TYPE (domain), POINTER ::  nest , parent
73       TYPE (grid_config_rec_type) config_flags
74     END SUBROUTINE med_nest_feedback
75     SUBROUTINE  blend_terrain ( ter_interpolated , ter_input , &
76                           ids , ide , jds , jde , kds , kde , &
77                           ims , ime , jms , jme , kms , kme , &
78                           ips , ipe , jps , jpe , kps , kpe )
79       INTEGER                           :: ids , ide , jds , jde , kds , kde , &
80                                            ims , ime , jms , jme , kms , kme , &
81                                            ips , ipe , jps , jpe , kps , kpe
82       REAL , DIMENSION(ims:ime,jms:jme) :: ter_interpolated
83       REAL , DIMENSION(ims:ime,jms:jme) :: ter_input
84     END SUBROUTINE blend_terrain
85     SUBROUTINE  store_terrain ( ter_interpolated , ter_input , &
86                           ids , ide , jds , jde , kds , kde , &
87                           ims , ime , jms , jme , kms , kme , &
88                           ips , ipe , jps , jpe , kps , kpe )
89       INTEGER                           :: ids , ide , jds , jde , kds , kde , &
90                                            ims , ime , jms , jme , kms , kme , &
91                                            ips , ipe , jps , jpe , kps , kpe
92       REAL , DIMENSION(ims:ime,jms:jme) :: ter_interpolated
93       REAL , DIMENSION(ims:ime,jms:jme) :: ter_input
94     END SUBROUTINE store_terrain
95   END INTERFACE
96
97  ! set grid pointer for code in deref_kludge (if used)
98   grid => nest
99
100   IF ( should_not_move( nest%id ) ) THEN
101      CALL wrf_message( 'Nest movement is disabled because of namelist settings' )
102      RETURN
103   ENDIF
104
105! if the nest has stopped don't do all this
106   IF ( WRFU_ClockIsStopTime(nest%domain_clock ,rc=ierr) ) RETURN
107
108! mask should be defined in nest domain
109
110   check_for_move: IF ( time_for_move ( parent , nest , dx, dy ) ) THEN
111
112     IF ( (dx .gt. 1 .or. dx .lt. -1 ) .or.  &
113          (dy .gt. 1 .or. dy .lt. -1 ) ) THEN
114       WRITE(mess,*)' invalid move: dx, dy ', dx, dy
115       CALL wrf_error_fatal( mess )
116     ENDIF
117
118     WRITE(mess,*)' moving ',grid%id,dx,dy
119     CALL wrf_message(mess)
120
121     CALL get_ijk_from_grid (  grid ,                   &
122                               ids, ide, jds, jde, kds, kde,    &
123                               ims, ime, jms, jme, kms, kme,    &
124                               ips, ipe, jps, jpe, kps, kpe    )
125
126     CALL wrf_dm_move_nest ( parent, nest%intermediate_grid, dx, dy )
127
128     CALL adjust_domain_dims_for_move( nest%intermediate_grid , dx, dy )
129
130     CALL get_ijk_from_grid (  grid ,                   &
131                               ids, ide, jds, jde, kds, kde,    &
132                               ims, ime, jms, jme, kms, kme,    &
133                               ips, ipe, jps, jpe, kps, kpe    )
134
135     grid => nest
136
137#if ( EM_CORE == 1 )
138     CALL shift_domain_em( grid, dx, dy  &
139!
140# include <actual_new_args.inc>
141!
142                           )
143#endif
144
145     px = grid%parent_grid_ratio*dx
146     py = grid%parent_grid_ratio*dy
147
148     grid%i_parent_start = grid%i_parent_start + px / grid%parent_grid_ratio
149     CALL nl_set_i_parent_start( grid%id, grid%i_parent_start )
150     grid%j_parent_start = grid%j_parent_start + py / grid%parent_grid_ratio
151     CALL nl_set_j_parent_start( grid%id, grid%j_parent_start )
152
153     IF ( wrf_dm_on_monitor() ) THEN
154       write(mess,*)  &
155         'Grid ',grid%id,' New SW corner (in parent x and y):',grid%i_parent_start, grid%j_parent_start
156       CALL wrf_message(TRIM(mess))
157     ENDIF
158
159     CALL med_interp_domain( parent, nest )
160
161     CALL nl_get_input_from_hires( nest%id , input_from_hires )
162     IF ( input_from_hires ) THEN
163
164! store horizontally interpolated terrain in temp location
165       CALL  store_terrain ( nest%ht_fine , nest%ht , &
166                             ids , ide , jds , jde , 1   , 1   , &
167                             ims , ime , jms , jme , 1   , 1   , &
168                             ips , ipe , jps , jpe , 1   , 1   )
169       CALL  store_terrain ( nest%mub_fine , nest%mub , &
170                             ids , ide , jds , jde , 1   , 1   , &
171                             ims , ime , jms , jme , 1   , 1   , &
172                             ips , ipe , jps , jpe , 1   , 1   )
173       CALL  store_terrain ( nest%phb_fine , nest%phb , &
174                             ids , ide , jds , jde , kds , kde , &
175                             ims , ime , jms , jme , kms , kme , &
176                             ips , ipe , jps , jpe , kps , kpe )
177
178       CALL  input_terrain_rsmas ( nest,                               &
179                                   ids , ide , jds , jde , 1   , 1   , &
180                                   ims , ime , jms , jme , 1   , 1   , &
181                                   ips , ipe , jps , jpe , 1   , 1   )
182
183       CALL  blend_terrain ( nest%ht_fine , nest%ht , &
184                             ids , ide , jds , jde , 1   , 1   , &
185                             ims , ime , jms , jme , 1   , 1   , &
186                             ips , ipe , jps , jpe , 1   , 1   )
187       CALL  blend_terrain ( nest%mub_fine , nest%mub , &
188                             ids , ide , jds , jde , 1   , 1   , &
189                             ims , ime , jms , jme , 1   , 1   , &
190                             ips , ipe , jps , jpe , 1   , 1   )
191       CALL  blend_terrain ( nest%phb_fine , nest%phb , &
192                             ids , ide , jds , jde , kds , kde , &
193                             ims , ime , jms , jme , kms , kme , &
194                             ips , ipe , jps , jpe , kps , kpe )
195!
196       CALL model_to_grid_config_rec ( parent%id , model_config_rec , config_flags )
197
198       CALL med_nest_feedback ( parent , nest , config_flags )
199       parent%imask_nostag = 1
200       parent%imask_xstag = 1
201       parent%imask_ystag = 1
202       parent%imask_xystag = 1
203
204! start_domain will key off "restart". Even if this is a restart run
205! we don't want it to here. Save the value, set it to false, and restore afterwards
206       saved_restart_value = config_flags%restart
207       config_flags%restart = .FALSE.
208       grid%restart = .FALSE.
209       CALL nl_set_restart ( 1, .FALSE. )
210       grid%press_adj = .FALSE.
211       CALL start_domain ( parent , .FALSE. )
212       config_flags%restart = saved_restart_value
213       grid%restart = saved_restart_value
214       CALL nl_set_restart ( 1,  saved_restart_value )
215
216     ENDIF
217
218!
219! masks associated with nest will have been set by shift_domain_em above
220     nest%moved = .true.
221! start_domain will key off "restart". Even if this is a restart run
222! we don't want it to here. Save the value, set it to false, and restore afterwards
223     saved_restart_value = config_flags%restart
224     config_flags%restart = .FALSE.
225     CALL nl_set_restart ( 1, .FALSE. )
226     grid%restart = .FALSE.
227     nest%press_adj = .FALSE.
228     CALL start_domain ( nest , .FALSE. )
229     config_flags%restart = saved_restart_value
230     grid%restart = saved_restart_value
231     CALL nl_set_restart ( 1,  saved_restart_value )
232     nest%moved = .false.
233     
234!
235! copy time level 2 to time level 1 in new regions of multi-time level fields
236! this should be registry generated.
237!
238#if ( EM_CORE == 1 )
239      do k = kms,kme
240        where ( nest%imask_xstag  .EQ. 1 ) nest%u_1(:,k,:)   = nest%u_2(:,k,:)
241        where ( nest%imask_ystag  .EQ. 1 ) nest%v_1(:,k,:)   = nest%v_2(:,k,:)
242        where ( nest%imask_nostag .EQ. 1 ) nest%t_1(:,k,:)   = nest%t_2(:,k,:)
243        where ( nest%imask_nostag .EQ. 1 ) nest%w_1(:,k,:)   = nest%w_2(:,k,:)
244        where ( nest%imask_nostag .EQ. 1 ) nest%ph_1(:,k,:)  = nest%ph_2(:,k,:)
245        where ( nest%imask_nostag .EQ. 1 ) nest%tp_1(:,k,:)  = nest%tp_2(:,k,:)
246        where ( nest%imask_nostag .EQ. 1 ) nest%tke_1(:,k,:) = nest%tke_2(:,k,:)
247      enddo
248      where ( nest%imask_nostag .EQ. 1 ) nest%mu_1(:,:)  = nest%mu_2(:,:)
249#endif
250!
251   ENDIF check_for_move
252#endif
253END SUBROUTINE med_nest_move
254
255LOGICAL FUNCTION time_for_move2 ( parent , grid , move_cd_x, move_cd_y )
256  ! Driver layer
257   USE module_domain
258   USE module_configure
259   USE module_compute_geop
260   USE module_dm
261   USE module_utility
262   IMPLICIT NONE
263! Arguments
264   TYPE(domain) , POINTER    :: parent, grid
265   INTEGER, INTENT(OUT)      :: move_cd_x , move_cd_y
266#ifdef MOVE_NESTS
267! Local
268   INTEGER  num_moves, rc
269   INTEGER  move_interval , move_id
270   TYPE(WRFU_Time) :: ct, st
271   TYPE(WRFU_TimeInterval) :: ti
272   CHARACTER*256 mess, timestr
273   INTEGER     :: ids, ide, jds, jde, kds, kde, &
274                  ims, ime, jms, jme, kms, kme, &
275                  ips, ipe, jps, jpe, kps, kpe
276   INTEGER :: is, ie, js, je, ierr
277   REAL    :: ipbar, pbar, jpbar, fact
278   REAL    :: last_vc_i , last_vc_j
279
280   REAL, ALLOCATABLE, DIMENSION(:,:) :: height_l, height
281   REAL, ALLOCATABLE, DIMENSION(:,:) :: psfc, xlat, xlong, terrain
282   REAL :: minh, maxh
283   INTEGER :: mini, minj, maxi, maxj, i, j, pgr, irad
284   REAL :: disp_x, disp_y, lag, radius, center_i, center_j, dx
285   REAL :: dijsmooth, vmax, vmin, a, b
286   REAL :: dc_i, dc_j   ! domain center
287   REAL :: maxws, ws
288   REAL :: pmin
289   INTEGER imploc, jmploc
290
291   INTEGER :: fje, fjs, fie, fis, fimloc, fjmloc, imloc, jmloc
292   INTEGER :: i_parent_start, j_parent_start
293   INTEGER :: max_vortex_speed, vortex_interval  ! meters per second and seconds
294   INTEGER :: track_level
295   REAL    :: rsmooth = 100000.  ! in meters
296
297   LOGICAL, EXTERNAL :: wrf_dm_on_monitor
298
299character*256 message, message2
300
301!#define MOVING_DIAGS
302# ifdef VORTEX_CENTER
303
304
305   CALL nl_get_parent_grid_ratio ( grid%id , pgr )
306   CALL nl_get_i_parent_start    ( grid%id , i_parent_start )
307   CALL nl_get_j_parent_start    ( grid%id , j_parent_start )
308   CALL nl_get_track_level       ( grid%id , track_level )
309
310!  WRITE(mess,*)'Vortex is tracked at ', track_level
311!  CALL wrf_message(mess)
312
313   CALL get_ijk_from_grid (  grid ,                        &
314                             ids, ide, jds, jde, kds, kde, &
315                             ims, ime, jms, jme, kms, kme, &
316                             ips, ipe, jps, jpe, kps, kpe  )
317
318! If the alarm is ringing, recompute the Vortex Center (VC); otherwise
319! use the previous position of VC.  VC is not recomputed ever step to
320! save on cost for global collection of height field and broadcast
321! of new center.
322
323#  ifdef MOVING_DIAGS
324write(0,*)'Check to see if COMPUTE_VORTEX_CENTER_ALARM is ringing? '
325#  endif
326   CALL nl_get_parent_grid_ratio ( grid%id , pgr )
327   CALL nl_get_dx ( grid%id , dx )
328
329   IF ( WRFU_AlarmIsRinging( grid%alarms( COMPUTE_VORTEX_CENTER_ALARM ), rc=rc ) ) THEN
330
331#  ifdef MOVING_DIAGS
332     write(0,*)'COMPUTE_VORTEX_CENTER_ALARM is ringing  '
333#  endif
334     CALL WRFU_AlarmRingerOff( grid%alarms( COMPUTE_VORTEX_CENTER_ALARM ), rc=rc )
335     CALL domain_clock_get( grid, current_timestr=timestr )
336
337     last_vc_i = grid%vc_i
338     last_vc_j = grid%vc_j
339
340     ALLOCATE ( height_l ( ims:ime , jms:jme ), STAT=ierr )
341     IF ( ierr .ne. 0 ) CALL wrf_error_fatal ('allocating height_l in time_for_move2')
342     IF ( wrf_dm_on_monitor() ) THEN
343       ALLOCATE ( height   ( ids:ide , jds:jde ), STAT=ierr )
344       IF ( ierr .ne. 0 ) CALL wrf_error_fatal ('allocating height in time_for_move2')
345       ALLOCATE ( psfc     ( ids:ide , jds:jde ), STAT=ierr )
346       IF ( ierr .ne. 0 ) CALL wrf_error_fatal ('allocating psfc in time_for_move2')
347       ALLOCATE ( xlat     ( ids:ide , jds:jde ), STAT=ierr )
348       IF ( ierr .ne. 0 ) CALL wrf_error_fatal ('allocating xlat in time_for_move2')
349       ALLOCATE ( xlong    ( ids:ide , jds:jde ), STAT=ierr )
350       IF ( ierr .ne. 0 ) CALL wrf_error_fatal ('allocating xlong in time_for_move2')
351       ALLOCATE ( terrain  ( ids:ide , jds:jde ), STAT=ierr )
352       IF ( ierr .ne. 0 ) CALL wrf_error_fatal ('allocating terrain in time_for_move2')
353     ELSE
354       ALLOCATE ( height   ( 1:1 , 1:1 ), STAT=ierr )
355       IF ( ierr .ne. 0 ) CALL wrf_error_fatal ('allocating height in time_for_move2')
356       ALLOCATE ( psfc     ( 1:1 , 1:1 ), STAT=ierr )
357       IF ( ierr .ne. 0 ) CALL wrf_error_fatal ('allocating psfc in time_for_move2')
358       ALLOCATE ( xlat     ( 1:1 , 1:1 ), STAT=ierr )
359       IF ( ierr .ne. 0 ) CALL wrf_error_fatal ('allocating xlat in time_for_move2')
360       ALLOCATE ( xlong    ( 1:1 , 1:1 ), STAT=ierr )
361       IF ( ierr .ne. 0 ) CALL wrf_error_fatal ('allocating xlong in time_for_move2')
362       ALLOCATE ( terrain  ( 1:1 , 1:1 ), STAT=ierr )
363       IF ( ierr .ne. 0 ) CALL wrf_error_fatal ('allocating terrain in time_for_move2')
364     ENDIF
365
366#  if (EM_CORE == 1)
367     CALL compute_500mb_height ( grid%ph_2 , grid%phb, grid%p, grid%pb, height_l , &
368                                 track_level,                  &
369                                 ids, ide, jds, jde, kds, kde, &
370                                 ims, ime, jms, jme, kms, kme, &
371                                 ips, ipe, jps, jpe, kps, kpe  )
372#  endif
373
374     CALL wrf_patch_to_global_real ( height_l , height , grid%domdesc, "z", "xy", &
375                                     ids, ide-1 , jds , jde-1 , 1 , 1 , &
376                                     ims, ime   , jms , jme   , 1 , 1 , &
377                                     ips, ipe   , jps , jpe   , 1 , 1   )
378     CALL wrf_patch_to_global_real ( grid%psfc , psfc , grid%domdesc, "z", "xy", &
379                                     ids, ide-1 , jds , jde-1 , 1 , 1 , &
380                                     ims, ime   , jms , jme   , 1 , 1 , &
381                                     ips, ipe   , jps , jpe   , 1 , 1   )
382     CALL wrf_patch_to_global_real ( grid%xlat , xlat , grid%domdesc, "z", "xy", &
383                                     ids, ide-1 , jds , jde-1 , 1 , 1 , &
384                                     ims, ime   , jms , jme   , 1 , 1 , &
385                                     ips, ipe   , jps , jpe   , 1 , 1   )
386     CALL wrf_patch_to_global_real ( grid%xlong , xlong , grid%domdesc, "z", "xy", &
387                                     ids, ide-1 , jds , jde-1 , 1 , 1 , &
388                                     ims, ime   , jms , jme   , 1 , 1 , &
389                                     ips, ipe   , jps , jpe   , 1 , 1   )
390     CALL wrf_patch_to_global_real ( grid%ht , terrain , grid%domdesc, "z", "xy", &
391                                     ids, ide-1 , jds , jde-1 , 1 , 1 , &
392                                     ims, ime   , jms , jme   , 1 , 1 , &
393                                     ips, ipe   , jps , jpe   , 1 , 1   )
394
395! calculate max wind speed
396     maxws = 0.
397     do j = jps, jpe
398       do i = ips, ipe
399         ws = grid%u10(i,j) * grid%u10(i,j) + grid%v10(i,j) * grid%v10(i,j)
400         if ( ws > maxws ) maxws = ws
401       enddo
402     enddo
403     maxws = sqrt ( maxws )
404     maxws = wrf_dm_max_real ( maxws )
405
406     monitor_only : IF ( wrf_dm_on_monitor() ) THEN
407
408!
409! This vortex center finding code adapted from the Hurricane version of MM5,
410! Courtesy:
411!
412!   Shuyi Chen et al., Rosenstiel School of Marine and Atmos. Sci., U. Miami.
413!   Spring, 2005
414!
415! Get the first guess vortex center about which we do our search
416! as mini and minh; minimum value is minh
417!
418
419       CALL nl_get_vortex_interval( grid%id , vortex_interval )
420       CALL nl_get_max_vortex_speed( grid%id , max_vortex_speed )
421
422       IF ( grid%vc_i < 0. .AND. grid%vc_j < 0. ) THEN
423          ! first time through
424          is = ids
425          ie = ide-1
426          js = jds
427          je = jde-1
428       ELSE
429          ! limit the search to an area around the vortex
430          ! that is limited by max_vortex_speed (default 40) meters per second from
431          ! previous location over vortex_interval (default 15 mins)
432
433          is = max( grid%vc_i - 60 * vortex_interval * max_vortex_speed / dx , 1.0 * ids )
434          js = max( grid%vc_j - 60 * vortex_interval * max_vortex_speed / dx , 1.0 * jds )
435          ie = min( grid%vc_i + 60 * vortex_interval * max_vortex_speed / dx , 1.0 * (ide-1) )
436          je = min( grid%vc_j + 60 * vortex_interval * max_vortex_speed / dx , 1.0 * (jde-1) )
437
438       ENDIF
439
440#  ifdef MOVING_DIAGS
441write(0,*)'search set around last position '
442write(0,*)'   is, ids-1,  ie,  ide-1 ', is, ids-1, ie, ide-1
443write(0,*)'   js, jds-1,  je,  jde-1 ', js, jds-1, je, jde-1
444#  endif
445
446       imploc = -1
447       jmploc = -1
448
449       ! find minimum psfc
450       pmin = 99999999.0     ! make this very large to be sure we find a minumum
451       DO j = js, je
452       DO i = is, ie
453         ! adjust approximately to sea level pressure (same as below: ATCF)
454         psfc(i,j)=psfc(i,j)+11.38*terrain(i,j)
455         IF ( psfc(i,j) .LT. pmin ) THEN
456           pmin = psfc(i,j)
457           imploc = i
458           jmploc = j
459         ENDIF
460       ENDDO
461       ENDDO
462
463       IF ( imploc .EQ. -1 .OR. jmploc .EQ. -1 ) THEN  ! if we fail to find a min there is something seriously wrong
464         WRITE(mess,*)'i,j,is,ie,js,je,imploc,jmploc ',i,j,is,ie,js,je,imploc,jmploc
465         CALL wrf_message(mess)
466         CALL wrf_error_fatal('time_for_move2: Method failure searching for minimum psfc.')
467       ENDIF
468
469       imloc = -1
470       jmloc = -1
471       maxi = -1
472       maxj = -1
473
474       ! find local min, max
475       vmin =  99999999.0
476       vmax = -99999999.0
477       DO j = js, je
478       DO i = is, ie
479         IF ( height(i,j) .LT. vmin ) THEN
480           vmin = height(i,j)
481           imloc = i
482           jmloc = j
483         ENDIF
484         IF ( height(i,j) .GT. vmax ) THEN
485           vmax = height(i,j)
486           maxi = i
487           maxj = j
488         ENDIF
489       ENDDO
490       ENDDO
491
492       IF ( imloc .EQ. -1 .OR. jmloc .EQ. -1 .OR. maxi .EQ. -1 .OR. maxj .EQ. -1 ) THEN
493         WRITE(mess,*)'i,j,is,ie,js,je,imloc,jmloc,maxi,maxj ',i,j,is,ie,js,je,imloc,jmloc,maxi,maxj
494         CALL wrf_message(mess)
495         CALL wrf_error_fatal('time_for_move2: Method failure searching max/min of height.')
496       ENDIF
497
498       fimloc = imloc
499       fjmloc = jmloc
500
501       if ( grid%xi .EQ. -1. ) grid%xi = fimloc
502       if ( grid%xj .EQ. -1. ) grid%xj = fjmloc
503
504       dijsmooth = rsmooth / dx
505
506       fjs = max(fjmloc-dijsmooth,1.0)
507       fje = min(fjmloc+dijsmooth,jde-2.0)
508       fis = max(fimloc-dijsmooth,1.0)
509       fie = min(fimloc+dijsmooth,ide-2.0)
510       js = fjs
511       je = fje
512       is = fis
513       ie = fie
514
515       vmin =  1000000.0
516       vmax = -1000000.0
517       DO j = js, je
518       DO i = is, ie
519         IF ( height(i,j) .GT. vmax ) THEN
520           vmax = height(i,j)
521         ENDIF
522       ENDDO
523       ENDDO
524
525       pbar  = 0.0
526       ipbar = 0.0
527       jpbar = 0.0
528
529       do j=js,je
530          do i=is,ie
531             fact = vmax - height(i,j)
532             pbar  = pbar + fact
533             ipbar = ipbar + fact*(i-is)
534             jpbar = jpbar + fact*(j-js)
535          enddo
536       enddo
537
538      IF ( pbar .NE. 0. ) THEN
539
540!     Compute an adjusted, smoothed, vortex center location in cross
541!     point index space.
542!     Time average. A is coef for old information; B is new
543!     If pbar is zero then just skip this, leave xi and xj alone,
544!     result will be no movement.
545         a = 0.0
546         b = 1.0
547         grid%xi =  (a * grid%xi + b * (is + ipbar / pbar))  / ( a + b )
548         grid%xj =  (a * grid%xj + b * (js + jpbar / pbar))  / ( a + b )
549
550         grid%vc_i = grid%xi + .5
551         grid%vc_j = grid%xj + .5
552
553
554      ENDIF
555
556#  ifdef MOVING_DIAGS
557write(0,*)'computed grid%vc_i, grid%vc_j ',grid%vc_i, grid%vc_j
558i = grid%vc_i ; j = grid%vc_j ; height( i,j ) = height(i,j) * 1.2   !mark the center
559CALL domain_clock_get( grid, current_timestr=message2 )
560WRITE ( message , FMT = '(A," on domain ",I3)' ) TRIM(message2), grid%id
561#  endif
562
563!
564        i = INT(grid%xi+.5)
565        j = INT(grid%xj+.5)
566        write(mess,'("ATCF"," ",A19," ",f8.2," ",f8.2," ",f6.1," ",f6.1)')                &
567                                       timestr(1:19),                               &
568                                       xlat(i,j),                                   &
569                                       xlong(i,j),                                  &
570                                       0.01*pmin,                                   &
571!already computed above                0.01*pmin+0.1138*terrain(imploc,jmploc),     &
572                                       maxws*1.94
573        CALL wrf_message(TRIM(mess))
574                           
575
576
577     ENDIF monitor_only
578
579     DEALLOCATE ( psfc )
580     DEALLOCATE ( xlat )
581     DEALLOCATE ( xlong )
582     DEALLOCATE ( terrain )
583     DEALLOCATE ( height )
584     DEALLOCATE ( height_l )
585
586     CALL wrf_dm_bcast_real( grid%vc_i , 1 )
587     CALL wrf_dm_bcast_real( grid%vc_j , 1 )
588
589     CALL wrf_dm_bcast_real( pmin , 1 )
590     CALL wrf_dm_bcast_integer( imploc , 1 )
591     CALL wrf_dm_bcast_integer( jmploc , 1 )
592
593#  ifdef MOVING_DIAGS
594write(0,*)'after bcast : grid%vc_i, grid%vc_j ',grid%vc_i, grid%vc_j
595#  endif
596
597
598   ENDIF   ! COMPUTE_VORTEX_CENTER_ALARM ringing
599
600#  ifdef MOVING_DIAGS
601write(0,*)'After ENDIF : grid%vc_i, grid%vc_j ',grid%vc_i, grid%vc_j
602#  endif
603
604   dc_i = (ide-ids+1)/2.    ! domain center
605   dc_j = (jde-jds+1)/2.
606
607   disp_x = grid%vc_i - dc_i * 1.0
608   disp_y = grid%vc_j - dc_j * 1.0
609
610#if 0
611! This appears to be an old, redundant, and perhaps even misnamed parameter.
612! Remove it from the namelist and Registry and just hard code it to
613! the default of 6. JM 20050721
614   CALL nl_get_vortex_search_radius( 1, irad )
615#else
616   irad = 6
617#endif
618
619   radius = irad
620
621   if ( disp_x .GT. 0 ) disp_x = min( disp_x , radius )
622   if ( disp_y .GT. 0 ) disp_y = min( disp_y , radius )
623
624   if ( disp_x .LT. 0 ) disp_x = max( disp_x , -radius )
625   if ( disp_y .LT. 0 ) disp_y = max( disp_y , -radius )
626
627   move_cd_x = int ( disp_x  / pgr )
628   move_cd_y = int ( disp_y  / pgr )
629
630   IF ( move_cd_x .GT. 0 ) move_cd_x = min ( move_cd_x , 1 )
631   IF ( move_cd_y .GT. 0 ) move_cd_y = min ( move_cd_y , 1 )
632   IF ( move_cd_x .LT. 0 ) move_cd_x = max ( move_cd_x , -1 )
633   IF ( move_cd_y .LT. 0 ) move_cd_y = max ( move_cd_y , -1 )
634
635   CALL domain_clock_get( grid, current_timestr=timestr )
636
637   WRITE(mess,*)timestr(1:19),' vortex center (in nest x and y): ',grid%vc_i, grid%vc_j
638   CALL wrf_message(TRIM(mess))
639   WRITE(mess,*)timestr(1:19),' grid   center (in nest x and y): ',     dc_i,      dc_j
640   CALL wrf_message(TRIM(mess))
641   WRITE(mess,*)timestr(1:19),' disp          : ',   disp_x,    disp_y
642   CALL wrf_message(TRIM(mess))
643   WRITE(mess,*)timestr(1:19),' move (rel cd) : ',move_cd_x, move_cd_y
644   CALL wrf_message(TRIM(mess))
645
646   grid%vc_i = grid%vc_i - move_cd_x * pgr
647   grid%vc_j = grid%vc_j - move_cd_y * pgr
648
649#  ifdef MOVING_DIAGS
650write(0,*)' changing grid%vc_i,  move_cd_x * pgr ', grid%vc_i, move_cd_x * pgr, move_cd_x, pgr
651#  endif
652
653   IF ( ( move_cd_x .NE. 0 ) .OR. ( move_cd_y .NE. 0 ) ) THEN
654     time_for_move2 = .TRUE.
655   ELSE
656     time_for_move2 = .FALSE.
657   ENDIF
658
659# else
660! from namelist
661   move_cd_x = 0
662   move_cd_y = 0
663   time_for_move2 = .FALSE.
664   CALL domain_clock_get( grid, current_time=ct, start_time=st )
665   CALL nl_get_num_moves( 1, num_moves )
666   IF ( num_moves .GT. max_moves ) THEN
667     WRITE(mess,*)'time_for_moves2: num_moves (',num_moves,') .GT. max_moves (',max_moves,')'
668     CALL wrf_error_fatal( TRIM(mess) )
669   ENDIF
670   DO i = 1, num_moves
671     CALL nl_get_move_id( i, move_id )
672     IF ( move_id .EQ. grid%id ) THEN
673       CALL nl_get_move_interval( i, move_interval )
674       IF ( move_interval .LT. 999999999 ) THEN
675         CALL WRFU_TimeIntervalSet ( ti, M=move_interval, rc=rc )
676         IF ( ct .GE. st + ti ) THEN
677           CALL nl_get_move_cd_x ( i, move_cd_x )
678           CALL nl_get_move_cd_y ( i, move_cd_y )
679           CALL nl_set_move_interval ( i, 999999999 )
680           time_for_move2 = .TRUE.
681           EXIT
682         ENDIF
683       ENDIF
684     ENDIF
685   ENDDO
686# endif
687   RETURN
688#else
689   time_for_move2 = .FALSE.
690#endif
691END FUNCTION time_for_move2
692
693LOGICAL FUNCTION time_for_move ( parent , grid , move_cd_x, move_cd_y )
694   USE module_domain
695   USE module_configure
696   USE module_dm
697USE module_timing
698   USE module_utility
699   IMPLICIT NONE
700! arguments
701   TYPE(domain) , POINTER    :: parent, grid, par, nst
702   INTEGER, INTENT(OUT)      :: move_cd_x , move_cd_y
703#ifdef MOVE_NESTS
704! local
705   INTEGER     :: corral_dist, kid
706   INTEGER     :: dw, de, ds, dn, pgr
707   INTEGER     :: would_move_x, would_move_y
708   INTEGER     :: cids, cide, cjds, cjde, ckds, ckde, &
709                  cims, cime, cjms, cjme, ckms, ckme, &
710                  cips, cipe, cjps, cjpe, ckps, ckpe, &
711                  nids, nide, njds, njde, nkds, nkde, &
712                  nims, nime, njms, njme, nkms, nkme, &
713                  nips, nipe, njps, njpe, nkps, nkpe
714! interface
715   INTERFACE
716     LOGICAL FUNCTION time_for_move2 ( parent , nest , dx , dy )
717        USE module_domain
718        USE module_utility
719        TYPE(domain) , POINTER    :: parent , nest
720        INTEGER, INTENT(OUT)      :: dx , dy
721     END FUNCTION time_for_move2
722   END INTERFACE
723! executable
724!
725! Simplifying assumption: domains in moving nest simulations have only
726! one parent and only one child.
727
728   IF   ( grid%num_nests .GT. 1 ) THEN
729     CALL wrf_error_fatal ( 'domains in moving nest simulations can have only 1 nest' )
730   ENDIF
731   kid = 1
732!
733! find out if this is the innermost nest (will not have kids)
734   IF   ( grid%num_nests .EQ. 0 ) THEN
735     ! code that executes on innermost nest
736     time_for_move = time_for_move2 ( parent , grid , move_cd_x, move_cd_y )
737
738     ! Make sure the parent can move before allowing the nest to approach
739     ! its boundary
740     par => grid%parents(1)%ptr
741     nst => grid
742
743     would_move_x = move_cd_x
744     would_move_y = move_cd_y
745
746     ! top of until loop
747100  CONTINUE
748       CALL nl_get_corral_dist ( nst%id , corral_dist )
749       CALL get_ijk_from_grid (  nst ,                               &
750                                 nids, nide, njds, njde, nkds, nkde, &
751                                 nims, nime, njms, njme, nkms, nkme, &
752                                 nips, nipe, njps, njpe, nkps, nkpe  )
753       CALL get_ijk_from_grid (  par ,                               &
754                                 cids, cide, cjds, cjde, ckds, ckde, &
755                                 cims, cime, cjms, cjme, ckms, ckme, &
756                                 cips, cipe, cjps, cjpe, ckps, ckpe  )
757       CALL nl_get_parent_grid_ratio ( nst%id , pgr )
758       ! perform measurements...
759       !  from western boundary
760       dw = nst%i_parent_start + would_move_x - cids
761       !  from southern boundary
762       ds = nst%j_parent_start + would_move_y - cjds
763       !  from eastern boundary
764       de = cide - ( nst%i_parent_start + (nide-nids+1)/pgr + would_move_x )
765       !  from northern boundary
766       dn = cjde - ( nst%j_parent_start + (njde-njds+1)/pgr + would_move_y )
767
768       ! would this generate a move on the parent?
769       would_move_x = 0
770       would_move_y = 0
771       if ( dw .LE. corral_dist ) would_move_x = would_move_x - 1
772       if ( de .LE. corral_dist ) would_move_x = would_move_x + 1
773       if ( ds .LE. corral_dist ) would_move_y = would_move_y - 1
774       if ( dn .LE. corral_dist ) would_move_y = would_move_y + 1
775
776     IF ( par%id .EQ. 1 ) THEN
777         IF ( would_move_x .NE. 0 .AND. move_cd_x .NE. 0 ) THEN
778           CALL wrf_message('MOAD can not move. Cancelling nest move in X')
779           if ( grid%num_nests .eq. 0 ) grid%vc_i = grid%vc_i + move_cd_x * pgr  ! cancel effect of move
780           move_cd_x = 0
781         ENDIF
782         IF ( would_move_y .NE. 0 .AND. move_cd_y .NE. 0 ) THEN
783           CALL wrf_message('MOAD can not move. Cancelling nest move in Y')
784           if ( grid%num_nests .eq. 0 ) grid%vc_j = grid%vc_j + move_cd_y * pgr  ! cancel effect of move
785           move_cd_y = 0
786         ENDIF
787     ELSE
788         nst => par
789         par => nst%parents(1)%ptr
790         GOTO 100
791     ENDIF
792
793! bottom of until loop
794     time_for_move = ( move_cd_x .NE. 0 ) .OR. ( move_cd_y .NE. 0 )
795
796   ELSE
797     ! code that executes on parent to see if parent needs to move
798     ! get closest number of cells we'll allow nest edge to approach parent bdy
799     CALL nl_get_corral_dist ( grid%nests(kid)%ptr%id , corral_dist )
800     ! get dims
801     CALL get_ijk_from_grid (  grid%nests(kid)%ptr ,               &
802                               nids, nide, njds, njde, nkds, nkde, &
803                               nims, nime, njms, njme, nkms, nkme, &
804                               nips, nipe, njps, njpe, nkps, nkpe  )
805     CALL get_ijk_from_grid (  grid ,                              &
806                               cids, cide, cjds, cjde, ckds, ckde, &
807                               cims, cime, cjms, cjme, ckms, ckme, &
808                               cips, cipe, cjps, cjpe, ckps, ckpe  )
809     CALL nl_get_parent_grid_ratio ( grid%nests(kid)%ptr%id , pgr )
810     ! perform measurements...
811     !  from western boundary
812     dw = grid%nests(kid)%ptr%i_parent_start - 1
813     !  from southern boundary
814     ds = grid%nests(kid)%ptr%j_parent_start - 1
815     !  from eastern boundary
816     de = cide - ( grid%nests(kid)%ptr%i_parent_start + (nide-nids+1)/pgr )
817     !  from northern boundary
818     dn = cjde - ( grid%nests(kid)%ptr%j_parent_start + (njde-njds+1)/pgr )
819
820     ! move this domain (the parent containing the moving nest)
821     ! in a direction that reestablishes the distance from
822     ! the boundary.
823     move_cd_x = 0
824     move_cd_y = 0
825     if ( dw .LE. corral_dist ) move_cd_x = move_cd_x - 1
826     if ( de .LE. corral_dist ) move_cd_x = move_cd_x + 1
827     if ( ds .LE. corral_dist ) move_cd_y = move_cd_y - 1
828     if ( dn .LE. corral_dist ) move_cd_y = move_cd_y + 1
829
830     time_for_move = ( move_cd_x .NE. 0 ) .OR. ( move_cd_y .NE. 0 )
831
832     IF ( time_for_move ) THEN
833       IF ( grid%id .EQ. 1 ) THEN
834
835         CALL wrf_message( 'DANGER: Nest has moved too close to boundary of outermost domain.' )
836         time_for_move = .FALSE.
837
838       ELSE
839         ! need to adjust the intermediate domain of the nest in relation to this
840         ! domain since we're moving
841
842         CALL wrf_dm_move_nest ( grid , grid%nests(kid)%ptr%intermediate_grid , -move_cd_x*pgr, -move_cd_y*pgr )
843         CALL adjust_domain_dims_for_move( grid%nests(kid)%ptr%intermediate_grid , -move_cd_x*pgr, -move_cd_y*pgr )
844         grid%nests(kid)%ptr%i_parent_start = grid%nests(kid)%ptr%i_parent_start - move_cd_x*pgr
845         CALL nl_set_i_parent_start( grid%nests(kid)%ptr%id, grid%nests(kid)%ptr%i_parent_start )
846         grid%nests(kid)%ptr%j_parent_start = grid%nests(kid)%ptr%j_parent_start - move_cd_y*pgr
847         CALL nl_set_j_parent_start( grid%nests(kid)%ptr%id, grid%nests(kid)%ptr%j_parent_start )
848
849       ENDIF
850     ENDIF
851
852   ENDIF
853
854   RETURN
855#else
856   time_for_move = .FALSE.
857#endif
858END FUNCTION time_for_move
859
860! Put any tests for non-moving options or conditions in here
861LOGICAL FUNCTION should_not_move ( id )
862  USE module_state_description
863  USE module_configure
864  IMPLICIT NONE
865  INTEGER, INTENT(IN) :: id
866 ! Local
867  LOGICAL retval
868  INTEGER cu_physics, ra_sw_physics, ra_lw_physics, ucmcall, obs_nudge_opt
869
870  retval = .FALSE.
871! check for GD ensemble cumulus, which can not move
872  CALL nl_get_cu_physics( id , cu_physics )
873  IF ( cu_physics .EQ. GDSCHEME ) THEN
874    CALL wrf_message('Grell cumulus can not be specified with moving nests. Movement disabled.')
875    retval = .TRUE.
876  ENDIF
877! check for CAM radiation scheme , which can not move
878  CALL nl_get_ra_sw_physics( id , ra_sw_physics )
879  IF ( ra_sw_physics .EQ. CAMSWSCHEME ) THEN
880    CALL wrf_message('CAM SW radiation can not be specified with moving nests. Movement disabled.')
881    retval = .TRUE.
882  ENDIF
883  CALL nl_get_ra_lw_physics( id , ra_lw_physics )
884  IF ( ra_lw_physics .EQ. CAMLWSCHEME ) THEN
885    CALL wrf_message('CAM LW radiation can not be specified with moving nests. Movement disabled.')
886    retval = .TRUE.
887  ENDIF
888! check for urban canopy Noah LSM, which can not move
889  CALL nl_get_ucmcall( id , ucmcall )
890  IF ( ucmcall .EQ. 1 ) THEN
891    CALL wrf_message('UCM Noah LSM can not be specified with moving nests. Movement disabled.')
892    retval = .TRUE.
893  ENDIF
894! check for observation nudging, which can not move
895  CALL nl_get_obs_nudge_opt( id , obs_nudge_opt )
896  IF ( obs_nudge_opt .EQ. 1 ) THEN
897    CALL wrf_message('Observation nudging can not be specified with moving nests. Movement disabled.')
898    retval = .TRUE.
899  ENDIF
900  should_not_move = retval
901END FUNCTION
902
Note: See TracBrowser for help on using the repository browser.