module_bc.F

References to this file elsewhere.
1 !WRF:MODEL_LAYER:BOUNDARY
2 !
3 
4 MODULE module_bc
5 
6    USE module_configure
7    USE module_wrf_error
8    IMPLICIT NONE
9 
10 !   TYPE bcs
11 !
12 !     LOGICAL                     :: periodic_x
13 !     LOGICAL                     :: symmetric_xs
14 !     LOGICAL                     :: symmetric_xe
15 !     LOGICAL                     :: open_xs
16 !     LOGICAL                     :: open_xe
17 !     LOGICAL                     :: periodic_y
18 !     LOGICAL                     :: symmetric_ys
19 !     LOGICAL                     :: symmetric_ye
20 !     LOGICAL                     :: open_ys
21 !     LOGICAL                     :: open_ye
22 !     LOGICAL                     :: nested
23 !     LOGICAL                     :: specified
24 !     LOGICAL                     :: top_radiation
25 !
26 !   END TYPE bcs
27 
28 !  set the bdyzone.  We are hardwiring this here and we'll
29 !  decide later where it should be set and stored
30 
31    INTEGER, PARAMETER            :: bdyzone = 4
32    INTEGER, PARAMETER            :: bdyzone_x = bdyzone
33    INTEGER, PARAMETER            :: bdyzone_y = bdyzone
34 
35    INTERFACE stuff_bdy
36      MODULE PROCEDURE stuff_bdy_new , stuff_bdy_old
37    END INTERFACE
38 
39    INTERFACE stuff_bdytend
40      MODULE PROCEDURE stuff_bdytend_new , stuff_bdytend_old
41    END INTERFACE
42 
43 CONTAINS
44 
45   SUBROUTINE boundary_condition_check ( config_flags, bzone, error, gn )
46 
47 !  this routine checks the boundary condition logicals 
48 !  to make sure that the boundary conditions are not over
49 !  or under specified.  The routine also checks that the
50 !  boundary zone is sufficiently sized for the specified
51 !  boundary conditions
52 
53   IMPLICIT NONE
54 
55   TYPE( grid_config_rec_type ) config_flags
56 
57   INTEGER, INTENT(IN   ) :: bzone, gn
58   INTEGER, INTENT(INOUT) :: error
59 
60 ! local variables
61 
62   INTEGER :: xs_bc, xe_bc, ys_bc, ye_bc, bzone_min
63   INTEGER :: nprocx, nprocy
64 
65   CALL wrf_debug( 100 , ' checking boundary conditions for grid ' )
66 
67   error = 0
68   xs_bc = 0
69   xe_bc = 0
70   ys_bc = 0
71   ye_bc = 0
72 
73 !  sum the number of conditions specified for each lateral boundary.
74 !  obviously, this number should be 1
75 
76   IF( config_flags%periodic_x ) THEN
77     xs_bc = xs_bc+1
78     xe_bc = xe_bc+1
79   ENDIF
80 
81   IF( config_flags%periodic_y ) THEN
82     ys_bc = ys_bc+1
83     ye_bc = ye_bc+1
84   ENDIF
85 
86   IF( config_flags%symmetric_xs ) xs_bc = xs_bc + 1
87   IF( config_flags%symmetric_xe ) xe_bc = xe_bc + 1
88   IF( config_flags%open_xs )      xs_bc = xs_bc + 1
89   IF( config_flags%open_xe )      xe_bc = xe_bc + 1
90 
91 
92   IF( config_flags%symmetric_ys ) ys_bc = ys_bc + 1
93   IF( config_flags%symmetric_ye ) ye_bc = ye_bc + 1
94   IF( config_flags%open_ys )      ys_bc = ys_bc + 1
95   IF( config_flags%open_ye )      ye_bc = ye_bc + 1
96 
97   IF( config_flags%nested ) THEN
98      xs_bc = xs_bc + 1
99      xe_bc = xe_bc + 1
100      ys_bc = ys_bc + 1
101      ye_bc = ye_bc + 1
102    ENDIF
103 
104   IF( config_flags%specified ) THEN
105      IF( .NOT. config_flags%periodic_x)xs_bc = xs_bc + 1
106      IF( .NOT. config_flags%periodic_x)xe_bc = xe_bc + 1
107      ys_bc = ys_bc + 1
108      ye_bc = ye_bc + 1
109    ENDIF
110 
111 !  check the number of conditions for each boundary
112 
113    IF( (xs_bc /= 1) .or. &
114        (xe_bc /= 1) .or. &
115        (ys_bc /= 1) .or. &
116        (ye_bc /= 1)         ) THEN
117 
118      error = 1
119 
120      write( wrf_err_message ,*) ' *** Error in boundary condition specification '
121      CALL wrf_message ( wrf_err_message )
122      write( wrf_err_message ,*) ' boundary conditions at xs ', xs_bc
123      CALL wrf_message ( wrf_err_message )
124      write( wrf_err_message ,*) ' boundary conditions at xe ', xe_bc
125      CALL wrf_message ( wrf_err_message )
126      write( wrf_err_message ,*) ' boundary conditions at ys ', ys_bc
127      CALL wrf_message ( wrf_err_message )
128      write( wrf_err_message ,*) ' boundary conditions at ye ', ye_bc
129      CALL wrf_message ( wrf_err_message )
130      write( wrf_err_message ,*) ' boundary conditions logicals are '
131      CALL wrf_message ( wrf_err_message )
132      write( wrf_err_message ,*) ' periodic_x   ',config_flags%periodic_x
133      CALL wrf_message ( wrf_err_message )
134      write( wrf_err_message ,*) ' periodic_y   ',config_flags%periodic_y
135      CALL wrf_message ( wrf_err_message )
136      write( wrf_err_message ,*) ' symmetric_xs ',config_flags%symmetric_xs
137      CALL wrf_message ( wrf_err_message )
138      write( wrf_err_message ,*) ' symmetric_xe ',config_flags%symmetric_xe
139      CALL wrf_message ( wrf_err_message )
140      write( wrf_err_message ,*) ' symmetric_ys ',config_flags%symmetric_ys
141      CALL wrf_message ( wrf_err_message )
142      write( wrf_err_message ,*) ' symmetric_ye ',config_flags%symmetric_ye
143      CALL wrf_message ( wrf_err_message )
144      write( wrf_err_message ,*) ' open_xs      ',config_flags%open_xs
145      CALL wrf_message ( wrf_err_message )
146      write( wrf_err_message ,*) ' open_xe      ',config_flags%open_xe
147      CALL wrf_message ( wrf_err_message )
148      write( wrf_err_message ,*) ' open_ys      ',config_flags%open_ys
149      CALL wrf_message ( wrf_err_message )
150      write( wrf_err_message ,*) ' open_ye      ',config_flags%open_ye
151      CALL wrf_message ( wrf_err_message )
152      write( wrf_err_message ,*) ' nested       ',config_flags%nested
153      CALL wrf_message ( wrf_err_message )
154      write( wrf_err_message ,*) ' specified    ',config_flags%specified
155      CALL wrf_message ( wrf_err_message )
156      CALL wrf_error_fatal( ' *** Error in boundary condition specification ' )
157 
158    ENDIF
159 
160 !  now check to see if boundary zone size is sufficient.
161 !  we could have the necessary boundary zone size be returned
162 !  to the calling routine.
163 
164    IF( config_flags%periodic_x   .or. &
165        config_flags%periodic_y   .or. &
166        config_flags%symmetric_xs .or. &
167        config_flags%symmetric_xe .or. &
168        config_flags%symmetric_ys .or. &
169        config_flags%symmetric_ye        )  THEN
170 
171        bzone_min = MAX( 1,                                  &
172                         (config_flags%h_mom_adv_order+1)/2, &
173                         (config_flags%h_sca_adv_order+1)/2 )
174 
175        IF( bzone < bzone_min) THEN  
176 
177          error = 2
178          WRITE ( wrf_err_message , * ) ' boundary zone not large enough '
179          CALL wrf_message ( wrf_err_message )
180          WRITE ( wrf_err_message , * ) ' boundary zone specified      ',bzone
181          CALL wrf_message ( wrf_err_message )
182          WRITE ( wrf_err_message , * ) ' minimum boundary zone needed ',bzone_min
183          CALL wrf_error_fatal ( wrf_err_message )
184 
185        ENDIF
186    ENDIF
187 
188    CALL wrf_debug ( 100 , ' boundary conditions OK for grid ' )
189 
190    END subroutine boundary_condition_check
191 
192 !--------------------------------------------------------------------------
193    SUBROUTINE set_physical_bc2d( dat, variable_in,  &
194                                  config_flags,           & 
195                                  ids,ide, jds,jde,   & ! domain dims
196                                  ims,ime, jms,jme,   & ! memory dims
197                                  ips,ipe, jps,jpe,   & ! patch  dims
198                                  its,ite, jts,jte   )      
199 
200 !  This subroutine sets the data in the boundary region, by direct
201 !  assignment if possible, for periodic and symmetric (wall)
202 !  boundary conditions.  Currently, we are only doing 1 variable
203 !  at a time - lots of overhead, so maybe this routine can be easily 
204 !  inlined later or we could pass multiple variables -
205 !  would probably want a largestep and smallstep version.
206 
207 !  15 Jan 99, Dave
208 !  Modified the incoming its,ite,jts,jte to truly be the tile size.
209 !  This required modifying the loop limits when the "istag" or "jstag"
210 !  is used, as this is only required at the end of the domain.
211 
212       IMPLICIT NONE
213 
214       INTEGER,      INTENT(IN   )    :: ids,ide, jds,jde
215       INTEGER,      INTENT(IN   )    :: ims,ime, jms,jme
216       INTEGER,      INTENT(IN   )    :: ips,ipe, jps,jpe
217       INTEGER,      INTENT(IN   )    :: its,ite, jts,jte
218       CHARACTER,    INTENT(IN   )    :: variable_in
219 
220       CHARACTER                      :: variable
221 
222       REAL,  DIMENSION( ims:ime , jms:jme ) :: dat
223       TYPE( grid_config_rec_type ) config_flags
224 
225       INTEGER  :: i, j, istag, jstag, itime
226 
227       LOGICAL  :: debug, open_bc_copy
228 
229 !------------
230 
231       debug = .false.
232 
233       open_bc_copy = .false.
234 
235       variable = variable_in
236       IF ( variable_in .ge. 'A' .and. variable_in .le. 'Z' ) THEN
237         variable = CHAR( ICHAR(variable_in) - ICHAR('A') + ICHAR('a') )
238       ENDIF
239       IF ((variable == 'u') .or. (variable == 'v') .or.  &
240           (variable == 'w') .or. (variable == 't') .or.  &
241           (variable == 'x') .or. (variable == 'y') .or.  &
242           (variable == 'r') .or. (variable == 'p') ) open_bc_copy = .true.
243 
244 !  begin, first set a staggering variable
245 
246       istag = -1
247       jstag = -1
248 
249       IF ((variable == 'u') .or. (variable == 'x')) istag = 0
250       IF ((variable == 'v') .or. (variable == 'y')) jstag = 0
251 
252       if(debug) then
253         write(6,*) ' in bc2d, var is ',variable, istag, jstag
254         write(6,*) ' b.cs are ',  &
255       config_flags%periodic_x,  &
256       config_flags%periodic_y
257       end if
258       
259 
260 
261 !  periodic conditions.
262 !  note, patch must cover full range in periodic dir, or else
263 !  its intra-patch communication that is handled elsewheres.
264 !  symmetry conditions can always be handled here, because no
265 !  outside patch communication is needed
266 
267       periodicity_x:  IF( ( config_flags%periodic_x ) ) THEN 
268         IF ( ( ids == ips ) .and.  ( ide == ipe ) ) THEN  ! test if east and west both on-processor 
269           IF ( its == ids ) THEN
270 
271             DO j = MAX(jds,jts-1), MIN(jte+1,jde+jstag)
272             DO i = 0,-(bdyzone-1),-1
273               dat(ids+i-1,j) = dat(ide+i-1,j)
274             ENDDO
275             ENDDO
276 
277           ENDIF
278 
279           IF ( ite == ide ) THEN
280 
281             DO j = MAX(jds,jts-1), MIN(jte+1,jde+jstag)
282 !!          DO i = 1 , bdyzone
283             DO i = -istag , bdyzone
284               dat(ide+i+istag,j) = dat(ids+i+istag,j)
285             ENDDO
286             ENDDO
287 
288           ENDIF
289         ENDIF
290 
291       ELSE 
292 
293         symmetry_xs: IF( ( config_flags%symmetric_xs ) .and.  &
294                          ( its == ids )                  )  THEN
295 
296           IF ( (variable /= 'u') .and. (variable /= 'x') ) THEN
297 
298             DO j = MAX(jds,jts-1), MIN(jte+1,jde+jstag)
299             DO i = 1, bdyzone
300               dat(ids-i,j) = dat(ids+i-1,j) !  here, dat(0) = dat(1), etc
301             ENDDO                             !  symmetry about dat(0.5) (u=0 pt)
302             ENDDO
303 
304           ELSE
305 
306             IF( variable == 'u' ) THEN
307 
308               DO j = MAX(jds,jts-1), MIN(jte+1,jde+jstag)
309               DO i = 0, bdyzone-1
310                 dat(ids-i,j) = - dat(ids+i,j) ! here, u(0) = - u(2), etc
311               ENDDO                             !  normal b.c symmetry at u(1)
312               ENDDO
313 
314             ELSE
315 
316               DO j = MAX(jds,jts-1), MIN(jte+1,jde+jstag)
317               DO i = 0, bdyzone-1
318                 dat(ids-i,j) =   dat(ids+i,j) ! here, phi(0) = phi(2), etc
319               ENDDO                             !  normal b.c symmetry at phi(1)
320               ENDDO
321 
322             END IF
323 
324           ENDIF
325 
326         ENDIF symmetry_xs
327 
328 
329 !  now the symmetry boundary at xe
330 
331         symmetry_xe: IF( ( config_flags%symmetric_xe ) .and.  &
332                          ( ite == ide )                  )  THEN
333 
334           IF ( (variable /= 'u') .and. (variable /= 'x') ) THEN
335 
336             DO j = MAX(jds,jts-1), MIN(jte+1,jde+jstag)
337             DO i = 1, bdyzone
338               dat(ide+i-1,j) = dat(ide-i,j)  !  sym. about dat(ide-0.5)
339             ENDDO
340             ENDDO
341 
342           ELSE
343 
344             IF (variable == 'u' ) THEN
345 
346               DO j = MAX(jds,jts-1), MIN(jte+1,jde+jstag)
347               DO i = 0, bdyzone-1
348                 dat(ide+i,j) = - dat(ide-i,j)  ! u(ide+1) = - u(ide-1), etc.
349               ENDDO
350               ENDDO
351 
352 
353             ELSE
354 
355               DO j = MAX(jds,jts-1), MIN(jte+1,jde+jstag)
356               DO i = 0, bdyzone-1
357                 dat(ide+i,j) = dat(ide-i,j)  !  phi(ide+1) = phi(ide-1), etc.
358               ENDDO
359               ENDDO
360 
361             END IF
362 
363           END IF 
364 
365         END IF symmetry_xe
366 
367 
368 !  set open b.c in X copy into boundary zone here.  WCS, 19 March 2000
369 
370         open_xs: IF( ( config_flags%open_xs   .or. &
371                        config_flags%specified .or. &
372                        config_flags%nested            ) .and.  &
373                          ( its == ids ) .and. open_bc_copy  )  THEN
374 
375             DO j = MAX(jds,jts-1), MIN(jte+1,jde+jstag)
376               dat(ids-1,j) = dat(ids,j) !  here, dat(0) = dat(1)
377               dat(ids-2,j) = dat(ids,j)
378               dat(ids-3,j) = dat(ids,j)
379             ENDDO
380 
381         ENDIF open_xs
382 
383 
384 !  now the open boundary copy at xe
385 
386         open_xe: IF( ( config_flags%open_xe   .or. &
387                        config_flags%specified .or. &
388                        config_flags%nested            ) .and.  &
389                           ( ite == ide ) .and. open_bc_copy  )  THEN
390 
391           IF ( variable /= 'u' .and. variable /= 'x') THEN
392 
393             DO j = MAX(jds,jts-1), MIN(jte+1,jde+jstag)
394               dat(ide  ,j) = dat(ide-1,j) 
395               dat(ide+1,j) = dat(ide-1,j) 
396               dat(ide+2,j) = dat(ide-1,j) 
397             ENDDO
398 
399           ELSE
400 
401             DO j = MAX(jds,jts-1), MIN(jte+1,jde+jstag)
402               dat(ide+1,j) = dat(ide,j)
403               dat(ide+2,j) = dat(ide,j)
404               dat(ide+3,j) = dat(ide,j)
405             ENDDO
406 
407           END IF 
408 
409         END IF open_xe
410 
411 !  end open b.c in X copy into boundary zone addition.  WCS, 19 March 2000
412 
413       END IF periodicity_x
414 
415 !  same procedure in y
416 
417       periodicity_y:  IF( ( config_flags%periodic_y ) ) THEN
418         IF ( ( jds == jps ) .and. ( jde == jpe ) )  THEN    ! test of both north and south on processor
419 
420           IF( jts == jds ) then
421 
422             DO j = 0, -(bdyzone-1), -1
423             DO i = MAX(ids,its-1), MIN(ite+1,ide+istag)
424               dat(i,jds+j-1) = dat(i,jde+j-1)
425             ENDDO
426             ENDDO
427 
428           END IF
429 
430           IF( jte == jde ) then
431 
432             DO j = -jstag, bdyzone
433             DO i = MAX(ids,its-1), MIN(ite+1,ide+istag)
434               dat(i,jde+j+jstag) = dat(i,jds+j+jstag)
435             ENDDO
436             ENDDO
437 
438           END IF
439 
440         END IF
441 
442       ELSE
443 
444         symmetry_ys: IF( ( config_flags%symmetric_ys ) .and.  &
445                          ( jts == jds)                   )  THEN
446 
447           IF ( (variable /= 'v') .and. (variable /= 'y') ) THEN
448 
449             DO j = 1, bdyzone
450             DO i = MAX(ids,its-1), MIN(ite+1,ide+istag)
451               dat(i,jds-j) = dat(i,jds+j-1) 
452             ENDDO
453             ENDDO
454 
455           ELSE
456 
457             IF (variable == 'v') THEN
458 
459               DO j = 1, bdyzone
460               DO i = MAX(ids,its-1), MIN(ite+1,ide+istag)
461                 dat(i,jds-j) = - dat(i,jds+j) 
462               ENDDO              
463               ENDDO
464 
465             ELSE
466 
467               DO j = 1, bdyzone
468               DO i = MAX(ids,its-1), MIN(ite+1,ide+istag)
469                 dat(i,jds-j) = dat(i,jds+j) 
470               ENDDO              
471               ENDDO
472 
473             END IF
474 
475           ENDIF
476 
477         ENDIF symmetry_ys
478 
479 !  now the symmetry boundary at ye
480 
481         symmetry_ye: IF( ( config_flags%symmetric_ye ) .and.  &
482                          ( jte == jde )                  )  THEN
483 
484           IF ( (variable /= 'v') .and. (variable /= 'y') ) THEN
485 
486             DO j = 1, bdyzone
487             DO i = MAX(ids,its-1), MIN(ite+1,ide+istag)
488               dat(i,jde+j-1) = dat(i,jde-j) 
489             ENDDO                               
490             ENDDO
491 
492           ELSE
493 
494             IF (variable == 'v' ) THEN
495 
496               DO j = 1, bdyzone
497               DO i = MAX(ids,its-1), MIN(ite+1,ide+istag)
498                 dat(i,jde+j) = - dat(i,jde-j)    ! bugfix: changed jds on rhs to jde , JM 20020410
499               ENDDO                               
500               ENDDO
501 
502             ELSE
503 
504               DO j = 1, bdyzone
505               DO i = MAX(ids,its-1), MIN(ite+1,ide+istag)
506                 dat(i,jde+j) = dat(i,jde-j)
507               ENDDO                               
508               ENDDO
509 
510             END IF
511 
512           ENDIF
513 
514         END IF symmetry_ye
515 
516 !  set open b.c in Y copy into boundary zone here.  WCS, 19 March 2000
517 
518         open_ys: IF( ( config_flags%open_ys   .or. &
519                        config_flags%specified .or. &
520                        config_flags%nested            ) .and.  &
521                          ( jts == jds) .and. open_bc_copy )  THEN
522 
523             DO i = MAX(ids,its-1), MIN(ite+1,ide+istag)
524               dat(i,jds-1) = dat(i,jds) 
525               dat(i,jds-2) = dat(i,jds) 
526               dat(i,jds-3) = dat(i,jds) 
527             ENDDO
528 
529         ENDIF open_ys
530 
531 !  now the open boundary copy at ye
532 
533         open_ye: IF( ( config_flags%open_ye   .or. &
534                        config_flags%specified .or. &
535                        config_flags%nested            ) .and.  &
536                          ( jte == jde ) .and. open_bc_copy )  THEN
537 
538           IF  (variable /= 'v' .and. variable /= 'y' ) THEN
539 
540             DO i = MAX(ids,its-1), MIN(ite+1,ide+istag)
541               dat(i,jde  ) = dat(i,jde-1) 
542               dat(i,jde+1) = dat(i,jde-1) 
543               dat(i,jde+2) = dat(i,jde-1) 
544             ENDDO                               
545 
546           ELSE
547 
548             DO i = MAX(ids,its-1), MIN(ite+1,ide+istag)
549               dat(i,jde+1) = dat(i,jde) 
550               dat(i,jde+2) = dat(i,jde) 
551               dat(i,jde+3) = dat(i,jde) 
552             ENDDO                               
553 
554           ENDIF
555 
556         END IF open_ye
557       
558 !  end open b.c in Y copy into boundary zone addition.  WCS, 19 March 2000
559 
560       END IF periodicity_y
561 
562 !  fix corners for doubly periodic domains
563 
564       IF ( config_flags%periodic_x .and. config_flags%periodic_y &
565            .and. (ids == ips) .and. (ide == ipe)                 &
566            .and. (jds == jps) .and. (jde == jpe)                   ) THEN
567 
568          IF ( (its == ids) .and. (jts == jds) ) THEN  ! lower left corner fill
569            DO j = 0, -(bdyzone-1), -1
570            DO i = 0, -(bdyzone-1), -1
571              dat(ids+i-1,jds+j-1) = dat(ide+i-1,jde+j-1)
572            ENDDO
573            ENDDO
574          END IF
575 
576          IF ( (ite == ide) .and. (jts == jds) ) THEN  ! lower right corner fill
577            DO j = 0, -(bdyzone-1), -1
578            DO i = 1, bdyzone
579              dat(ide+i+istag,jds+j-1) = dat(ids+i+istag,jde+j-1)
580            ENDDO
581            ENDDO
582          END IF
583 
584          IF ( (ite == ide) .and. (jte == jde) ) THEN  ! upper right corner fill
585            DO j = 1, bdyzone
586            DO i = 1, bdyzone
587              dat(ide+i+istag,jde+j+jstag) = dat(ids+i+istag,jds+j+jstag)
588            ENDDO
589            ENDDO
590          END IF
591 
592          IF ( (its == ids) .and. (jte == jde) ) THEN  ! upper left corner fill
593            DO j = 1, bdyzone
594            DO i = 0, -(bdyzone-1), -1
595              dat(ids+i-1,jde+j+jstag) = dat(ide+i-1,jds+j+jstag)
596            ENDDO
597            ENDDO
598          END IF
599 
600        END IF
601 
602    END SUBROUTINE set_physical_bc2d
603 
604 !-----------------------------------
605 
606    SUBROUTINE set_physical_bc3d( dat, variable_in,        &
607                                config_flags,                   & 
608                                ids,ide, jds,jde, kds,kde,  & ! domain dims
609                                ims,ime, jms,jme, kms,kme,  & ! memory dims
610                                ips,ipe, jps,jpe, kps,kpe,  & ! patch  dims
611                                its,ite, jts,jte, kts,kte )
612 
613 !  This subroutine sets the data in the boundary region, by direct
614 !  assignment if possible, for periodic and symmetric (wall)
615 !  boundary conditions.  Currently, we are only doing 1 variable
616 !  at a time - lots of overhead, so maybe this routine can be easily 
617 !  inlined later or we could pass multiple variables -
618 !  would probably want a largestep and smallstep version.
619 
620 !  15 Jan 99, Dave
621 !  Modified the incoming its,ite,jts,jte to truly be the tile size.
622 !  This required modifying the loop limits when the "istag" or "jstag"
623 !  is used, as this is only required at the end of the domain.
624 
625       IMPLICIT NONE
626 
627       INTEGER,      INTENT(IN   )    :: ids,ide, jds,jde, kds,kde
628       INTEGER,      INTENT(IN   )    :: ims,ime, jms,jme, kms,kme
629       INTEGER,      INTENT(IN   )    :: ips,ipe, jps,jpe, kps,kpe
630       INTEGER,      INTENT(IN   )    :: its,ite, jts,jte, kts,kte
631       CHARACTER,    INTENT(IN   )    :: variable_in
632 
633       CHARACTER                      :: variable
634 
635       REAL,  DIMENSION( ims:ime , kms:kme , jms:jme ) :: dat
636       TYPE( grid_config_rec_type ) config_flags
637 
638       INTEGER  :: i, j, k, istag, jstag, itime, k_end
639 
640       LOGICAL  :: debug, open_bc_copy
641 
642 !------------
643 
644       debug = .false.
645 
646       open_bc_copy = .false.
647 
648       variable = variable_in
649       IF ( variable_in .ge. 'A' .and. variable_in .le. 'Z' ) THEN
650         variable = CHAR( ICHAR(variable_in) - ICHAR('A') + ICHAR('a') )
651       ENDIF
652 
653       IF ((variable == 'u') .or. (variable == 'v') .or.     &
654           (variable == 'w') .or. (variable == 't') .or.     &
655           (variable == 'd') .or. (variable == 'e') .or. &
656           (variable == 'x') .or. (variable == 'y') .or. &
657           (variable == 'f') .or. (variable == 'r') .or. &
658           (variable == 'p')                        ) open_bc_copy = .true.
659 
660 !  begin, first set a staggering variable
661 
662       istag = -1
663       jstag = -1
664       k_end = max(1,min(kde-1,kte))
665 
666 
667       IF ((variable == 'u') .or. (variable == 'x')) istag = 0
668       IF ((variable == 'v') .or. (variable == 'y')) jstag = 0
669       IF ((variable == 'd') .or. (variable == 'xy')) then
670          istag = 0
671          jstag = 0
672       ENDIF
673       IF ((variable == 'e') ) then
674          istag = 0
675          k_end = min(kde,kte)
676       ENDIF
677 
678       IF ((variable == 'f') ) then
679          jstag = 0
680          k_end = min(kde,kte)
681       ENDIF
682 
683       IF ( variable == 'w')  k_end = min(kde,kte)
684 
685 !      k_end = kte
686 
687       if(debug) then
688         write(6,*) ' in bc, var is ',variable, istag, jstag, kte, k_end
689         write(6,*) ' b.cs are ',  &
690       config_flags%periodic_x,  &
691       config_flags%periodic_y
692       end if
693       
694 
695 
696 !  periodic conditions.
697 !  note, patch must cover full range in periodic dir, or else
698 !  its intra-patch communication that is handled elsewheres.
699 !  symmetry conditions can always be handled here, because no
700 !  outside patch communication is needed
701 
702       periodicity_x:  IF( ( config_flags%periodic_x ) ) THEN
703 
704         IF ( ( ids == ips ) .and. ( ide == ipe ) ) THEN  ! test if both east and west on-processor
705           IF ( its == ids ) THEN
706 
707             DO j = MAX(jds,jts-1), MIN(jte+1,jde+jstag)
708             DO k = kts, k_end
709             DO i = 0,-(bdyzone-1),-1
710               dat(ids+i-1,k,j) = dat(ide+i-1,k,j)
711             ENDDO
712             ENDDO
713             ENDDO
714 
715           ENDIF
716 
717 
718           IF ( ite == ide ) THEN
719 
720             DO j = MAX(jds,jts-1), MIN(jte+1,jde+jstag)
721             DO k = kts, k_end
722             DO i = -istag , bdyzone
723               dat(ide+i+istag,k,j) = dat(ids+i+istag,k,j)
724             ENDDO
725             ENDDO
726             ENDDO
727 
728           ENDIF
729 
730         ENDIF
731 
732       ELSE 
733 
734         symmetry_xs: IF( ( config_flags%symmetric_xs ) .and.  &
735                          ( its == ids )                  )  THEN
736 
737           IF ( (variable /= 'u') .and. (variable /= 'x') ) THEN
738 
739             DO j = MAX(jds,jts-1), MIN(jte+1,jde+jstag)
740             DO k = kts, k_end
741             DO i = 1, bdyzone
742               dat(ids-i,k,j) = dat(ids+i-1,k,j) !  here, dat(0) = dat(1), etc
743             ENDDO                                 !  symmetry about dat(0.5) (u = 0 pt)
744             ENDDO
745             ENDDO
746 
747           ELSE
748 
749             IF ( variable == 'u' ) THEN
750 
751               DO j = MAX(jds,jts-1), MIN(jte+1,jde+jstag)
752               DO k = kts, k_end
753               DO i = 1, bdyzone
754                 dat(ids-i,k,j) = - dat(ids+i,k,j) ! here, u(0) = - u(2), etc
755               ENDDO                                 !  normal b.c symmetry at u(1)
756               ENDDO
757               ENDDO
758 
759             ELSE
760 
761               DO j = MAX(jds,jts-1), MIN(jte+1,jde+jstag)
762               DO k = kts, k_end
763               DO i = 1, bdyzone
764                 dat(ids-i,k,j) = dat(ids+i,k,j) ! here, phi(0) = phi(2), etc
765               ENDDO                               !  normal b.c symmetry at phi(1)
766               ENDDO
767               ENDDO
768 
769             END IF
770 
771           ENDIF
772 
773         ENDIF symmetry_xs
774 
775 
776 !  now the symmetry boundary at xe
777 
778         symmetry_xe: IF( ( config_flags%symmetric_xe ) .and.  &
779                          ( ite == ide )                  )  THEN
780 
781           IF ( (variable /= 'u') .and. (variable /= 'x') ) THEN
782 
783             DO j = MAX(jds,jts-1), MIN(jte+1,jde+jstag)
784             DO k = kts, k_end
785             DO i = 1, bdyzone
786               dat(ide+i-1,k,j) = dat(ide-i,k,j)  !  sym. about dat(ide-0.5)
787             ENDDO
788             ENDDO
789             ENDDO
790 
791           ELSE
792 
793             IF (variable == 'u') THEN
794 
795               DO j = MAX(jds,jts-1), MIN(jte+1,jde+jstag)
796               DO k = kts, k_end
797               DO i = 1, bdyzone
798                 dat(ide+i,k,j) = - dat(ide-i,k,j)  ! u(ide+1) = - u(ide-1), etc.
799               ENDDO
800               ENDDO
801               ENDDO
802 
803             ELSE
804 
805               DO j = MAX(jds,jts-1), MIN(jte+1,jde+jstag)
806               DO k = kts, k_end
807               DO i = 1, bdyzone
808                 dat(ide+i,k,j) = dat(ide-i,k,j)  ! phi(ide+1) = - phi(ide-1), etc.
809               ENDDO
810               ENDDO
811               ENDDO
812 
813              END IF
814 
815           END IF 
816 
817         END IF symmetry_xe
818 
819 !  set open b.c in X copy into boundary zone here.  WCS, 19 March 2000
820 
821         open_xs: IF( ( config_flags%open_xs   .or. &
822                        config_flags%specified .or. &
823                        config_flags%nested            ) .and.  &
824                          ( its == ids ) .and. open_bc_copy  )  THEN
825 
826             DO j = jts-bdyzone, MIN(jte,jde+jstag)+bdyzone
827             DO k = kts, k_end
828               dat(ids-1,k,j) = dat(ids,k,j) !  here, dat(0) = dat(1), etc
829               dat(ids-2,k,j) = dat(ids,k,j)
830               dat(ids-3,k,j) = dat(ids,k,j)
831             ENDDO
832             ENDDO
833 
834         ENDIF open_xs
835 
836 
837 !  now the open_xe boundary copy 
838 
839         open_xe: IF( ( config_flags%open_xe   .or. &
840                        config_flags%specified .or. &
841                        config_flags%nested            ) .and.  &
842                          ( ite == ide ) .and. open_bc_copy )  THEN
843 
844           IF (variable /= 'u' .and. variable /= 'x' ) THEN
845 
846             DO j = jts-bdyzone, MIN(jte,jde+jstag)+bdyzone
847             DO k = kts, k_end
848               dat(ide  ,k,j) = dat(ide-1,k,j)
849               dat(ide+1,k,j) = dat(ide-1,k,j)
850               dat(ide+2,k,j) = dat(ide-1,k,j)
851             ENDDO
852             ENDDO
853 
854           ELSE
855 
856 !!!!!!! I am not sure about this one!  JM 20020402
857             DO j = MAX(jds,jts-1)-bdyzone, MIN(jte+1,jde+jstag)+bdyzone
858             DO k = kts, k_end
859               dat(ide+1,k,j) = dat(ide,k,j)
860               dat(ide+2,k,j) = dat(ide,k,j)
861               dat(ide+3,k,j) = dat(ide,k,j)
862             ENDDO
863             ENDDO
864 
865           END IF 
866 
867         END IF open_xe
868 
869 !  end open b.c in X copy into boundary zone addition.  WCS, 19 March 2000
870 
871       END IF periodicity_x
872 
873 !  same procedure in y
874 
875       periodicity_y:  IF( ( config_flags%periodic_y ) ) THEN
876         IF ( ( jds == jps ) .and. ( jde == jpe ) )  THEN      ! test if both north and south on processor
877           IF( jts == jds ) then
878 
879             DO j = 0, -(bdyzone-1), -1
880             DO k = kts, k_end
881             DO i = MAX(ids,its-1), MIN(ite+1,ide+istag)
882               dat(i,k,jds+j-1) = dat(i,k,jde+j-1)
883             ENDDO
884             ENDDO
885             ENDDO
886 
887           END IF
888 
889           IF( jte == jde ) then
890 
891             DO j = -jstag, bdyzone
892             DO k = kts, k_end
893             DO i = MAX(ids,its-1), MIN(ite+1,ide+istag)
894               dat(i,k,jde+j+jstag) = dat(i,k,jds+j+jstag)
895             ENDDO
896             ENDDO
897             ENDDO
898 
899           END IF
900 
901         END IF
902 
903       ELSE
904 
905         symmetry_ys: IF( ( config_flags%symmetric_ys ) .and.  &
906                          ( jts == jds)                   )  THEN
907 
908           IF ( (variable /= 'v') .and. (variable /= 'y') ) THEN
909 
910             DO j = 1, bdyzone
911             DO k = kts, k_end
912             DO i = MAX(ids,its-1), MIN(ite+1,ide+istag)
913               dat(i,k,jds-j) = dat(i,k,jds+j-1) 
914             ENDDO                               
915             ENDDO
916             ENDDO
917 
918           ELSE
919 
920             IF (variable == 'v') THEN
921 
922               DO j = 1, bdyzone
923               DO k = kts, k_end
924               DO i = MAX(ids,its-1), MIN(ite+1,ide+istag)
925                 dat(i,k,jds-j) = - dat(i,k,jds+j) 
926               ENDDO              
927               ENDDO
928               ENDDO
929 
930             ELSE
931 
932               DO j = 1, bdyzone
933               DO k = kts, k_end
934               DO i = MAX(ids,its-1), MIN(ite+1,ide+istag)
935                 dat(i,k,jds-j) = dat(i,k,jds+j) 
936               ENDDO              
937               ENDDO
938               ENDDO
939 
940             END IF
941 
942           ENDIF
943 
944         ENDIF symmetry_ys
945 
946 !  now the symmetry boundary at ye
947 
948         symmetry_ye: IF( ( config_flags%symmetric_ye ) .and.  &
949                          ( jte == jde )                  )  THEN
950 
951           IF ( (variable /= 'v') .and. (variable /= 'y') ) THEN
952 
953             DO j = 1, bdyzone
954             DO k = kts, k_end
955             DO i = MAX(ids,its-1), MIN(ite+1,ide+istag)
956               dat(i,k,jde+j-1) = dat(i,k,jde-j) 
957             ENDDO                               
958             ENDDO
959             ENDDO
960 
961           ELSE
962 
963             IF ( variable == 'v' ) THEN
964 
965               DO j = 1, bdyzone
966               DO k = kts, k_end
967               DO i = MAX(ids,its-1), MIN(ite+1,ide+istag)
968                 dat(i,k,jde+j) = - dat(i,k,jde-j) 
969               ENDDO                               
970               ENDDO
971               ENDDO
972 
973             ELSE
974 
975               DO j = 1, bdyzone
976               DO k = kts, k_end
977               DO i = MAX(ids,its-1), MIN(ite+1,ide+istag)
978                 dat(i,k,jde+j) = dat(i,k,jde-j) 
979               ENDDO                               
980               ENDDO
981               ENDDO
982 
983             END IF
984 
985           ENDIF
986 
987         END IF symmetry_ye
988       
989 !  set open b.c in Y copy into boundary zone here.  WCS, 19 March 2000
990 
991         open_ys: IF( ( config_flags%open_ys   .or. &
992                        config_flags%specified .or. &
993                        config_flags%nested            ) .and.  &
994                          ( jts == jds) .and. open_bc_copy )  THEN
995 
996             DO k = kts, k_end
997             DO i = MAX(ids,its-1), MIN(ite+1,ide+istag)
998               dat(i,k,jds-1) = dat(i,k,jds) 
999               dat(i,k,jds-2) = dat(i,k,jds) 
1000               dat(i,k,jds-3) = dat(i,k,jds) 
1001             ENDDO
1002             ENDDO
1003 
1004         ENDIF open_ys
1005 
1006 !  now the open boundary copy at ye
1007 
1008         open_ye: IF( ( config_flags%open_ye   .or. &
1009                        config_flags%specified .or. &
1010                        config_flags%nested            ) .and.  &
1011                          ( jte == jde ) .and. open_bc_copy )  THEN
1012 
1013           IF (variable /= 'v' .and. variable /= 'y' ) THEN
1014 
1015             DO k = kts, k_end
1016             DO i = MAX(ids,its-1), MIN(ite+1,ide+istag)
1017               dat(i,k,jde  ) = dat(i,k,jde-1) 
1018               dat(i,k,jde+1) = dat(i,k,jde-1) 
1019               dat(i,k,jde+2) = dat(i,k,jde-1) 
1020             ENDDO                               
1021             ENDDO
1022 
1023           ELSE
1024 
1025             DO k = kts, k_end
1026             DO i = MAX(ids,its-1), MIN(ite+1,ide+istag)
1027               dat(i,k,jde+1) = dat(i,k,jde) 
1028               dat(i,k,jde+2) = dat(i,k,jde) 
1029               dat(i,k,jde+3) = dat(i,k,jde) 
1030             ENDDO                               
1031             ENDDO
1032 
1033           ENDIF
1034 
1035       END IF open_ye
1036 
1037 !  end open b.c in Y copy into boundary zone addition.  WCS, 19 March 2000
1038 
1039       END IF periodicity_y
1040 
1041 !  fix corners for doubly periodic domains
1042 
1043       IF ( config_flags%periodic_x .and. config_flags%periodic_y &
1044            .and. (ids == ips) .and. (ide == ipe)                 &
1045            .and. (jds == jps) .and. (jde == jpe)                   ) THEN
1046 
1047          IF ( (its == ids) .and. (jts == jds) ) THEN  ! lower left corner fill
1048            DO j = 0, -(bdyzone-1), -1
1049            DO k = kts, k_end
1050            DO i = 0, -(bdyzone-1), -1
1051              dat(ids+i-1,k,jds+j-1) = dat(ide+i-1,k,jde+j-1)
1052            ENDDO
1053            ENDDO
1054            ENDDO
1055          END IF
1056 
1057          IF ( (ite == ide) .and. (jts == jds) ) THEN  ! lower right corner fill
1058            DO j = 0, -(bdyzone-1), -1
1059            DO k = kts, k_end
1060            DO i = 1, bdyzone
1061              dat(ide+i+istag,k,jds+j-1) = dat(ids+i+istag,k,jde+j-1)
1062            ENDDO
1063            ENDDO
1064            ENDDO
1065          END IF
1066 
1067          IF ( (ite == ide) .and. (jte == jde) ) THEN  ! upper right corner fill
1068            DO j = 1, bdyzone
1069            DO k = kts, k_end
1070            DO i = 1, bdyzone
1071              dat(ide+i+istag,k,jde+j+jstag) = dat(ids+i+istag,k,jds+j+jstag)
1072            ENDDO
1073            ENDDO
1074            ENDDO
1075          END IF
1076 
1077          IF ( (its == ids) .and. (jte == jde) ) THEN  ! upper left corner fill
1078            DO j = 1, bdyzone
1079            DO k = kts, k_end
1080            DO i = 0, -(bdyzone-1), -1
1081              dat(ids+i-1,k,jde+j+jstag) = dat(ide+i-1,k,jds+j+jstag)
1082            ENDDO
1083            ENDDO
1084            ENDDO
1085          END IF
1086 
1087        END IF
1088 
1089    END SUBROUTINE set_physical_bc3d
1090 
1091    SUBROUTINE init_module_bc
1092    END SUBROUTINE init_module_bc
1093 
1094 !------------------------------------------------------------------------
1095    SUBROUTINE relax_bdytend   ( field, field_tend,                     &
1096                                 field_bdy_xs, field_bdy_xe,            &
1097                                 field_bdy_ys, field_bdy_ye,            &
1098                                 field_bdy_tend_xs, field_bdy_tend_xe,  &
1099                                 field_bdy_tend_ys, field_bdy_tend_ye,  &
1100                                 variable_in, config_flags,             &
1101                                 spec_bdy_width, spec_zone, relax_zone, &
1102                                 dtbc, fcx, gcx,             &
1103                                 ids,ide, jds,jde, kds,kde,  & ! domain dims
1104                                 ims,ime, jms,jme, kms,kme,  & ! memory dims
1105                                 ips,ipe, jps,jpe, kps,kpe,  & ! patch  dims
1106                                 its,ite, jts,jte, kts,kte )
1107 
1108 !  This subroutine adds the tendencies in the boundary relaxation region, for specified
1109 !  boundary conditions.  
1110 !  spec_bdy_width is only used to dimension the boundary arrays.
1111 !  relax_zone is the inner edge of the boundary relaxation zone treated here.
1112 !  spec_zone is the width of the outer specified b.c.s that are not changed here.
1113 !  (JD July 2000)
1114 
1115       IMPLICIT NONE
1116 
1117       INTEGER,      INTENT(IN   )    :: ids,ide, jds,jde, kds,kde
1118       INTEGER,      INTENT(IN   )    :: ims,ime, jms,jme, kms,kme
1119       INTEGER,      INTENT(IN   )    :: ips,ipe, jps,jpe, kps,kpe
1120       INTEGER,      INTENT(IN   )    :: its,ite, jts,jte, kts,kte
1121       INTEGER,      INTENT(IN   )    :: spec_bdy_width, spec_zone, relax_zone
1122       REAL,         INTENT(IN   )    :: dtbc
1123       CHARACTER,    INTENT(IN   )    :: variable_in
1124 
1125 
1126       REAL,  DIMENSION( ims:ime , kms:kme , jms:jme ), INTENT(IN   ) :: field
1127       REAL,  DIMENSION( ims:ime , kms:kme , jms:jme ), INTENT(INOUT) :: field_tend
1128       REAL,  DIMENSION( jms:jme , kds:kde , spec_bdy_width ), INTENT(IN   ) :: field_bdy_xs, field_bdy_xe
1129       REAL,  DIMENSION( ims:ime , kds:kde , spec_bdy_width ), INTENT(IN   ) :: field_bdy_ys, field_bdy_ye
1130       REAL,  DIMENSION( jms:jme , kds:kde , spec_bdy_width ), INTENT(IN   ) :: field_bdy_tend_xs, field_bdy_tend_xe
1131       REAL,  DIMENSION( ims:ime , kds:kde , spec_bdy_width ), INTENT(IN   ) :: field_bdy_tend_ys, field_bdy_tend_ye
1132       REAL,  DIMENSION( spec_bdy_width ), INTENT(IN   ) :: fcx, gcx
1133       TYPE( grid_config_rec_type ) config_flags
1134 
1135       CHARACTER  :: variable
1136       INTEGER    :: i, j, k, ibs, ibe, jbs, jbe, itf, jtf, ktf, im1, ip1
1137       INTEGER    :: b_dist, b_limit
1138       REAL       :: fls0, fls1, fls2, fls3, fls4
1139       LOGICAL    :: periodic_x
1140 
1141       periodic_x = config_flags%periodic_x
1142       variable = variable_in
1143 
1144       IF (variable == 'U') variable = 'u'
1145       IF (variable == 'V') variable = 'v'
1146       IF (variable == 'M') variable = 'm'
1147       IF (variable == 'H') variable = 'h'
1148 
1149       ibs = ids
1150       ibe = ide-1
1151       itf = min(ite,ide-1)
1152       jbs = jds
1153       jbe = jde-1
1154       jtf = min(jte,jde-1)
1155       ktf = kde-1
1156       IF (variable == 'u') ibe = ide
1157       IF (variable == 'u') itf = min(ite,ide)
1158       IF (variable == 'v') jbe = jde
1159       IF (variable == 'v') jtf = min(jte,jde)
1160       IF (variable == 'm') ktf = kte
1161       IF (variable == 'h') ktf = kte
1162 
1163       IF (jts - jbs .lt. relax_zone) THEN
1164 ! Y-start boundary
1165         DO j = max(jts,jbs+spec_zone), min(jtf,jbs+relax_zone-1)
1166           b_dist = j - jbs
1167           b_limit = b_dist
1168           IF(periodic_x)b_limit = 0
1169           DO k = kts, ktf
1170             DO i = max(its,b_limit+ibs), min(itf,ibe-b_limit)
1171               im1 = max(i-1,ibs)
1172               ip1 = min(i+1,ibe)
1173               fls0 = field_bdy_ys(i, k, b_dist+1) &
1174                    + dtbc * field_bdy_tend_ys(i, k, b_dist+1) &
1175                    - field(i,k,j)
1176               fls1 = field_bdy_ys(im1, k, b_dist+1) &
1177                    + dtbc * field_bdy_tend_ys(im1, k, b_dist+1) &
1178                    - field(im1,k,j)
1179               fls2 = field_bdy_ys(ip1, k, b_dist+1) &
1180                    + dtbc * field_bdy_tend_ys(ip1, k, b_dist+1) &
1181                    - field(ip1,k,j)
1182               fls3 = field_bdy_ys(i, k, b_dist) &
1183                    + dtbc * field_bdy_tend_ys(i, k, b_dist) &
1184                    - field(i,k,j-1)
1185               fls4 = field_bdy_ys(i, k, b_dist+2) &
1186                    + dtbc * field_bdy_tend_ys(i, k, b_dist+2) &
1187                    - field(i,k,j+1)
1188               field_tend(i,k,j) = field_tend(i,k,j) &
1189                                 + fcx(b_dist+1)*fls0 &
1190                                 - gcx(b_dist+1)*(fls1+fls2+fls3+fls4-4.*fls0)
1191             ENDDO
1192           ENDDO
1193         ENDDO
1194       ENDIF
1195 
1196       IF (jbe - jtf .lt. relax_zone) THEN
1197 ! Y-end boundary
1198         DO j = max(jts,jbe-relax_zone+1), min(jtf,jbe-spec_zone)
1199           b_dist = jbe - j
1200           b_limit = b_dist
1201           IF(periodic_x)b_limit = 0
1202           DO k = kts, ktf
1203             DO i = max(its,b_limit+ibs), min(itf,ibe-b_limit)
1204               im1 = max(i-1,ibs)
1205               ip1 = min(i+1,ibe)
1206               fls0 = field_bdy_ye(i, k, b_dist+1) &
1207                    + dtbc * field_bdy_tend_ye(i, k, b_dist+1) &
1208                    - field(i,k,j)
1209               fls1 = field_bdy_ye(im1, k, b_dist+1) &
1210                    + dtbc * field_bdy_tend_ye(im1, k, b_dist+1) &
1211                    - field(im1,k,j)
1212               fls2 = field_bdy_ye(ip1, k, b_dist+1) &
1213                    + dtbc * field_bdy_tend_ye(ip1, k, b_dist+1) &
1214                    - field(ip1,k,j)
1215               fls3 = field_bdy_ye(i, k, b_dist) &
1216                    + dtbc * field_bdy_tend_ye(i, k, b_dist) &
1217                    - field(i,k,j+1)
1218               fls4 = field_bdy_ye(i, k, b_dist+2) &
1219                    + dtbc * field_bdy_tend_ye(i, k, b_dist+2) &
1220                    - field(i,k,j-1)
1221               field_tend(i,k,j) = field_tend(i,k,j) &
1222                                 + fcx(b_dist+1)*fls0 &
1223                                 - gcx(b_dist+1)*(fls1+fls2+fls3+fls4-4.*fls0)
1224 
1225             ENDDO
1226           ENDDO
1227         ENDDO
1228       ENDIF
1229 
1230     IF(.NOT.periodic_x)THEN
1231       IF (its - ibs .lt. relax_zone) THEN
1232 ! X-start boundary
1233         DO i = max(its,ibs+spec_zone), min(itf,ibs+relax_zone-1)
1234           b_dist = i - ibs
1235           DO k = kts, ktf
1236             DO j = max(jts,b_dist+jbs+1), min(jtf,jbe-b_dist-1)
1237               fls0 = field_bdy_xs(j, k, b_dist+1) &
1238                    + dtbc * field_bdy_tend_xs(j, k, b_dist+1) &
1239                    - field(i,k,j)
1240               fls1 = field_bdy_xs(j-1, k, b_dist+1) &
1241                    + dtbc * field_bdy_tend_xs(j-1, k, b_dist+1) &
1242                    - field(i,k,j-1)
1243               fls2 = field_bdy_xs(j+1, k, b_dist+1) &
1244                    + dtbc * field_bdy_tend_xs(j+1, k, b_dist+1) &
1245                    - field(i,k,j+1)
1246               fls3 = field_bdy_xs(j, k, b_dist) &
1247                    + dtbc * field_bdy_tend_xs(j, k, b_dist) &
1248                    - field(i-1,k,j)
1249               fls4 = field_bdy_xs(j, k, b_dist+2) &
1250                    + dtbc * field_bdy_tend_xs(j, k, b_dist+2) &
1251                    - field(i+1,k,j)
1252               field_tend(i,k,j) = field_tend(i,k,j) &
1253                                 + fcx(b_dist+1)*fls0 &
1254                                 - gcx(b_dist+1)*(fls1+fls2+fls3+fls4-4.*fls0)
1255 
1256             ENDDO
1257           ENDDO
1258         ENDDO
1259       ENDIF
1260 
1261       IF (ibe - itf .lt. relax_zone) THEN
1262 ! X-end boundary
1263         DO i = max(its,ibe-relax_zone+1), min(itf,ibe-spec_zone)
1264           b_dist = ibe - i
1265           DO k = kts, ktf
1266             DO j = max(jts,b_dist+jbs+1), min(jtf,jbe-b_dist-1)
1267               fls0 = field_bdy_xe(j, k, b_dist+1) &
1268                    + dtbc * field_bdy_tend_xe(j, k, b_dist+1) &
1269                    - field(i,k,j)
1270               fls1 = field_bdy_xe(j-1, k, b_dist+1) &
1271                    + dtbc * field_bdy_tend_xe(j-1, k, b_dist+1) &
1272                    - field(i,k,j-1)
1273               fls2 = field_bdy_xe(j+1, k, b_dist+1) &
1274                    + dtbc * field_bdy_tend_xe(j+1, k, b_dist+1) &
1275                    - field(i,k,j+1)
1276               fls3 = field_bdy_xe(j, k, b_dist) &
1277                    + dtbc * field_bdy_tend_xe(j, k, b_dist) &
1278                    - field(i+1,k,j)
1279               fls4 = field_bdy_xe(j, k, b_dist+2) &
1280                    + dtbc * field_bdy_tend_xe(j, k, b_dist+2) &
1281                    - field(i-1,k,j)
1282               field_tend(i,k,j) = field_tend(i,k,j) &
1283                                 + fcx(b_dist+1)*fls0 &
1284                                 - gcx(b_dist+1)*(fls1+fls2+fls3+fls4-4.*fls0)
1285             ENDDO
1286           ENDDO
1287         ENDDO
1288       ENDIF
1289     ENDIF
1290 
1291 
1292    END SUBROUTINE relax_bdytend
1293 !------------------------------------------------------------------------
1294 
1295    SUBROUTINE spec_bdytend   ( field_tend,                           &
1296                                field_bdy_xs, field_bdy_xe,           &
1297                                field_bdy_ys, field_bdy_ye,           &
1298                                field_bdy_tend_xs, field_bdy_tend_xe, &
1299                                field_bdy_tend_ys, field_bdy_tend_ye, &
1300                                variable_in, config_flags, & 
1301                                spec_bdy_width, spec_zone, &
1302                                ids,ide, jds,jde, kds,kde,  & ! domain dims
1303                                ims,ime, jms,jme, kms,kme,  & ! memory dims
1304                                ips,ipe, jps,jpe, kps,kpe,  & ! patch  dims
1305                                its,ite, jts,jte, kts,kte )
1306 
1307 !  This subroutine sets the tendencies in the boundary specified region.
1308 !  spec_bdy_width is only used to dimension the boundary arrays.
1309 !  spec_zone is the width of the outer specified b.c.s that are set here.
1310 !  (JD July 2000)
1311 
1312       IMPLICIT NONE
1313 
1314       INTEGER,      INTENT(IN   )    :: ids,ide, jds,jde, kds,kde
1315       INTEGER,      INTENT(IN   )    :: ims,ime, jms,jme, kms,kme
1316       INTEGER,      INTENT(IN   )    :: ips,ipe, jps,jpe, kps,kpe
1317       INTEGER,      INTENT(IN   )    :: its,ite, jts,jte, kts,kte
1318       INTEGER,      INTENT(IN   )    :: spec_bdy_width, spec_zone
1319       CHARACTER,    INTENT(IN   )    :: variable_in
1320 
1321 
1322       REAL,  DIMENSION( ims:ime , kms:kme , jms:jme ), INTENT(OUT  ) :: field_tend
1323       REAL,  DIMENSION( jms:jme , kds:kde , spec_bdy_width ), INTENT(IN   ) :: field_bdy_xs, field_bdy_xe
1324       REAL,  DIMENSION( ims:ime , kds:kde , spec_bdy_width ), INTENT(IN   ) :: field_bdy_ys, field_bdy_ye 
1325       REAL,  DIMENSION( jms:jme , kds:kde , spec_bdy_width ), INTENT(IN   ) :: field_bdy_tend_xs, field_bdy_tend_xe
1326       REAL,  DIMENSION( ims:ime , kds:kde , spec_bdy_width ), INTENT(IN   ) :: field_bdy_tend_ys, field_bdy_tend_ye 
1327       TYPE( grid_config_rec_type ) config_flags
1328 
1329       CHARACTER  :: variable
1330       INTEGER    :: i, j, k, ibs, ibe, jbs, jbe, itf, jtf, ktf
1331       INTEGER    :: b_dist, b_limit
1332       LOGICAL    :: periodic_x
1333 
1334       periodic_x = config_flags%periodic_x
1335 
1336       variable = variable_in
1337 
1338       IF (variable == 'U') variable = 'u'
1339       IF (variable == 'V') variable = 'v'
1340       IF (variable == 'M') variable = 'm'
1341       IF (variable == 'H') variable = 'h'
1342 
1343       ibs = ids
1344       ibe = ide-1
1345       itf = min(ite,ide-1)
1346       jbs = jds
1347       jbe = jde-1
1348       jtf = min(jte,jde-1)
1349       ktf = kde-1
1350       IF (variable == 'u') ibe = ide
1351       IF (variable == 'u') itf = min(ite,ide)
1352       IF (variable == 'v') jbe = jde
1353       IF (variable == 'v') jtf = min(jte,jde)
1354       IF (variable == 'm') ktf = kte
1355       IF (variable == 'h') ktf = kte
1356 
1357       IF (jts - jbs .lt. spec_zone) THEN
1358 ! Y-start boundary
1359         DO j = jts, min(jtf,jbs+spec_zone-1)
1360           b_dist = j - jbs
1361           b_limit = b_dist
1362           IF(periodic_x)b_limit = 0
1363           DO k = kts, ktf
1364             DO i = max(its,b_limit+ibs), min(itf,ibe-b_limit)
1365               field_tend(i,k,j) = field_bdy_tend_ys(i, k, b_dist+1)
1366             ENDDO
1367           ENDDO
1368         ENDDO
1369       ENDIF 
1370       IF (jbe - jtf .lt. spec_zone) THEN 
1371 ! Y-end boundary 
1372         DO j = max(jts,jbe-spec_zone+1), jtf 
1373           b_dist = jbe - j 
1374           b_limit = b_dist
1375           IF(periodic_x)b_limit = 0
1376           DO k = kts, ktf 
1377             DO i = max(its,b_limit+ibs), min(itf,ibe-b_limit)
1378               field_tend(i,k,j) = field_bdy_tend_ye(i, k, b_dist+1)
1379             ENDDO
1380           ENDDO
1381         ENDDO
1382       ENDIF 
1383 
1384     IF(.NOT.periodic_x)THEN
1385       IF (its - ibs .lt. spec_zone) THEN
1386 ! X-start boundary
1387         DO i = its, min(itf,ibs+spec_zone-1)
1388           b_dist = i - ibs
1389           DO k = kts, ktf
1390             DO j = max(jts,b_dist+jbs+1), min(jtf,jbe-b_dist-1)
1391               field_tend(i,k,j) = field_bdy_tend_xs(j, k, b_dist+1)
1392             ENDDO
1393           ENDDO
1394         ENDDO
1395       ENDIF 
1396 
1397       IF (ibe - itf .lt. spec_zone) THEN
1398 ! X-end boundary
1399         DO i = max(its,ibe-spec_zone+1), itf
1400           b_dist = ibe - i
1401           DO k = kts, ktf
1402             DO j = max(jts,b_dist+jbs+1), min(jtf,jbe-b_dist-1)
1403               field_tend(i,k,j) = field_bdy_tend_xe(j, k, b_dist+1)
1404             ENDDO
1405           ENDDO
1406         ENDDO
1407       ENDIF 
1408     ENDIF
1409 
1410    END SUBROUTINE spec_bdytend
1411 !------------------------------------------------------------------------
1412 
1413    SUBROUTINE spec_bdyupdate(  field,      &
1414                                field_tend, dt,            &
1415                                variable_in, config_flags, & 
1416                                spec_zone,                  &
1417                                ids,ide, jds,jde, kds,kde,  & ! domain dims
1418                                ims,ime, jms,jme, kms,kme,  & ! memory dims
1419                                ips,ipe, jps,jpe, kps,kpe,  & ! patch  dims
1420                                its,ite, jts,jte, kts,kte )
1421 
1422 !  This subroutine adds the tendencies in the boundary specified region.
1423 !  spec_zone is the width of the outer specified b.c.s that are set here.
1424 !  (JD August 2000)
1425 
1426       IMPLICIT NONE
1427 
1428       INTEGER,      INTENT(IN   )    :: ids,ide, jds,jde, kds,kde
1429       INTEGER,      INTENT(IN   )    :: ims,ime, jms,jme, kms,kme
1430       INTEGER,      INTENT(IN   )    :: ips,ipe, jps,jpe, kps,kpe
1431       INTEGER,      INTENT(IN   )    :: its,ite, jts,jte, kts,kte
1432       INTEGER,      INTENT(IN   )    :: spec_zone
1433       CHARACTER,    INTENT(IN   )    :: variable_in
1434       REAL,         INTENT(IN   )    :: dt
1435 
1436 
1437       REAL,  DIMENSION( ims:ime , kms:kme , jms:jme ), INTENT(INOUT) :: field
1438       REAL,  DIMENSION( ims:ime , kms:kme , jms:jme ), INTENT(IN   ) :: field_tend
1439       TYPE( grid_config_rec_type ) config_flags
1440 
1441       CHARACTER  :: variable
1442       INTEGER    :: i, j, k, ibs, ibe, jbs, jbe, itf, jtf, ktf
1443       INTEGER    :: b_dist, b_limit
1444       LOGICAL    :: periodic_x
1445 
1446       periodic_x = config_flags%periodic_x
1447 
1448       variable = variable_in
1449 
1450       IF (variable == 'U') variable = 'u'
1451       IF (variable == 'V') variable = 'v'
1452       IF (variable == 'M') variable = 'm'
1453       IF (variable == 'H') variable = 'h'
1454 
1455       ibs = ids
1456       ibe = ide-1
1457       itf = min(ite,ide-1)
1458       jbs = jds
1459       jbe = jde-1
1460       jtf = min(jte,jde-1)
1461       ktf = kde-1
1462       IF (variable == 'u') ibe = ide
1463       IF (variable == 'u') itf = min(ite,ide)
1464       IF (variable == 'v') jbe = jde
1465       IF (variable == 'v') jtf = min(jte,jde)
1466       IF (variable == 'm') ktf = kte
1467       IF (variable == 'h') ktf = kte
1468 
1469       IF (jts - jbs .lt. spec_zone) THEN
1470 ! Y-start boundary
1471         DO j = jts, min(jtf,jbs+spec_zone-1)
1472           b_dist = j - jbs
1473           b_limit = b_dist
1474           IF(periodic_x)b_limit = 0
1475           DO k = kts, ktf
1476             DO i = max(its,b_limit+ibs), min(itf,ibe-b_limit)
1477               field(i,k,j) = field(i,k,j) + dt*field_tend(i,k,j) 
1478             ENDDO
1479           ENDDO
1480         ENDDO
1481       ENDIF 
1482       IF (jbe - jtf .lt. spec_zone) THEN 
1483 ! Y-end boundary 
1484         DO j = max(jts,jbe-spec_zone+1), jtf 
1485           b_dist = jbe - j 
1486           b_limit = b_dist
1487           IF(periodic_x)b_limit = 0
1488           DO k = kts, ktf 
1489             DO i = max(its,b_limit+ibs), min(itf,ibe-b_limit)
1490               field(i,k,j) = field(i,k,j) + dt*field_tend(i,k,j) 
1491             ENDDO
1492           ENDDO
1493         ENDDO
1494       ENDIF 
1495 
1496     IF(.NOT.periodic_x)THEN
1497       IF (its - ibs .lt. spec_zone) THEN
1498 ! X-start boundary
1499         DO i = its, min(itf,ibs+spec_zone-1)
1500           b_dist = i - ibs
1501           DO k = kts, ktf
1502             DO j = max(jts,b_dist+jbs+1), min(jtf,jbe-b_dist-1)
1503               field(i,k,j) = field(i,k,j) + dt*field_tend(i,k,j) 
1504             ENDDO
1505           ENDDO
1506         ENDDO
1507       ENDIF 
1508 
1509       IF (ibe - itf .lt. spec_zone) THEN
1510 ! X-end boundary
1511         DO i = max(its,ibe-spec_zone+1), itf
1512           b_dist = ibe - i
1513           DO k = kts, ktf
1514             DO j = max(jts,b_dist+jbs+1), min(jtf,jbe-b_dist-1)
1515               field(i,k,j) = field(i,k,j) + dt*field_tend(i,k,j) 
1516             ENDDO
1517           ENDDO
1518         ENDDO
1519       ENDIF 
1520     ENDIF
1521 
1522    END SUBROUTINE spec_bdyupdate
1523 !------------------------------------------------------------------------
1524 
1525    SUBROUTINE zero_grad_bdy (  field,                     &
1526                                variable_in, config_flags, & 
1527                                spec_zone,                  &
1528                                ids,ide, jds,jde, kds,kde,  & ! domain dims
1529                                ims,ime, jms,jme, kms,kme,  & ! memory dims
1530                                ips,ipe, jps,jpe, kps,kpe,  & ! patch  dims
1531                                its,ite, jts,jte, kts,kte )
1532 
1533 !  This subroutine sets zero gradient conditions in the boundary specified region.
1534 !  spec_zone is the width of the outer specified b.c.s that are set here.
1535 !  (JD August 2000)
1536 
1537       IMPLICIT NONE
1538 
1539       INTEGER,      INTENT(IN   )    :: ids,ide, jds,jde, kds,kde
1540       INTEGER,      INTENT(IN   )    :: ims,ime, jms,jme, kms,kme
1541       INTEGER,      INTENT(IN   )    :: ips,ipe, jps,jpe, kps,kpe
1542       INTEGER,      INTENT(IN   )    :: its,ite, jts,jte, kts,kte
1543       INTEGER,      INTENT(IN   )    :: spec_zone
1544       CHARACTER,    INTENT(IN   )    :: variable_in
1545 
1546 
1547       REAL,  DIMENSION( ims:ime , kms:kme , jms:jme ), INTENT(INOUT) :: field
1548       TYPE( grid_config_rec_type ) config_flags
1549 
1550       CHARACTER  :: variable
1551       INTEGER    :: i, j, k, ibs, ibe, jbs, jbe, itf, jtf, ktf, i_inner, j_inner
1552       INTEGER    :: b_dist, b_limit
1553       LOGICAL    :: periodic_x
1554 
1555       periodic_x = config_flags%periodic_x
1556 
1557       variable = variable_in
1558 
1559       IF (variable == 'U') variable = 'u'
1560       IF (variable == 'V') variable = 'v'
1561 
1562       ibs = ids
1563       ibe = ide-1
1564       itf = min(ite,ide-1)
1565       jbs = jds
1566       jbe = jde-1
1567       jtf = min(jte,jde-1)
1568       ktf = kde-1
1569       IF (variable == 'u') ibe = ide
1570       IF (variable == 'u') itf = min(ite,ide)
1571       IF (variable == 'v') jbe = jde
1572       IF (variable == 'v') jtf = min(jte,jde)
1573       IF (variable == 'w') ktf = kde
1574 
1575       IF (jts - jbs .lt. spec_zone) THEN
1576 ! Y-start boundary
1577         DO j = jts, min(jtf,jbs+spec_zone-1)
1578           b_dist = j - jbs
1579           b_limit = b_dist
1580           IF(periodic_x)b_limit = 0
1581           DO k = kts, ktf
1582             DO i = max(its,b_limit+ibs), min(itf,ibe-b_limit)
1583               i_inner = max(i,ibs+spec_zone)
1584               i_inner = min(i_inner,ibe-spec_zone)
1585               IF(periodic_x)i_inner = i
1586               field(i,k,j) = field(i_inner,k,jbs+spec_zone)
1587             ENDDO
1588           ENDDO
1589         ENDDO
1590       ENDIF 
1591       IF (jbe - jtf .lt. spec_zone) THEN 
1592 ! Y-end boundary 
1593         DO j = max(jts,jbe-spec_zone+1), jtf 
1594           b_dist = jbe - j 
1595           b_limit = b_dist
1596           IF(periodic_x)b_limit = 0
1597           DO k = kts, ktf 
1598             DO i = max(its,b_limit+ibs), min(itf,ibe-b_limit)
1599               i_inner = max(i,ibs+spec_zone)
1600               i_inner = min(i_inner,ibe-spec_zone)
1601               IF(periodic_x)i_inner = i
1602               field(i,k,j) = field(i_inner,k,jbe-spec_zone)
1603             ENDDO
1604           ENDDO
1605         ENDDO
1606       ENDIF 
1607 
1608     IF(.NOT.periodic_x)THEN
1609       IF (its - ibs .lt. spec_zone) THEN
1610 ! X-start boundary
1611         DO i = its, min(itf,ibs+spec_zone-1)
1612           b_dist = i - ibs
1613           DO k = kts, ktf
1614             DO j = max(jts,b_dist+jbs+1), min(jtf,jbe-b_dist-1)
1615               j_inner = max(j,jbs+spec_zone)
1616               j_inner = min(j_inner,jbe-spec_zone)
1617               field(i,k,j) = field(ibs+spec_zone,k,j_inner)
1618             ENDDO
1619           ENDDO
1620         ENDDO
1621       ENDIF 
1622 
1623       IF (ibe - itf .lt. spec_zone) THEN
1624 ! X-end boundary
1625         DO i = max(its,ibe-spec_zone+1), itf
1626           b_dist = ibe - i
1627           DO k = kts, ktf
1628             DO j = max(jts,b_dist+jbs+1), min(jtf,jbe-b_dist-1)
1629               j_inner = max(j,jbs+spec_zone)
1630               j_inner = min(j_inner,jbe-spec_zone)
1631               field(i,k,j) = field(ibe-spec_zone,k,j_inner)
1632             ENDDO
1633           ENDDO
1634         ENDDO
1635       ENDIF 
1636     ENDIF
1637 
1638    END SUBROUTINE zero_grad_bdy
1639 !------------------------------------------------------------------------
1640 
1641    SUBROUTINE flow_dep_bdy  (  field,                     &
1642                                u, v, config_flags, & 
1643                                spec_zone,                  &
1644                                ids,ide, jds,jde, kds,kde,  & ! domain dims
1645                                ims,ime, jms,jme, kms,kme,  & ! memory dims
1646                                ips,ipe, jps,jpe, kps,kpe,  & ! patch  dims
1647                                its,ite, jts,jte, kts,kte )
1648 
1649 !  This subroutine sets zero gradient conditions for outflow and zero value
1650 !  for inflow in the boundary specified region. Note that field must be unstaggered.
1651 !  The velocities, u and v, will only be used to check their sign (coupled vels OK)
1652 !  spec_zone is the width of the outer specified b.c.s that are set here.
1653 !  (JD August 2000)
1654 
1655       IMPLICIT NONE
1656 
1657       INTEGER,      INTENT(IN   )    :: ids,ide, jds,jde, kds,kde
1658       INTEGER,      INTENT(IN   )    :: ims,ime, jms,jme, kms,kme
1659       INTEGER,      INTENT(IN   )    :: ips,ipe, jps,jpe, kps,kpe
1660       INTEGER,      INTENT(IN   )    :: its,ite, jts,jte, kts,kte
1661       INTEGER,      INTENT(IN   )    :: spec_zone
1662 
1663 
1664       REAL,  DIMENSION( ims:ime , kms:kme , jms:jme ), INTENT(INOUT) :: field
1665       REAL,  DIMENSION( ims:ime , kms:kme , jms:jme ), INTENT(IN   ) :: u
1666       REAL,  DIMENSION( ims:ime , kms:kme , jms:jme ), INTENT(IN   ) :: v
1667       TYPE( grid_config_rec_type ) config_flags
1668 
1669       INTEGER    :: i, j, k, ibs, ibe, jbs, jbe, itf, jtf, ktf, i_inner, j_inner
1670       INTEGER    :: b_dist, b_limit
1671       LOGICAL    :: periodic_x
1672 
1673       periodic_x = config_flags%periodic_x
1674 
1675       ibs = ids
1676       ibe = ide-1
1677       itf = min(ite,ide-1)
1678       jbs = jds
1679       jbe = jde-1
1680       jtf = min(jte,jde-1)
1681       ktf = kde-1
1682 
1683       IF (jts - jbs .lt. spec_zone) THEN
1684 ! Y-start boundary
1685         DO j = jts, min(jtf,jbs+spec_zone-1)
1686           b_dist = j - jbs
1687           b_limit = b_dist
1688           IF(periodic_x)b_limit = 0
1689           DO k = kts, ktf
1690             DO i = max(its,b_limit+ibs), min(itf,ibe-b_limit)
1691               i_inner = max(i,ibs+spec_zone)
1692               i_inner = min(i_inner,ibe-spec_zone)
1693               IF(periodic_x)i_inner = i
1694               IF(v(i,k,j) .lt. 0.)THEN
1695                 field(i,k,j) = field(i_inner,k,jbs+spec_zone)
1696               ELSE
1697                 field(i,k,j) = 0.
1698               ENDIF
1699             ENDDO
1700           ENDDO
1701         ENDDO
1702       ENDIF 
1703       IF (jbe - jtf .lt. spec_zone) THEN 
1704 ! Y-end boundary 
1705         DO j = max(jts,jbe-spec_zone+1), jtf 
1706           b_dist = jbe - j 
1707           b_limit = b_dist
1708           IF(periodic_x)b_limit = 0
1709           DO k = kts, ktf 
1710             DO i = max(its,b_limit+ibs), min(itf,ibe-b_limit)
1711               i_inner = max(i,ibs+spec_zone)
1712               i_inner = min(i_inner,ibe-spec_zone)
1713               IF(periodic_x)i_inner = i
1714               IF(v(i,k,j+1) .gt. 0.)THEN
1715                 field(i,k,j) = field(i_inner,k,jbe-spec_zone)
1716               ELSE
1717                 field(i,k,j) = 0.
1718               ENDIF
1719             ENDDO
1720           ENDDO
1721         ENDDO
1722       ENDIF 
1723 
1724     IF(.NOT.periodic_x)THEN
1725       IF (its - ibs .lt. spec_zone) THEN
1726 ! X-start boundary
1727         DO i = its, min(itf,ibs+spec_zone-1)
1728           b_dist = i - ibs
1729           DO k = kts, ktf
1730             DO j = max(jts,b_dist+jbs+1), min(jtf,jbe-b_dist-1)
1731               j_inner = max(j,jbs+spec_zone)
1732               j_inner = min(j_inner,jbe-spec_zone)
1733               IF(u(i,k,j) .lt. 0.)THEN
1734                 field(i,k,j) = field(ibs+spec_zone,k,j_inner)
1735               ELSE
1736                 field(i,k,j) = 0.
1737               ENDIF
1738             ENDDO
1739           ENDDO
1740         ENDDO
1741       ENDIF 
1742 
1743       IF (ibe - itf .lt. spec_zone) THEN
1744 ! X-end boundary
1745         DO i = max(its,ibe-spec_zone+1), itf
1746           b_dist = ibe - i
1747           DO k = kts, ktf
1748             DO j = max(jts,b_dist+jbs+1), min(jtf,jbe-b_dist-1)
1749               j_inner = max(j,jbs+spec_zone)
1750               j_inner = min(j_inner,jbe-spec_zone)
1751               IF(u(i+1,k,j) .gt. 0.)THEN
1752                 field(i,k,j) = field(ibe-spec_zone,k,j_inner)
1753               ELSE
1754                 field(i,k,j) = 0.
1755               ENDIF
1756             ENDDO
1757           ENDDO
1758         ENDDO
1759       ENDIF 
1760     ENDIF
1761 
1762    END SUBROUTINE flow_dep_bdy
1763 
1764 !------------------------------------------------------------------------------
1765 
1766  SUBROUTINE stuff_bdy_new ( data3d , space_bdy_xs, space_bdy_xe, space_bdy_ys, space_bdy_ye, &
1767                              char_stagger , &
1768                              spec_bdy_width , &
1769                              ids, ide, jds, jde, kds, kde , &
1770                              ims, ime, jms, jme, kms, kme , & 
1771                              its, ite, jts, jte, kts, kte )
1772  
1773  !  This routine puts the data in the 3d arrays into the proper locations
1774  !  for the lateral boundary arrays.
1775  
1776     USE module_state_description
1777     
1778     IMPLICIT NONE
1779  
1780     INTEGER , INTENT(IN) :: ids, ide, jds, jde, kds, kde
1781     INTEGER , INTENT(IN) :: ims, ime, jms, jme, kms, kme
1782     INTEGER , INTENT(IN) :: its, ite, jts, jte, kts, kte
1783     INTEGER , INTENT(IN) :: spec_bdy_width
1784     REAL , DIMENSION(ims:ime,kms:kme,jms:jme) , INTENT(IN) :: data3d
1785     REAL , DIMENSION(jms:jme,kds:kde,spec_bdy_width) , INTENT(OUT) :: space_bdy_xs, space_bdy_xe
1786     REAL , DIMENSION(ims:ime,kds:kde,spec_bdy_width) , INTENT(OUT) :: space_bdy_ys, space_bdy_ye
1787     CHARACTER (LEN=1) , INTENT(IN) :: char_stagger
1788  
1789     INTEGER :: i , ii , j , jj , k
1790  
1791     !  There are four lateral boundary locations that are stored.
1792  
1793     !  X start boundary
1794  
1795     IF ( char_stagger .EQ. 'W' ) THEN
1796        DO j = MAX(jds,jts) , MIN(jde-1,jte)
1797        DO k = kds , kde
1798        DO i = MAX(ids,its) , MIN(ids + spec_bdy_width - 1,ite)
1799           space_bdy_xs(j,k,i) = data3d(i,k,j)
1800        END DO
1801        END DO
1802        END DO
1803     ELSE IF ( char_stagger .EQ. 'M' ) THEN
1804        DO j = MAX(jds,jts) , MIN(jde-1,jte)
1805        DO k = kds , kde
1806        DO i = MAX(ids,its) , MIN(ids + spec_bdy_width - 1,ite)
1807           space_bdy_xs(j,k,i) = data3d(i,k,j)
1808        END DO
1809        END DO
1810        END DO
1811     ELSE IF ( char_stagger .EQ. 'V' ) THEN
1812        DO j = MAX(jds,jts) , MIN(jde,jte)
1813        DO k = kds , kde - 1
1814        DO i = MAX(ids,its) , MIN(ids + spec_bdy_width - 1,ite)
1815           space_bdy_xs(j,k,i) = data3d(i,k,j)
1816        END DO
1817        END DO
1818        END DO
1819     ELSE
1820        DO j = MAX(jds,jts) , MIN(jde-1,jte)
1821        DO k = kds , kde - 1
1822        DO i = MAX(ids,its) , MIN(ids + spec_bdy_width - 1,ite)
1823           space_bdy_xs(j,k,i) = data3d(i,k,j)
1824        END DO
1825        END DO
1826        END DO
1827     END IF
1828  
1829     !  X end boundary
1830  
1831     IF      ( char_stagger .EQ. 'U' ) THEN
1832        DO j = MAX(jds,jts) , MIN(jde-1,jte)
1833        DO k = kds , kde - 1
1834        DO i = MIN(ide,ite) , MAX(ide - spec_bdy_width + 1,its) , -1
1835           ii = ide - i + 1
1836           space_bdy_xe(j,k,ii) = data3d(i,k,j)
1837        END DO
1838        END DO
1839        END DO
1840     ELSE IF ( char_stagger .EQ. 'V' ) THEN
1841        DO j = MAX(jds,jts) , MIN(jde,jte)
1842        DO k = kds , kde - 1
1843        DO i = MIN(ide - 1,ite) , MAX(ide - spec_bdy_width,its) , -1
1844           ii = ide - i
1845           space_bdy_xe(j,k,ii) = data3d(i,k,j)
1846        END DO
1847        END DO
1848        END DO
1849     ELSE IF ( char_stagger .EQ. 'W' ) THEN
1850        DO j = MAX(jds,jts) , MIN(jde-1,jte)
1851        DO k = kds , kde
1852        DO i = MIN(ide - 1,ite) , MAX(ide - spec_bdy_width,its) , -1
1853           ii = ide - i
1854           space_bdy_xe(j,k,ii) = data3d(i,k,j)
1855        END DO
1856        END DO
1857        END DO
1858     ELSE IF ( char_stagger .EQ. 'M' ) THEN
1859        DO j = MAX(jds,jts) , MIN(jde-1,jte)
1860        DO k = kds , kde
1861        DO i = MIN(ide - 1,ite) , MAX(ide - spec_bdy_width,its) , -1
1862           ii = ide - i
1863           space_bdy_xe(j,k,ii) = data3d(i,k,j)
1864        END DO
1865        END DO
1866        END DO
1867     ELSE
1868        DO j = MAX(jds,jts) , MIN(jde-1,jte)
1869        DO k = kds , kde - 1
1870        DO i = MIN(ide - 1,ite) , MAX(ide - spec_bdy_width,its) , -1
1871           ii = ide - i
1872           space_bdy_xe(j,k,ii) = data3d(i,k,j)
1873        END DO
1874        END DO
1875        END DO
1876     END IF
1877  
1878     !  Y start boundary
1879  
1880     IF ( char_stagger .EQ. 'W' ) THEN
1881        DO j = MAX(jds,jts) , MIN(jds + spec_bdy_width - 1,jte)
1882        DO k = kds , kde
1883        DO i = MAX(ids,its) , MIN(ide-1,ite)
1884           space_bdy_ys(i,k,j) = data3d(i,k,j)
1885        END DO
1886        END DO
1887        END DO
1888     ELSE IF ( char_stagger .EQ. 'M' ) THEN
1889        DO j = MAX(jds,jts) , MIN(jds + spec_bdy_width - 1,jte)
1890        DO k = kds , kde
1891        DO i = MAX(ids,its) , MIN(ide-1,ite)
1892           space_bdy_ys(i,k,j) = data3d(i,k,j)
1893        END DO
1894        END DO
1895        END DO
1896     ELSE IF ( char_stagger .EQ. 'U' ) THEN
1897        DO j = MAX(jds,jts) , MIN(jds + spec_bdy_width - 1,jte)
1898        DO k = kds , kde - 1
1899        DO i = MAX(ids,its) , MIN(ide,ite)
1900           space_bdy_ys(i,k,j) = data3d(i,k,j)
1901        END DO
1902        END DO
1903        END DO
1904     ELSE
1905        DO j = MAX(jds,jts) , MIN(jds + spec_bdy_width - 1,jte)
1906        DO k = kds , kde - 1
1907        DO i = MAX(ids,its) , MIN(ide-1,ite)
1908           space_bdy_ys(i,k,j) = data3d(i,k,j)
1909        END DO
1910        END DO
1911        END DO
1912     END IF
1913  
1914     !  Y end boundary
1915  
1916     IF      ( char_stagger .EQ. 'V' ) THEN
1917        DO j = MIN(jde,jte) , MAX(jde - spec_bdy_width + 1,jts) , -1
1918        DO k = kds , kde - 1
1919        DO i = MAX(ids,its) , MIN(ide-1,ite)
1920           jj = jde - j + 1
1921           space_bdy_ye(i,k,jj) = data3d(i,k,j)
1922        END DO
1923        END DO
1924        END DO
1925     ELSE IF ( char_stagger .EQ. 'U' ) THEN
1926        DO j = MIN(jde-1,jte) , MAX(jde - spec_bdy_width,jts) , -1
1927        DO k = kds , kde - 1
1928        DO i = MAX(ids,its) , MIN(ide,ite)
1929           jj = jde - j
1930           space_bdy_ye(i,k,jj) = data3d(i,k,j)
1931        END DO
1932        END DO
1933        END DO
1934     ELSE IF ( char_stagger .EQ. 'W' ) THEN
1935        DO j = MIN(jde-1,jte) , MAX(jde - spec_bdy_width,jts) , -1
1936        DO k = kds , kde
1937        DO i = MAX(ids,its) , MIN(ide-1,ite)
1938           jj = jde - j
1939           space_bdy_ye(i,k,jj) = data3d(i,k,j)
1940        END DO
1941        END DO
1942        END DO
1943     ELSE IF ( char_stagger .EQ. 'M' ) THEN
1944        DO j = MIN(jde-1,jte) , MAX(jde - spec_bdy_width,jts) , -1
1945        DO k = kds , kde
1946        DO i = MAX(ids,its) , MIN(ide-1,ite)
1947           jj = jde - j
1948           space_bdy_ye(i,k,jj) = data3d(i,k,j)
1949        END DO
1950        END DO
1951        END DO
1952     ELSE
1953        DO j = MIN(jde-1,jte) , MAX(jde - spec_bdy_width,jts) , -1
1954        DO k = kds , kde - 1
1955        DO i = MAX(ids,its) , MIN(ide-1,ite)
1956           jj = jde - j
1957           space_bdy_ye(i,k,jj) = data3d(i,k,j)
1958        END DO
1959        END DO
1960        END DO
1961     END IF
1962     
1963  END SUBROUTINE stuff_bdy_new
1964  
1965  SUBROUTINE stuff_bdytend_new ( data3dnew , data3dold , time_diff , &
1966                              space_bdy_xs, space_bdy_xe, space_bdy_ys, space_bdy_ye, &
1967                              char_stagger , &
1968                              spec_bdy_width , &
1969                              ids, ide, jds, jde, kds, kde , &
1970                              ims, ime, jms, jme, kms, kme , & 
1971                              its, ite, jts, jte, kts, kte )
1972  
1973  !  This routine puts the tendency data into the proper locations
1974  !  for the lateral boundary arrays.
1975  
1976     USE module_state_description
1977     
1978     IMPLICIT NONE
1979  
1980     INTEGER , INTENT(IN) :: ids, ide, jds, jde, kds, kde
1981     INTEGER , INTENT(IN) :: ims, ime, jms, jme, kms, kme
1982     INTEGER , INTENT(IN) :: its, ite, jts, jte, kts, kte
1983     INTEGER , INTENT(IN) :: spec_bdy_width
1984     REAL , DIMENSION(ims:ime,kms:kme,jms:jme) , INTENT(IN) :: data3dnew , data3dold
1985     REAL , DIMENSION(jms:jme,kds:kde,spec_bdy_width) , INTENT(OUT) :: space_bdy_xs, space_bdy_xe
1986     REAL , DIMENSION(ims:ime,kds:kde,spec_bdy_width) , INTENT(OUT) :: space_bdy_ys, space_bdy_ye
1987     CHARACTER (LEN=1) , INTENT(IN) :: char_stagger
1988     REAL , INTENT(IN) :: time_diff ! seconds
1989  
1990     INTEGER :: i , ii , j , jj , k
1991  
1992     !  There are four lateral boundary locations that are stored.
1993  
1994     !  X start boundary
1995  
1996     IF ( char_stagger .EQ. 'W' ) THEN
1997        DO j = MAX(jds,jts) , MIN(jde-1,jte)
1998        DO k = kds , kde
1999        DO i = MAX(ids,its) , MIN(ids + spec_bdy_width - 1,ite)
2000           space_bdy_xs(j,k,i) = ( data3dnew(i,k,j) - data3dold(i,k,j) ) / time_diff
2001        END DO
2002        END DO
2003        END DO
2004     ELSE IF ( char_stagger .EQ. 'M' ) THEN
2005        DO j = MAX(jds,jts) , MIN(jde-1,jte)
2006        DO k = kds , kde
2007        DO i = MAX(ids,its) , MIN(ids + spec_bdy_width - 1,ite)
2008           space_bdy_xs(j,k,i) = ( data3dnew(i,k,j) - data3dold(i,k,j) ) / time_diff
2009        END DO
2010        END DO
2011        END DO
2012     ELSE IF ( char_stagger .EQ. 'V' ) THEN
2013        DO j = MAX(jds,jts) , MIN(jde,jte)
2014        DO k = kds , kde - 1
2015        DO i = MAX(ids,its) , MIN(ids + spec_bdy_width - 1,ite)
2016           space_bdy_xs(j,k,i) = ( data3dnew(i,k,j) - data3dold(i,k,j) ) / time_diff
2017        END DO
2018        END DO
2019        END DO
2020     ELSE
2021        DO j = MAX(jds,jts) , MIN(jde-1,jte)
2022        DO k = kds , kde - 1
2023        DO i = MAX(ids,its) , MIN(ids + spec_bdy_width - 1,ite)
2024           space_bdy_xs(j,k,i) = ( data3dnew(i,k,j) - data3dold(i,k,j) ) / time_diff
2025        END DO
2026        END DO
2027        END DO
2028     END IF
2029  
2030     !  X end boundary
2031  
2032     IF      ( char_stagger .EQ. 'U' ) THEN
2033        DO j = MAX(jds,jts) , MIN(jde-1,jte)
2034        DO k = kds , kde - 1
2035        DO i = MIN(ide,ite) , MAX(ide - spec_bdy_width + 1,its) , -1
2036           ii = ide - i + 1
2037           space_bdy_xe(j,k,ii) = ( data3dnew(i,k,j) - data3dold(i,k,j) ) / time_diff
2038        END DO
2039        END DO
2040        END DO
2041     ELSE IF ( char_stagger .EQ. 'V' ) THEN
2042        DO j = MAX(jds,jts) , MIN(jde,jte)
2043        DO k = kds , kde - 1
2044        DO i = MIN(ide - 1,ite) , MAX(ide - spec_bdy_width,its) , -1
2045           ii = ide - i
2046           space_bdy_xe(j,k,ii) = ( data3dnew(i,k,j) - data3dold(i,k,j) ) / time_diff
2047        END DO
2048        END DO
2049        END DO
2050     ELSE IF ( char_stagger .EQ. 'W' ) THEN
2051        DO j = MAX(jds,jts) , MIN(jde-1,jte)
2052        DO k = kds , kde
2053        DO i = MIN(ide - 1,ite) , MAX(ide - spec_bdy_width,its) , -1
2054           ii = ide - i
2055           space_bdy_xe(j,k,ii) = ( data3dnew(i,k,j) - data3dold(i,k,j) ) / time_diff
2056        END DO
2057        END DO
2058        END DO
2059     ELSE IF ( char_stagger .EQ. 'M' ) THEN
2060        DO j = MAX(jds,jts) , MIN(jde-1,jte)
2061        DO k = kds , kde
2062        DO i = MIN(ide - 1,ite) , MAX(ide - spec_bdy_width,its) , -1
2063           ii = ide - i
2064           space_bdy_xe(j,k,ii) = ( data3dnew(i,k,j) - data3dold(i,k,j) ) / time_diff
2065        END DO
2066        END DO
2067        END DO
2068     ELSE
2069        DO j = MAX(jds,jts) , MIN(jde-1,jte)
2070        DO k = kds , kde - 1
2071        DO i = MIN(ide - 1,ite) , MAX(ide - spec_bdy_width,its) , -1
2072           ii = ide - i
2073           space_bdy_xe(j,k,ii) = ( data3dnew(i,k,j) - data3dold(i,k,j) ) / time_diff
2074        END DO
2075        END DO
2076        END DO
2077     END IF
2078  
2079     !  Y start boundary
2080  
2081     IF ( char_stagger .EQ. 'W' ) THEN
2082        DO j = MAX(jds,jts) , MIN(jds + spec_bdy_width - 1,jte)
2083        DO k = kds , kde
2084        DO i = MAX(ids,its) , MIN(ide-1,ite)
2085           space_bdy_ys(i,k,j) = ( data3dnew(i,k,j) - data3dold(i,k,j) ) / time_diff
2086        END DO
2087        END DO
2088        END DO
2089     ELSE IF ( char_stagger .EQ. 'M' ) THEN
2090        DO j = MAX(jds,jts) , MIN(jds + spec_bdy_width - 1,jte)
2091        DO k = kds , kde
2092        DO i = MAX(ids,its) , MIN(ide-1,ite)
2093           space_bdy_ys(i,k,j) = ( data3dnew(i,k,j) - data3dold(i,k,j) ) / time_diff
2094        END DO
2095        END DO
2096        END DO
2097     ELSE IF ( char_stagger .EQ. 'U' ) THEN
2098        DO j = MAX(jds,jts) , MIN(jds + spec_bdy_width - 1,jte)
2099        DO k = kds , kde - 1
2100        DO i = MAX(ids,its) , MIN(ide,ite)
2101           space_bdy_ys(i,k,j) = ( data3dnew(i,k,j) - data3dold(i,k,j) ) / time_diff
2102        END DO
2103        END DO
2104        END DO
2105     ELSE
2106        DO j = MAX(jds,jts) , MIN(jds + spec_bdy_width - 1,jte)
2107        DO k = kds , kde - 1
2108        DO i = MAX(ids,its) , MIN(ide-1,ite)
2109           space_bdy_ys(i,k,j) = ( data3dnew(i,k,j) - data3dold(i,k,j) ) / time_diff
2110        END DO
2111        END DO
2112        END DO
2113     END IF
2114  
2115     !  Y end boundary
2116  
2117     IF      ( char_stagger .EQ. 'V' ) THEN
2118        DO j = MIN(jde,jte) , MAX(jde - spec_bdy_width + 1,jts) , -1
2119        DO k = kds , kde - 1
2120        DO i = MAX(ids,its) , MIN(ide-1,ite)
2121           jj = jde - j + 1
2122           space_bdy_ye(i,k,jj) = ( data3dnew(i,k,j) - data3dold(i,k,j) ) / time_diff
2123        END DO
2124        END DO
2125        END DO
2126     ELSE IF ( char_stagger .EQ. 'U' ) THEN
2127        DO j = MIN(jde-1,jte) , MAX(jde - spec_bdy_width,jts) , -1
2128        DO k = kds , kde - 1
2129        DO i = MAX(ids,its) , MIN(ide,ite)
2130           jj = jde - j
2131           space_bdy_ye(i,k,jj) = ( data3dnew(i,k,j) - data3dold(i,k,j) ) / time_diff
2132        END DO
2133        END DO
2134        END DO
2135     ELSE IF ( char_stagger .EQ. 'W' ) THEN
2136        DO j = MIN(jde-1,jte) , MAX(jde - spec_bdy_width,jts) , -1
2137        DO k = kds , kde
2138        DO i = MAX(ids,its) , MIN(ide-1,ite)
2139           jj = jde - j
2140           space_bdy_ye(i,k,jj) = ( data3dnew(i,k,j) - data3dold(i,k,j) ) / time_diff
2141        END DO
2142        END DO
2143        END DO
2144     ELSE IF ( char_stagger .EQ. 'M' ) THEN
2145        DO j = MIN(jde-1,jte) , MAX(jde - spec_bdy_width,jts) , -1
2146        DO k = kds , kde
2147        DO i = MAX(ids,its) , MIN(ide-1,ite)
2148           jj = jde - j
2149           space_bdy_ye(i,k,jj) = ( data3dnew(i,k,j) - data3dold(i,k,j) ) / time_diff
2150        END DO
2151        END DO
2152        END DO
2153     ELSE
2154        DO j = MIN(jde-1,jte) , MAX(jde - spec_bdy_width,jts) , -1
2155        DO k = kds , kde - 1
2156        DO i = MAX(ids,its) , MIN(ide-1,ite)
2157           jj = jde - j
2158           space_bdy_ye(i,k,jj) = ( data3dnew(i,k,j) - data3dold(i,k,j) ) / time_diff
2159        END DO
2160        END DO
2161        END DO
2162     END IF
2163     
2164  END SUBROUTINE stuff_bdytend_new
2165 
2166 !--- old versions for use with modules that use the old bdy data structures ---
2167 
2168  SUBROUTINE stuff_bdy_old ( data3d , space_bdy , char_stagger , &
2169                              ijds , ijde , spec_bdy_width , &
2170                              ids, ide, jds, jde, kds, kde , &
2171                              ims, ime, jms, jme, kms, kme , & 
2172                              its, ite, jts, jte, kts, kte )
2173  
2174  !  This routine puts the data in the 3d arrays into the proper locations
2175  !  for the lateral boundary arrays.
2176  
2177     USE module_state_description
2178     
2179     IMPLICIT NONE
2180  
2181     INTEGER , INTENT(IN) :: ids, ide, jds, jde, kds, kde
2182     INTEGER , INTENT(IN) :: ims, ime, jms, jme, kms, kme
2183     INTEGER , INTENT(IN) :: its, ite, jts, jte, kts, kte
2184     INTEGER , INTENT(IN) :: ijds , ijde , spec_bdy_width
2185     REAL , DIMENSION(ims:ime,kms:kme,jms:jme) , INTENT(IN) :: data3d
2186     REAL , DIMENSION(ijds:ijde,kds:kde,spec_bdy_width,4) , INTENT(OUT) :: space_bdy
2187     CHARACTER (LEN=1) , INTENT(IN) :: char_stagger
2188  
2189     INTEGER :: i , ii , j , jj , k
2190  
2191     !  There are four lateral boundary locations that are stored.
2192  
2193     !  X start boundary
2194  
2195     IF ( char_stagger .EQ. 'W' ) THEN
2196        DO j = MAX(jds,jts) , MIN(jde-1,jte)
2197        DO k = kds , kde
2198        DO i = MAX(ids,its) , MIN(ids + spec_bdy_width - 1,ite)
2199           space_bdy(j,k,i,P_XSB) = data3d(i,k,j)
2200        END DO
2201        END DO
2202        END DO
2203     ELSE IF ( char_stagger .EQ. 'M' ) THEN
2204        DO j = MAX(jds,jts) , MIN(jde-1,jte)
2205        DO k = kds , kde
2206        DO i = MAX(ids,its) , MIN(ids + spec_bdy_width - 1,ite)
2207           space_bdy(j,k,i,P_XSB) = data3d(i,k,j)
2208        END DO
2209        END DO
2210        END DO
2211     ELSE IF ( char_stagger .EQ. 'V' ) THEN
2212        DO j = MAX(jds,jts) , MIN(jde,jte)
2213        DO k = kds , kde - 1
2214        DO i = MAX(ids,its) , MIN(ids + spec_bdy_width - 1,ite)
2215           space_bdy(j,k,i,P_XSB) = data3d(i,k,j)
2216        END DO
2217        END DO
2218        END DO
2219     ELSE
2220        DO j = MAX(jds,jts) , MIN(jde-1,jte)
2221        DO k = kds , kde - 1
2222        DO i = MAX(ids,its) , MIN(ids + spec_bdy_width - 1,ite)
2223           space_bdy(j,k,i,P_XSB) = data3d(i,k,j)
2224        END DO
2225        END DO
2226        END DO
2227     END IF
2228  
2229     !  X end boundary
2230  
2231     IF      ( char_stagger .EQ. 'U' ) THEN
2232        DO j = MAX(jds,jts) , MIN(jde-1,jte)
2233        DO k = kds , kde - 1
2234        DO i = MIN(ide,ite) , MAX(ide - spec_bdy_width + 1,its) , -1
2235           ii = ide - i + 1
2236           space_bdy(j,k,ii,P_XEB) = data3d(i,k,j)
2237        END DO
2238        END DO
2239        END DO
2240     ELSE IF ( char_stagger .EQ. 'V' ) THEN
2241        DO j = MAX(jds,jts) , MIN(jde,jte)
2242        DO k = kds , kde - 1
2243        DO i = MIN(ide - 1,ite) , MAX(ide - spec_bdy_width,its) , -1
2244           ii = ide - i
2245           space_bdy(j,k,ii,P_XEB) = data3d(i,k,j)
2246        END DO
2247        END DO
2248        END DO
2249     ELSE IF ( char_stagger .EQ. 'W' ) THEN
2250        DO j = MAX(jds,jts) , MIN(jde-1,jte)
2251        DO k = kds , kde
2252        DO i = MIN(ide - 1,ite) , MAX(ide - spec_bdy_width,its) , -1
2253           ii = ide - i
2254           space_bdy(j,k,ii,P_XEB) = data3d(i,k,j)
2255        END DO
2256        END DO
2257        END DO
2258     ELSE IF ( char_stagger .EQ. 'M' ) THEN
2259        DO j = MAX(jds,jts) , MIN(jde-1,jte)
2260        DO k = kds , kde
2261        DO i = MIN(ide - 1,ite) , MAX(ide - spec_bdy_width,its) , -1
2262           ii = ide - i
2263           space_bdy(j,k,ii,P_XEB) = data3d(i,k,j)
2264        END DO
2265        END DO
2266        END DO
2267     ELSE
2268        DO j = MAX(jds,jts) , MIN(jde-1,jte)
2269        DO k = kds , kde - 1
2270        DO i = MIN(ide - 1,ite) , MAX(ide - spec_bdy_width,its) , -1
2271           ii = ide - i
2272           space_bdy(j,k,ii,P_XEB) = data3d(i,k,j)
2273        END DO
2274        END DO
2275        END DO
2276     END IF
2277  
2278     !  Y start boundary
2279  
2280     IF ( char_stagger .EQ. 'W' ) THEN
2281        DO j = MAX(jds,jts) , MIN(jds + spec_bdy_width - 1,jte)
2282        DO k = kds , kde
2283        DO i = MAX(ids,its) , MIN(ide-1,ite)
2284           space_bdy(i,k,j,P_YSB) = data3d(i,k,j)
2285        END DO
2286        END DO
2287        END DO
2288     ELSE IF ( char_stagger .EQ. 'M' ) THEN
2289        DO j = MAX(jds,jts) , MIN(jds + spec_bdy_width - 1,jte)
2290        DO k = kds , kde
2291        DO i = MAX(ids,its) , MIN(ide-1,ite)
2292           space_bdy(i,k,j,P_YSB) = data3d(i,k,j)
2293        END DO
2294        END DO
2295        END DO
2296     ELSE IF ( char_stagger .EQ. 'U' ) THEN
2297        DO j = MAX(jds,jts) , MIN(jds + spec_bdy_width - 1,jte)
2298        DO k = kds , kde - 1
2299        DO i = MAX(ids,its) , MIN(ide,ite)
2300           space_bdy(i,k,j,P_YSB) = data3d(i,k,j)
2301        END DO
2302        END DO
2303        END DO
2304     ELSE
2305        DO j = MAX(jds,jts) , MIN(jds + spec_bdy_width - 1,jte)
2306        DO k = kds , kde - 1
2307        DO i = MAX(ids,its) , MIN(ide-1,ite)
2308           space_bdy(i,k,j,P_YSB) = data3d(i,k,j)
2309        END DO
2310        END DO
2311        END DO
2312     END IF
2313  
2314     !  Y end boundary
2315  
2316     IF      ( char_stagger .EQ. 'V' ) THEN
2317        DO j = MIN(jde,jte) , MAX(jde - spec_bdy_width + 1,jts) , -1
2318        DO k = kds , kde - 1
2319        DO i = MAX(ids,its) , MIN(ide-1,ite)
2320           jj = jde - j + 1
2321           space_bdy(i,k,jj,P_YEB) = data3d(i,k,j)
2322        END DO
2323        END DO
2324        END DO
2325     ELSE IF ( char_stagger .EQ. 'U' ) THEN
2326        DO j = MIN(jde-1,jte) , MAX(jde - spec_bdy_width,jts) , -1
2327        DO k = kds , kde - 1
2328        DO i = MAX(ids,its) , MIN(ide,ite)
2329           jj = jde - j
2330           space_bdy(i,k,jj,P_YEB) = data3d(i,k,j)
2331        END DO
2332        END DO
2333        END DO
2334     ELSE IF ( char_stagger .EQ. 'W' ) THEN
2335        DO j = MIN(jde-1,jte) , MAX(jde - spec_bdy_width,jts) , -1
2336        DO k = kds , kde
2337        DO i = MAX(ids,its) , MIN(ide-1,ite)
2338           jj = jde - j
2339           space_bdy(i,k,jj,P_YEB) = data3d(i,k,j)
2340        END DO
2341        END DO
2342        END DO
2343     ELSE IF ( char_stagger .EQ. 'M' ) THEN
2344        DO j = MIN(jde-1,jte) , MAX(jde - spec_bdy_width,jts) , -1
2345        DO k = kds , kde
2346        DO i = MAX(ids,its) , MIN(ide-1,ite)
2347           jj = jde - j
2348           space_bdy(i,k,jj,P_YEB) = data3d(i,k,j)
2349        END DO
2350        END DO
2351        END DO
2352     ELSE
2353        DO j = MIN(jde-1,jte) , MAX(jde - spec_bdy_width,jts) , -1
2354        DO k = kds , kde - 1
2355        DO i = MAX(ids,its) , MIN(ide-1,ite)
2356           jj = jde - j
2357           space_bdy(i,k,jj,P_YEB) = data3d(i,k,j)
2358        END DO
2359        END DO
2360        END DO
2361     END IF
2362     
2363  END SUBROUTINE stuff_bdy_old
2364  
2365  SUBROUTINE stuff_bdytend_old ( data3dnew , data3dold , time_diff , space_bdy , char_stagger , &
2366                              ijds , ijde , spec_bdy_width , &
2367                              ids, ide, jds, jde, kds, kde , &
2368                              ims, ime, jms, jme, kms, kme , & 
2369                              its, ite, jts, jte, kts, kte )
2370  
2371  !  This routine puts the tendency data into the proper locations
2372  !  for the lateral boundary arrays.
2373  
2374     USE module_state_description
2375     
2376     IMPLICIT NONE
2377  
2378     INTEGER , INTENT(IN) :: ids, ide, jds, jde, kds, kde
2379     INTEGER , INTENT(IN) :: ims, ime, jms, jme, kms, kme
2380     INTEGER , INTENT(IN) :: its, ite, jts, jte, kts, kte
2381     INTEGER , INTENT(IN) :: ijds , ijde , spec_bdy_width
2382     REAL , DIMENSION(ims:ime,kms:kme,jms:jme) , INTENT(IN) :: data3dnew , data3dold
2383 !    REAL , DIMENSION(:,:,:,:) , INTENT(OUT) :: space_bdy
2384     REAL , DIMENSION(ijds:ijde,kds:kde,spec_bdy_width,4) , INTENT(OUT) :: space_bdy
2385     CHARACTER (LEN=1) , INTENT(IN) :: char_stagger
2386     REAL , INTENT(IN) :: time_diff ! seconds
2387  
2388     INTEGER :: i , ii , j , jj , k
2389  
2390     !  There are four lateral boundary locations that are stored.
2391  
2392     !  X start boundary
2393  
2394     IF ( char_stagger .EQ. 'W' ) THEN
2395        DO j = MAX(jds,jts) , MIN(jde-1,jte)
2396        DO k = kds , kde
2397        DO i = MAX(ids,its) , MIN(ids + spec_bdy_width - 1,ite)
2398           space_bdy(j,k,i,P_XSB) = ( data3dnew(i,k,j) - data3dold(i,k,j) ) / time_diff
2399 !         space_bdy(j,k,i,P_XSB) = 0. ! zeroout
2400        END DO
2401        END DO
2402        END DO
2403     ELSE IF ( char_stagger .EQ. 'M' ) THEN
2404        DO j = MAX(jds,jts) , MIN(jde-1,jte)
2405        DO k = kds , kde
2406        DO i = MAX(ids,its) , MIN(ids + spec_bdy_width - 1,ite)
2407           space_bdy(j,k,i,P_XSB) = ( data3dnew(i,k,j) - data3dold(i,k,j) ) / time_diff
2408 !         space_bdy(j,k,i,P_XSB) = 0. ! zeroout
2409        END DO
2410        END DO
2411        END DO
2412     ELSE IF ( char_stagger .EQ. 'V' ) THEN
2413        DO j = MAX(jds,jts) , MIN(jde,jte)
2414        DO k = kds , kde - 1
2415        DO i = MAX(ids,its) , MIN(ids + spec_bdy_width - 1,ite)
2416           space_bdy(j,k,i,P_XSB) = ( data3dnew(i,k,j) - data3dold(i,k,j) ) / time_diff
2417 !         space_bdy(j,k,i,P_XSB) = 0. ! zeroout
2418        END DO
2419        END DO
2420        END DO
2421     ELSE
2422        DO j = MAX(jds,jts) , MIN(jde-1,jte)
2423        DO k = kds , kde - 1
2424        DO i = MAX(ids,its) , MIN(ids + spec_bdy_width - 1,ite)
2425           space_bdy(j,k,i,P_XSB) = ( data3dnew(i,k,j) - data3dold(i,k,j) ) / time_diff
2426 !         space_bdy(j,k,i,P_XSB) = 0. ! zeroout
2427        END DO
2428        END DO
2429        END DO
2430     END IF
2431  
2432     !  X end boundary
2433  
2434     IF      ( char_stagger .EQ. 'U' ) THEN
2435        DO j = MAX(jds,jts) , MIN(jde-1,jte)
2436        DO k = kds , kde - 1
2437        DO i = MIN(ide,ite) , MAX(ide - spec_bdy_width + 1,its) , -1
2438           ii = ide - i + 1
2439           space_bdy(j,k,ii,P_XEB) = ( data3dnew(i,k,j) - data3dold(i,k,j) ) / time_diff
2440 !         space_bdy(j,k,ii,P_XEB) = 0. ! zeroout
2441        END DO
2442        END DO
2443        END DO
2444     ELSE IF ( char_stagger .EQ. 'V' ) THEN
2445        DO j = MAX(jds,jts) , MIN(jde,jte)
2446        DO k = kds , kde - 1
2447        DO i = MIN(ide - 1,ite) , MAX(ide - spec_bdy_width,its) , -1
2448           ii = ide - i
2449           space_bdy(j,k,ii,P_XEB) = ( data3dnew(i,k,j) - data3dold(i,k,j) ) / time_diff
2450 !         space_bdy(j,k,ii,P_XEB) = 0. ! zeroout
2451        END DO
2452        END DO
2453        END DO
2454     ELSE IF ( char_stagger .EQ. 'W' ) THEN
2455        DO j = MAX(jds,jts) , MIN(jde-1,jte)
2456        DO k = kds , kde
2457        DO i = MIN(ide - 1,ite) , MAX(ide - spec_bdy_width,its) , -1
2458           ii = ide - i
2459           space_bdy(j,k,ii,P_XEB) = ( data3dnew(i,k,j) - data3dold(i,k,j) ) / time_diff
2460 !         space_bdy(j,k,ii,P_XEB) = 0. ! zeroout
2461        END DO
2462        END DO
2463        END DO
2464     ELSE IF ( char_stagger .EQ. 'M' ) THEN
2465        DO j = MAX(jds,jts) , MIN(jde-1,jte)
2466        DO k = kds , kde
2467        DO i = MIN(ide - 1,ite) , MAX(ide - spec_bdy_width,its) , -1
2468           ii = ide - i
2469           space_bdy(j,k,ii,P_XEB) = ( data3dnew(i,k,j) - data3dold(i,k,j) ) / time_diff
2470 !         space_bdy(j,k,ii,P_XEB) = 0. ! zeroout
2471        END DO
2472        END DO
2473        END DO
2474     ELSE
2475        DO j = MAX(jds,jts) , MIN(jde-1,jte)
2476        DO k = kds , kde - 1
2477        DO i = MIN(ide - 1,ite) , MAX(ide - spec_bdy_width,its) , -1
2478           ii = ide - i
2479           space_bdy(j,k,ii,P_XEB) = ( data3dnew(i,k,j) - data3dold(i,k,j) ) / time_diff
2480 !         space_bdy(j,k,ii,P_XEB) = 0. ! zeroout
2481        END DO
2482        END DO
2483        END DO
2484     END IF
2485  
2486     !  Y start boundary
2487  
2488     IF ( char_stagger .EQ. 'W' ) THEN
2489        DO j = MAX(jds,jts) , MIN(jds + spec_bdy_width - 1,jte)
2490        DO k = kds , kde
2491        DO i = MAX(ids,its) , MIN(ide-1,ite)
2492           space_bdy(i,k,j,P_YSB) = ( data3dnew(i,k,j) - data3dold(i,k,j) ) / time_diff
2493 !         space_bdy(i,k,j,P_YSB) = 0. ! zeroout
2494        END DO
2495        END DO
2496        END DO
2497     ELSE IF ( char_stagger .EQ. 'M' ) THEN
2498        DO j = MAX(jds,jts) , MIN(jds + spec_bdy_width - 1,jte)
2499        DO k = kds , kde
2500        DO i = MAX(ids,its) , MIN(ide-1,ite)
2501           space_bdy(i,k,j,P_YSB) = ( data3dnew(i,k,j) - data3dold(i,k,j) ) / time_diff
2502 !         space_bdy(i,k,j,P_YSB) = 0. ! zeroout
2503        END DO
2504        END DO
2505        END DO
2506     ELSE IF ( char_stagger .EQ. 'U' ) THEN
2507        DO j = MAX(jds,jts) , MIN(jds + spec_bdy_width - 1,jte)
2508        DO k = kds , kde - 1
2509        DO i = MAX(ids,its) , MIN(ide,ite)
2510           space_bdy(i,k,j,P_YSB) = ( data3dnew(i,k,j) - data3dold(i,k,j) ) / time_diff
2511 !         space_bdy(i,k,j,P_YSB) = 0. ! zeroout
2512        END DO
2513        END DO
2514        END DO
2515     ELSE
2516        DO j = MAX(jds,jts) , MIN(jds + spec_bdy_width - 1,jte)
2517        DO k = kds , kde - 1
2518        DO i = MAX(ids,its) , MIN(ide-1,ite)
2519           space_bdy(i,k,j,P_YSB) = ( data3dnew(i,k,j) - data3dold(i,k,j) ) / time_diff
2520 !         space_bdy(i,k,j,P_YSB) = 0. ! zeroout
2521        END DO
2522        END DO
2523        END DO
2524     END IF
2525  
2526     !  Y end boundary
2527  
2528     IF      ( char_stagger .EQ. 'V' ) THEN
2529        DO j = MIN(jde,jte) , MAX(jde - spec_bdy_width + 1,jts) , -1
2530        DO k = kds , kde - 1
2531        DO i = MAX(ids,its) , MIN(ide-1,ite)
2532           jj = jde - j + 1
2533           space_bdy(i,k,jj,P_YEB) = ( data3dnew(i,k,j) - data3dold(i,k,j) ) / time_diff
2534 !         space_bdy(i,k,jj,P_YEB) = 0. ! zeroout
2535        END DO
2536        END DO
2537        END DO
2538     ELSE IF ( char_stagger .EQ. 'U' ) THEN
2539        DO j = MIN(jde-1,jte) , MAX(jde - spec_bdy_width,jts) , -1
2540        DO k = kds , kde - 1
2541        DO i = MAX(ids,its) , MIN(ide,ite)
2542           jj = jde - j
2543           space_bdy(i,k,jj,P_YEB) = ( data3dnew(i,k,j) - data3dold(i,k,j) ) / time_diff
2544 !         space_bdy(i,k,jj,P_YEB) = 0. ! zeroout
2545        END DO
2546        END DO
2547        END DO
2548     ELSE IF ( char_stagger .EQ. 'W' ) THEN
2549        DO j = MIN(jde-1,jte) , MAX(jde - spec_bdy_width,jts) , -1
2550        DO k = kds , kde
2551        DO i = MAX(ids,its) , MIN(ide-1,ite)
2552           jj = jde - j
2553           space_bdy(i,k,jj,P_YEB) = ( data3dnew(i,k,j) - data3dold(i,k,j) ) / time_diff
2554 !         space_bdy(i,k,jj,P_YEB) = 0. ! zeroout
2555        END DO
2556        END DO
2557        END DO
2558     ELSE IF ( char_stagger .EQ. 'M' ) THEN
2559        DO j = MIN(jde-1,jte) , MAX(jde - spec_bdy_width,jts) , -1
2560        DO k = kds , kde
2561        DO i = MAX(ids,its) , MIN(ide-1,ite)
2562           jj = jde - j
2563           space_bdy(i,k,jj,P_YEB) = ( data3dnew(i,k,j) - data3dold(i,k,j) ) / time_diff
2564 !         space_bdy(i,k,jj,P_YEB) = 0. ! zeroout
2565        END DO
2566        END DO
2567        END DO
2568     ELSE
2569        DO j = MIN(jde-1,jte) , MAX(jde - spec_bdy_width,jts) , -1
2570        DO k = kds , kde - 1
2571        DO i = MAX(ids,its) , MIN(ide-1,ite)
2572           jj = jde - j
2573           space_bdy(i,k,jj,P_YEB) = ( data3dnew(i,k,j) - data3dold(i,k,j) ) / time_diff
2574 !         space_bdy(i,k,jj,P_YEB) = 0. ! zeroout
2575        END DO
2576        END DO
2577        END DO
2578     END IF
2579     
2580  END SUBROUTINE stuff_bdytend_old
2581 
2582 END MODULE module_bc
2583 
2584 SUBROUTINE get_bdyzone_x ( bzx )
2585   USE module_bc
2586   IMPLICIT NONE
2587   INTEGER bzx
2588   bzx = bdyzone_x
2589 END SUBROUTINE get_bdyzone_x
2590 
2591 SUBROUTINE get_bdyzone_y ( bzy)
2592   USE module_bc
2593   IMPLICIT NONE
2594   INTEGER bzy
2595   bzy = bdyzone_y
2596 END SUBROUTINE get_bdyzone_y
2597 
2598 SUBROUTINE get_bdyzone ( bz)
2599   USE module_bc
2600   IMPLICIT NONE
2601   INTEGER bz
2602   bz = bdyzone
2603 END SUBROUTINE get_bdyzone
2604