bdy_2d.F90 29.5 KB
Newer Older
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15
#include "cppdefs.h"
!-----------------------------------------------------------------------
!BOP
!
! !MODULE:  bdy_2d - 2D boundary conditions \label{bdy-2d}
!
! !INTERFACE:
   module bdy_2d
!
! !DESCRIPTION:
!
! Here, the two-dimensional boundary
! conditions for sea surface elevation and transports are handled.
!
! !USES:
Knut's avatar
Knut committed
16
   use parameters,only: g
17 18
   use halo_zones, only : z_TAG,H_TAG,U_TAG,V_TAG
   use domain, only: imin,jmin,imax,jmax,kmax,H,az,au,av
Knut's avatar
Knut committed
19
   use domain, only: nsbv,nsbvl,nbdy,NWB,NNB,NEB,NSB,bdy_index,bdy_index_l
20
   use domain, only: bdy_2d_desc,bdy_2d_type
Knut's avatar
Knut committed
21
   use domain, only: need_2d_bdy_elev,need_2d_bdy_u,need_2d_bdy_v
22 23
   use domain, only: wi,wfj,wlj,nj,nfi,nli,ei,efj,elj,sj,sfi,sli
   use domain, only: min_depth
Knut's avatar
Knut committed
24
   use time, only: write_time_string,timestr
25
   use variables_2d, only: dtm,z,zo,D,U,DU,V,DV
26
   use domain, only: dxu,dyv
27
   use exceptions
28 29 30 31 32 33 34

   IMPLICIT NONE
!
   private
!
! !PUBLIC DATA MEMBERS:
   public init_bdy_2d, do_bdy_2d
Knut's avatar
Knut committed
35
   character(len=PATH_MAX),public       :: bdyfile_2d
Knut's avatar
Knut committed
36 37 38
   integer,public                       :: bdyfmt_2d
   integer,public                       :: bdy2d_ramp=-1
   integer,public                       :: bdy2d_sponge_size=0
Knut's avatar
Knut committed
39
   REALTYPE,dimension(:),pointer,public :: bdy_data,bdy_data_u,bdy_data_v
40
!
Knut's avatar
Knut committed
41 42
! !PRIVATE DATA MEMBERS:
   private bdy2d_active,bdy2d_need_elev,bdy2d_need_vel
Knut's avatar
Knut committed
43 44
   REALTYPE                             :: ramp=_ONE_
   logical                              :: ramp_is_active=.false.
Knut's avatar
Knut committed
45
   REALTYPE,dimension(:),allocatable    :: sp(:)
Knut's avatar
Knut committed
46
!
47 48 49 50 51 52 53 54 55 56 57 58 59 60 61
! !REVISION HISTORY:
!  Original author(s): Karsten Bolding & Hans Burchard
!
!EOP
!-----------------------------------------------------------------------

   contains

!-----------------------------------------------------------------------
!BOP
!
! !IROUTINE: init_bdy_2d - initialising 2D boundary conditions
! \label{sec-init-bdy-2d}
!
! !INTERFACE:
Knut's avatar
Knut committed
62
   subroutine init_bdy_2d(bdy2d,hotstart)
63 64 65 66 67 68
!
! !DESCRIPTION:
!
! !USES:
   IMPLICIT NONE
!
Knut's avatar
Knut committed
69 70 71 72 73 74
! !INPUT PARAMETERS:
   logical, intent(in)                 :: hotstart
!
! !INPUT/OUTPUT PARAMETERS:
   logical, intent(inout)              :: bdy2d
!
75
! !LOCAL VARIABLES:
76
   integer :: i,j,n,l,shift, rc
77 78 79 80 81 82 83 84 85
!EOP
!-------------------------------------------------------------------------
!BOC
#ifdef DEBUG
   integer, save :: Ncall = 0
   Ncall = Ncall+1
   write(debug,*) 'init_bdy_2d() # ',Ncall
#endif

Knut's avatar
Knut committed
86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 126 127 128 129 130 131 132 133 134 135 136 137
   LEVEL2 'init_bdy_2d'

   if (bdy2d) then

      do l=1,nbdy
         if (bdy2d_need_elev(bdy_2d_type(l))) then
            need_2d_bdy_elev = .true.
            exit
         end if
      end do

      l = 0
      do n = 1,NWB
         l = l+1
         if (bdy2d_need_vel(bdy_2d_type(l))) then
            need_2d_bdy_u = .true.
            exit
         end if
      end do
      if (.not. need_2d_bdy_u) then
         l = l + NNB
         do n = 1,NEB
            l = l+1
            if (bdy2d_need_vel(bdy_2d_type(l))) then
               need_2d_bdy_u = .true.
               exit
            end if
         end do
      end if

      l = NWB
      do n = 1,NNB
         l = l+1
         if (bdy2d_need_vel(bdy_2d_type(l))) then
            need_2d_bdy_v = .true.
            exit
         end if
      end do
      if (.not. need_2d_bdy_v) then
         l = l + NEB
         do n = 1,NSB
            l = l+1
            if (bdy2d_need_vel(bdy_2d_type(l))) then
               need_2d_bdy_v = .true.
               exit
            end if
         end do
      end if

      if (need_2d_bdy_elev .or. need_2d_bdy_u .or. need_2d_bdy_v) then
         LEVEL3 'bdyfile_2d=',TRIM(bdyfile_2d)
         LEVEL3 'bdyfmt_2d=',bdyfmt_2d
Knut's avatar
Knut committed
138 139
         if (bdy2d_ramp .gt. 1) then
            LEVEL3 'bdy2d_ramp=',bdy2d_ramp
Knut's avatar
Knut committed
140
            ramp_is_active = .true.
Knut's avatar
Knut committed
141
            if (hotstart) then
Knut's avatar
Knut committed
142
               LEVEL4 'WARNING: hotstart is .true. AND bdy2d_ramp .gt. 1'
Knut's avatar
Knut committed
143 144
               LEVEL4 'WARNING: .. be sure you know what you are doing ..'
            end if
Knut's avatar
Knut committed
145
         end if
146
         if (need_2d_bdy_elev) then
Knut's avatar
Knut committed
147
            allocate(bdy_data(nsbvl),stat=rc)
148 149 150
            if (rc /= 0) stop 'init_bdy_2d: Error allocating memory (bdy_data)'
         end if
         if (need_2d_bdy_u) then
Knut's avatar
Knut committed
151
            allocate(bdy_data_u(nsbvl),stat=rc)
152 153 154
            if (rc /= 0) stop 'init_bdy_2d: Error allocating memory (bdy_data_u)'
         end if
         if (need_2d_bdy_v) then
Knut's avatar
Knut committed
155
            allocate(bdy_data_v(nsbvl),stat=rc)
156 157
            if (rc /= 0) stop 'init_bdy_2d: Error allocating memory (bdy_data_v)'
         end if
Knut's avatar
Knut committed
158 159 160 161 162 163 164
      else
         bdy2d = .false.
      end if

   else

      do l=1,nbdy
Knut's avatar
Knut committed
165
         if (bdy2d_active(bdy_2d_type(l))) then
166 167
            LEVEL3 'bdy2d=F resets local 2D bdy #',l
            LEVEL4 'old: ',trim(bdy_2d_desc(bdy_2d_type(l)))
Knut's avatar
Knut committed
168
            bdy_2d_type(l) = CONSTANT
169
            LEVEL4 'new: ',trim(bdy_2d_desc(bdy_2d_type(l)))
Knut's avatar
Knut committed
170 171 172 173
         end if
      end do

   end if
174

Knut's avatar
Knut committed
175 176 177 178 179 180 181 182 183 184 185 186
   if (bdy2d_sponge_size .gt. 0) then
      allocate(sp(bdy2d_sponge_size),stat=rc)
      if (rc /= 0) stop 'init_bdy_2d: Error allocating memory (sp)'

!        Sponge layer factors according to Martinsen and Engedahl, 1987.
!        Note (KK): factor=1 (bdy cell) does not count for sponge size
!                   (in contrast to earlier GETM)
      LEVEL3 "sponge layer factors for surface elevation:"
      do i=1,bdy2d_sponge_size
         sp(i) = ((_ONE_+bdy2d_sponge_size-i)/(_ONE_+bdy2d_sponge_size))**2
         LEVEL4 "sp(",i,")=",real(sp(i))
      end do
187 188
   else
      bdy2d_sponge_size = 0
Knut's avatar
Knut committed
189 190
   end if

191 192 193 194 195 196 197 198 199 200 201 202 203 204 205 206 207 208 209 210 211 212 213 214 215 216 217 218 219 220 221 222 223 224 225 226 227 228 229 230 231 232 233 234 235 236 237 238 239 240 241 242 243 244 245 246 247 248 249 250 251 252 253 254 255 256 257 258 259 260 261 262 263 264 265 266 267 268 269 270 271 272 273 274 275 276 277 278 279 280 281 282 283 284
   l = 0
   do n = 1,NWB
      l = l+1
      i = wi(n)
      do j = wfj(n),wlj(n)
         if (az(i,j) .eq. 2) then
            select case (bdy_2d_type(l))
               case (CONSTANT,CLAMPED_ELEV,CLAMPED)
                  shift = bdy2d_sponge_size
               case (ZERO_GRADIENT,SOMMERFELD,CLAMPED_VEL,FLATHER_VEL)
                  shift = 1
               case default
                  shift = 0
            end select
            if ( i+shift .gt. imax ) then
               FATAL 'local western bdy #',n,'too close to eastern subdomain edge'
               call getm_error('init_bdy_2d()', &
                               'western open bdy too close to eastern subdomain edge')
            else
               exit
            end if
         end if
      end do
   end do
   do n = 1,NNB
      l = l+1
      j = nj(n)
      do i = nfi(n),nli(n)
         if (az(i,j) .eq. 2) then
            select case (bdy_2d_type(l))
               case (CONSTANT,CLAMPED_ELEV,CLAMPED)
                  shift = bdy2d_sponge_size
               case (ZERO_GRADIENT,SOMMERFELD,CLAMPED_VEL,FLATHER_VEL)
                  shift = 1
               case default
                  shift = 0
            end select
            if ( j-shift .lt. jmin ) then
               FATAL 'local northern bdy #',n,'too close to southern subdomain edge'
               call getm_error('init_bdy_2d()', &
                               'northern open bdy too close to southern subdomain edge')
            else
               exit
            end if
         end if
      end do
   end do
   do n = 1,NEB
      l = l+1
      i = ei(n)
      do j = efj(n),elj(n)
         if (az(i,j) .eq. 2) then
            select case (bdy_2d_type(l))
               case (CONSTANT,CLAMPED_ELEV,CLAMPED)
                  shift = bdy2d_sponge_size
               case (ZERO_GRADIENT,SOMMERFELD,CLAMPED_VEL,FLATHER_VEL)
                  shift = 1
               case default
                  shift = 0
            end select
            if ( i-shift .lt. imin ) then
               FATAL 'local eastern bdy #',n,'too close to western subdomain edge'
               call getm_error('init_bdy_2d()', &
                               'eastern open bdy too close to western subdomain edge')
            else
               exit
            end if
         end if
      end do
   end do
   do n = 1,NSB
      l = l+1
      j = sj(n)
      do i = sfi(n),sli(n)
         if (az(i,j) .eq. 2) then
            select case (bdy_2d_type(l))
               case (CONSTANT,CLAMPED_ELEV,CLAMPED)
                  shift = bdy2d_sponge_size
               case (ZERO_GRADIENT,SOMMERFELD,CLAMPED_VEL,FLATHER_VEL)
                  shift = 1
               case default
                  shift = 0
            end select
            if ( j+shift .gt. jmax ) then
               FATAL 'local southern bdy #',n,'too close to northern subdomain edge'
               call getm_error('init_bdy_2d()', &
                               'southern open bdy too close to northern subdomain edge')
            else
               exit
            end if
         end if
      end do
   end do

Knut's avatar
Knut committed
285

286 287 288 289 290 291 292 293 294 295 296 297 298 299 300 301 302 303 304 305 306 307 308 309 310 311 312 313
#ifdef DEBUG
   write(debug,*) 'Leaving init_bdy_2d()'
   write(debug,*)
#endif
   return
   end subroutine init_bdy_2d
!EOC

!-----------------------------------------------------------------------
!BOP
!
! !IROUTINE:  do_bdy_2d  - updating 2D boundary conditions
! \label{sec-do-bdy-2d}
!
! !INTERFACE:
   subroutine do_bdy_2d(loop,tag)
!
! !DESCRIPTION:
!
! !USES:
   IMPLICIT NONE
!
! !INPUT PARAMETERS:
   integer, intent(in)                 :: loop,tag
!
! !INPUT/OUTPUT PARAMETERS:
!
! !LOCAL VARIABLES:
Knut's avatar
Knut committed
314
   REALTYPE                  :: cfl,depth,a
Knut's avatar
Knut committed
315
   integer                   :: i,j,k,ii,jj,kl,l,n
Knut's avatar
Knut committed
316
   REALTYPE, parameter       :: theta = _HALF_
317 318 319 320 321 322 323 324 325 326
!
!EOP
!-----------------------------------------------------------------------
!BOC
#ifdef DEBUG
   integer, save :: Ncall = 0
   Ncall = Ncall+1
   write(debug,*) 'do_bdy_2d() # ',Ncall
#endif

Knut's avatar
Knut committed
327
#if 0
328 329 330 331 332 333 334
   select case (bdyfmt_2d)
      case (NO_DATA)
      case (ANALYTICAL)
      case (ASCII)
      case (NETCDF)
!        Read in get_2d_bdy() via get_2d_bdy_ncdf()
      case default
Knut's avatar
Knut committed
335
         stop 'do_bdy_2d(): invalid bdyfmt_2d'
336
   end select
Knut's avatar
Knut committed
337
#endif
338 339 340

!  Data read - do time interpolation

Knut's avatar
Knut committed
341 342 343 344 345 346 347 348 349 350 351
   if (ramp_is_active) then
      if (loop .ge. bdy2d_ramp) then
         ramp = _ONE_
         ramp_is_active = .false.
         STDERR LINE
         call write_time_string()
         LEVEL3 timestr,': finished bdy2d_ramp=',bdy2d_ramp
         STDERR LINE
      else
         ramp = _ONE_*loop/bdy2d_ramp
      end if
Knut's avatar
Knut committed
352
   end if
353

Knut's avatar
Knut committed
354 355 356 357 358 359 360 361
   select case (tag)

      case (z_TAG,H_TAG)

         l = 0
         do n = 1,NWB
            l = l+1
            k = bdy_index(l)
Knut's avatar
Knut committed
362
            kl = bdy_index_l(l)
Knut's avatar
Knut committed
363
            i = wi(n)
Knut's avatar
Knut committed
364 365 366
            select case (bdy_2d_type(l))
               case (ZERO_GRADIENT,CLAMPED_VEL,FLATHER_VEL)
                  do j = wfj(n),wlj(n)
Knut's avatar
Knut committed
367
                     if ( az(i+1,j) .ne. 0 ) z(i,j) = z(i+1,j)
Knut's avatar
Knut committed
368 369 370
                  end do
               case (SOMMERFELD)
                  do j = wfj(n),wlj(n)
Knut's avatar
Knut committed
371
                     if ( az(i+1,j) .ne. 0 ) then
Knut's avatar
Knut committed
372 373 374 375 376 377
                     cfl = sqrt(g*_HALF_*(D(i,j)+D(i+1,j)))*dtm/DXU
                     z(i,j) = (                                             &
                                (_ONE_ - _TWO_*cfl*(_ONE_-theta))*z (i  ,j) &
                               +(_ONE_ + _TWO_*cfl*(_ONE_-theta))*zo(i+1,j) &
                               -(_ONE_ - _TWO_*cfl*theta        )*z (i+1,j) &
                              )/(_ONE_ + _TWO_*cfl*theta        )
Knut's avatar
Knut committed
378
                     end if
Knut's avatar
Knut committed
379
                  end do
Knut's avatar
Knut committed
380
               case (CLAMPED_ELEV,CLAMPED)
Knut's avatar
Knut committed
381
                  do j = wfj(n),wlj(n)
Knut's avatar
Knut committed
382
                     z(i,j) = max(ramp*bdy_data(kl),-H(i,j)+min_depth)
Knut's avatar
Knut committed
383
                     k = k+1
Knut's avatar
Knut committed
384
                     kl = kl + 1
Knut's avatar
Knut committed
385 386 387
                  end do
               case (FLATHER_ELEV)
                  do j = wfj(n),wlj(n)
Knut's avatar
Knut committed
388
                     if ( az(i+1,j) .ne. 0 ) then
Knut's avatar
Knut committed
389 390 391 392
!                    Note (KK): approximate interface depths at vel-time stage
!                               by spatial mean at last sse-time stage
                     depth = _HALF_*(D(i,j)+D(i+1,j))
!                    Note (KK): note approximation of sse at vel-time stage
Knut's avatar
Knut committed
393 394 395 396 397
                     a = _TWO_/sqrt(g*depth)*(U(i,j)-ramp*bdy_data_u(kl)*depth)
                     else
                     a = _ZERO_
                     end if
                     z(i,j) = max(ramp*bdy_data(kl)-a,-H(i,j)+min_depth)
Knut's avatar
Knut committed
398
                     k = k+1
Knut's avatar
Knut committed
399
                     kl = kl + 1
Knut's avatar
Knut committed
400 401
                  end do
            end select
Knut's avatar
Knut committed
402 403 404 405
         end do
         do n = 1,NNB
            l = l+1
            k = bdy_index(l)
Knut's avatar
Knut committed
406
            kl = bdy_index_l(l)
Knut's avatar
Knut committed
407
            j = nj(n)
Knut's avatar
Knut committed
408 409 410
            select case (bdy_2d_type(l))
               case (ZERO_GRADIENT,CLAMPED_VEL,FLATHER_VEL)
                  do i = nfi(n),nli(n)
Knut's avatar
Knut committed
411
                     if ( az(i,j-1) .ne. 0 ) z(i,j) = z(i,j-1)
Knut's avatar
Knut committed
412 413 414
                  end do
               case (SOMMERFELD)
                  do i = nfi(n),nli(n)
Knut's avatar
Knut committed
415
                     if ( az(i,j-1) .ne. 0 ) then
Knut's avatar
Knut committed
416 417 418 419 420 421
                     cfl = sqrt(g*_HALF_*(D(i,j-1)+D(i,j)))*dtm/DYVJM1
                     z(i,j) = (                                             &
                                (_ONE_ - _TWO_*cfl*(_ONE_-theta))*z (i,j  ) &
                               +(_ONE_ + _TWO_*cfl*(_ONE_-theta))*zo(i,j-1) &
                               -(_ONE_ - _TWO_*cfl*theta        )*z (i,j-1) &
                              )/(_ONE_ + _TWO_*cfl*theta        )
Knut's avatar
Knut committed
422
                     end if
Knut's avatar
Knut committed
423
                  end do
Knut's avatar
Knut committed
424
               case (CLAMPED_ELEV,CLAMPED)
Knut's avatar
Knut committed
425
                  do i = nfi(n),nli(n)
Knut's avatar
Knut committed
426
                     z(i,j) = max(ramp*bdy_data(kl),-H(i,j)+min_depth)
Knut's avatar
Knut committed
427
                     k = k+1
Knut's avatar
Knut committed
428
                     kl = kl + 1
Knut's avatar
Knut committed
429 430 431
                  end do
               case (FLATHER_ELEV)
                  do i = nfi(n),nli(n)
Knut's avatar
Knut committed
432
                     if ( az(i,j-1) .ne. 0 ) then
Knut's avatar
Knut committed
433 434 435 436
!                    Note (KK): approximate interface depths at vel-time stage
!                               by spatial mean at last sse-time stage
                     depth = _HALF_*(D(i,j-1)+D(i,j))
!                    Note (KK): note approximation of sse at vel-time stage
Knut's avatar
Knut committed
437 438 439 440 441
                     a = _TWO_/sqrt(g*depth)*(V(i,j-1)-ramp*bdy_data_v(kl)*depth)
                     else
                     a = _ZERO_
                     end if
                     z(i,j) = max(ramp*bdy_data(kl)+a,-H(i,j)+min_depth)
Knut's avatar
Knut committed
442
                     k = k+1
Knut's avatar
Knut committed
443
                     kl = kl + 1
Knut's avatar
Knut committed
444 445
                  end do
            end select
Knut's avatar
Knut committed
446 447 448 449
         end do
         do n = 1,NEB
            l = l+1
            k = bdy_index(l)
Knut's avatar
Knut committed
450
            kl = bdy_index_l(l)
Knut's avatar
Knut committed
451
            i = ei(n)
Knut's avatar
Knut committed
452 453 454
            select case (bdy_2d_type(l))
               case (ZERO_GRADIENT,CLAMPED_VEL,FLATHER_VEL)
                  do j = efj(n),elj(n)
Knut's avatar
Knut committed
455
                     if ( az(i-1,j) .ne. 0 ) z(i,j) = z(i-1,j)
Knut's avatar
Knut committed
456 457 458
                  end do
               case (SOMMERFELD)
                  do j = efj(n),elj(n)
Knut's avatar
Knut committed
459
                     if ( az(i-1,j) .ne. 0 ) then
Knut's avatar
Knut committed
460 461 462 463 464 465
                     cfl = sqrt(g*_HALF_*(D(i-1,j)+D(i,j)))*dtm/DXUIM1
                     z(i,j) = (                                             &
                                (_ONE_ - _TWO_*cfl*(_ONE_-theta))*z (i  ,j) &
                               +(_ONE_ + _TWO_*cfl*(_ONE_-theta))*zo(i-1,j) &
                               -(_ONE_ - _TWO_*cfl*theta        )*z (i-1,j) &
                              )/(_ONE_ + _TWO_*cfl*theta        )
Knut's avatar
Knut committed
466
                     end if
Knut's avatar
Knut committed
467
                  end do
Knut's avatar
Knut committed
468
               case (CLAMPED_ELEV,CLAMPED)
Knut's avatar
Knut committed
469
                  do j = efj(n),elj(n)
Knut's avatar
Knut committed
470
                     z(i,j) = max(ramp*bdy_data(kl),-H(i,j)+min_depth)
Knut's avatar
Knut committed
471
                     k = k+1
Knut's avatar
Knut committed
472
                     kl = kl + 1
Knut's avatar
Knut committed
473 474 475
                  end do
               case (FLATHER_ELEV)
                  do j = efj(n),elj(n)
Knut's avatar
Knut committed
476
                     if ( az(i-1,j) .ne. 0 ) then
Knut's avatar
Knut committed
477 478 479 480
!                    Note (KK): approximate interface depths at vel-time stage
!                               by spatial mean at last sse-time stage
                     depth = _HALF_*(D(i-1,j)+D(i,j))
!                    Note (KK): note approximation of sse at vel-time stage
Knut's avatar
Knut committed
481 482 483 484 485
                     a = _TWO_/sqrt(g*depth)*(U(i-1,j)-ramp*bdy_data_u(kl)*depth)
                     else
                     a = _ZERO_
                     end if
                     z(i,j) = max(ramp*bdy_data(kl)+a,-H(i,j)+min_depth)
Knut's avatar
Knut committed
486
                     k = k+1
Knut's avatar
Knut committed
487
                     kl = kl + 1
Knut's avatar
Knut committed
488 489
                  end do
            end select
Knut's avatar
Knut committed
490 491 492 493
         end do
         do n = 1,NSB
            l = l+1
            k = bdy_index(l)
Knut's avatar
Knut committed
494
            kl = bdy_index_l(l)
Knut's avatar
Knut committed
495
            j = sj(n)
Knut's avatar
Knut committed
496 497 498
            select case (bdy_2d_type(l))
               case (ZERO_GRADIENT,CLAMPED_VEL,FLATHER_VEL)
                  do i = sfi(n),sli(n)
Knut's avatar
Knut committed
499
                     if ( az(i,j+1) .ne. 0 ) z(i,j) = z(i,j+1)
Knut's avatar
Knut committed
500 501 502
                  end do
               case (SOMMERFELD)
                  do i = sfi(n),sli(n)
Knut's avatar
Knut committed
503
                     if ( az(i,j+1) .ne. 0 ) then
Knut's avatar
Knut committed
504 505 506 507 508 509
                     cfl = sqrt(g*_HALF_*(D(i,j)+D(i,j+1)))*dtm/DYV
                     z(i,j) = (                                             &
                                (_ONE_ - _TWO_*cfl*(_ONE_-theta))*z (i,j  ) &
                               +(_ONE_ + _TWO_*cfl*(_ONE_-theta))*zo(i,j+1) &
                               -(_ONE_ - _TWO_*cfl*theta        )*z (i,j+1) &
                              )/(_ONE_ + _TWO_*cfl*theta        )
Knut's avatar
Knut committed
510
                     end if
Knut's avatar
Knut committed
511
                  end do
Knut's avatar
Knut committed
512
               case (CLAMPED_ELEV,CLAMPED)
Knut's avatar
Knut committed
513
                  do i = sfi(n),sli(n)
Knut's avatar
Knut committed
514
                     z(i,j) = max(ramp*bdy_data(kl),-H(i,j)+min_depth)
Knut's avatar
Knut committed
515
                     k = k+1
Knut's avatar
Knut committed
516
                     kl = kl + 1
Knut's avatar
Knut committed
517 518 519
                  end do
               case (FLATHER_ELEV)
                  do i = sfi(n),sli(n)
Knut's avatar
Knut committed
520
                     if ( az(i,j+1) .ne. 0 ) then
Knut's avatar
Knut committed
521 522 523 524
!                    Note (KK): approximate interface depths at vel-time stage
!                               by spatial mean at last sse-time stage
                     depth = _HALF_*(D(i,j)+D(i,j+1))
!                    Note (KK): note approximation of sse at vel-time stage
Knut's avatar
Knut committed
525 526 527 528 529
                     a = _TWO_/sqrt(g*depth)*(V(i,j)-ramp*bdy_data_v(kl)*depth)
                     else
                     a = _ZERO_
                     end if
                     z(i,j) = max(ramp*bdy_data(kl)-a,-H(i,j)+min_depth)
Knut's avatar
Knut committed
530
                     k = k+1
Knut's avatar
Knut committed
531
                     kl = kl + 1
Knut's avatar
Knut committed
532 533
                  end do
            end select
Knut's avatar
Knut committed
534 535
         end do

536 537 538 539 540 541 542 543
         if (bdy2d_sponge_size .gt. 0) then
            l = 0
            do n = 1,NWB
               l = l+1
               select case (bdy_2d_type(l))
                  case (CONSTANT,CLAMPED_ELEV,CLAMPED)
                     i = wi(n)
                     do j = wfj(n),wlj(n)
544 545 546 547 548 549 550 551 552
                        if (az(i,j) .eq. 2) then
                           do ii=1,bdy2d_sponge_size
                              if (az(i+ii,j) .eq. 1) then
                                 z(i+ii,j) = sp(ii)*z(i,j)+(_ONE_-sp(ii))*z(i+ii,j)
                              else
                                 exit
                              end if
                           end do
                        end if
553 554 555 556 557 558 559 560 561
                     end do
               end select
            end do
            do n = 1,NNB
               l = l+1
               select case (bdy_2d_type(l))
                  case (CONSTANT,CLAMPED_ELEV,CLAMPED)
                     j = nj(n)
                     do i = nfi(n),nli(n)
562 563 564 565 566 567 568 569 570
                        if (az(i,j) .eq. 2) then
                           do jj=1,bdy2d_sponge_size
                              if (az(i,j-jj) .eq. 1) then
                                 z(i,j-jj) = sp(jj)*z(i,j)+(_ONE_-sp(jj))*z(i,j-jj)
                              else
                                 exit
                              end if
                           end do
                        end if
571 572 573 574 575 576 577 578 579
                     end do
               end select
            end do
            do n = 1,NEB
               l = l+1
               select case (bdy_2d_type(l))
                  case (CONSTANT,CLAMPED_ELEV,CLAMPED)
                     i = ei(n)
                     do j = efj(n),elj(n)
580 581 582 583 584 585 586 587 588
                        if (az(i,j) .eq. 2) then
                           do ii=1,bdy2d_sponge_size
                              if (az(i-ii,j) .eq. 1) then
                                 z(i-ii,j) = sp(ii)*z(i,j)+(_ONE_-sp(ii))*z(i-ii,j)
                              else
                                 exit
                              end if
                           end do
                        end if
589 590 591 592 593 594 595 596 597
                     end do
               end select
            end do
            do n = 1,NSB
               l = l+1
               select case (bdy_2d_type(l))
                  case (CONSTANT,CLAMPED_ELEV,CLAMPED)
                     j = sj(n)
                     do i = sfi(n),sli(n)
598 599 600 601 602 603 604 605 606
                        if (az(i,j) .eq. 2) then
                           do jj=1,bdy2d_sponge_size
                              if (az(i,j+jj) .eq. 1) then
                                 z(i,j+jj) = sp(jj)*z(i,j)+(_ONE_-sp(jj))*z(i,j+jj)
                              else
                                 exit
                              end if
                           end do
                        end if
607 608 609 610 611 612
                     end do
               end select
            end do
         end if


Knut's avatar
Knut committed
613 614 615 616 617 618
      case (U_TAG)

         l = 0
         do n = 1,NWB
            l = l+1
            k = bdy_index(l)
Knut's avatar
Knut committed
619
            kl = bdy_index_l(l)
Knut's avatar
Knut committed
620 621 622 623
            i = wi(n)
            select case (bdy_2d_type(l))
               case (FLATHER_VEL)
                  do j = wfj(n),wlj(n)
Knut's avatar
Knut committed
624
                     if ( az(i+1,j) .ne. 0 ) then
Knut's avatar
Knut committed
625 626 627 628
!                    Note (KK): approximate interface depths at vel-time stage
!                               by spatial mean at last sse-time stage
                     depth = _HALF_*(D(i,j)+D(i+1,j))
!                    Note (KK): note approximation of sse at vel-time stage
Knut's avatar
Knut committed
629 630
                     U(i,j) = ramp*bdy_data_u(kl)*depth &
                              - _HALF_*sqrt(g*depth)*(z(i,j)-ramp*bdy_data(kl))
Knut's avatar
Knut committed
631
                     end if
Knut's avatar
Knut committed
632
                     k = k+1
Knut's avatar
Knut committed
633
                     kl = kl + 1
Knut's avatar
Knut committed
634 635 636
                  end do
               case (CLAMPED_VEL,CLAMPED)
                  do j = wfj(n),wlj(n)
Knut's avatar
Knut committed
637
                     if ( az(i+1,j) .ne. 0 ) then
Knut's avatar
Knut committed
638 639 640
!                    Note (KK): approximate interface depths at vel-time stage
!                               by spatial mean at last sse-time stage
                     depth = _HALF_*(D(i,j)+D(i+1,j))
Knut's avatar
Knut committed
641
                     U(i,j) = ramp*bdy_data_u(kl)*depth
Knut's avatar
Knut committed
642
                     end if
Knut's avatar
Knut committed
643
                     k = k+1
Knut's avatar
Knut committed
644
                     kl = kl + 1
Knut's avatar
Knut committed
645 646 647 648 649 650 651
                  end do
            end select
         end do
         l = l + NNB
         do n = 1,NEB
            l = l+1
            k = bdy_index(l)
Knut's avatar
Knut committed
652
            kl = bdy_index_l(l)
Knut's avatar
Knut committed
653 654 655 656
            i = ei(n)
            select case (bdy_2d_type(l))
               case (FLATHER_VEL)
                  do j = efj(n),elj(n)
Knut's avatar
Knut committed
657
                     if ( az(i-1,j) .ne. 0 ) then
Knut's avatar
Knut committed
658 659 660 661
!                    Note (KK): approximate interface depths at vel-time stage
!                               by spatial mean at last sse-time stage
                     depth = _HALF_*(D(i-1,j)+D(i,j))
!                    Note (KK): note approximation of sse at vel-time stage
Knut's avatar
Knut committed
662 663
                     U(i-1,j) = ramp*bdy_data_u(kl)*depth &
                                + _HALF_*sqrt(g*depth)*(z(i,j)-ramp*bdy_data(kl))
Knut's avatar
Knut committed
664
                     end if
Knut's avatar
Knut committed
665
                     k = k+1
Knut's avatar
Knut committed
666
                     kl = kl + 1
Knut's avatar
Knut committed
667 668 669
                  end do
               case (CLAMPED_VEL,CLAMPED)
                  do j = efj(n),elj(n)
Knut's avatar
Knut committed
670
                     if ( az(i-1,j) .ne. 0 ) then
Knut's avatar
Knut committed
671 672 673
!                    Note (KK): approximate interface depths at vel-time stage
!                               by spatial mean at last sse-time stage
                     depth = _HALF_*(D(i-1,j)+D(i,j))
Knut's avatar
Knut committed
674
                     U(i-1,j) = ramp*bdy_data_u(kl)*depth
Knut's avatar
Knut committed
675
                     end if
Knut's avatar
Knut committed
676
                     k = k+1
Knut's avatar
Knut committed
677
                     kl = kl + 1
Knut's avatar
Knut committed
678 679 680 681 682 683 684 685 686 687
                  end do
            end select
         end do

      case (V_TAG)

         l = NWB
         do n = 1,NNB
            l = l+1
            k = bdy_index(l)
Knut's avatar
Knut committed
688
            kl = bdy_index_l(l)
Knut's avatar
Knut committed
689 690 691 692
            j = nj(n)
            select case (bdy_2d_type(l))
               case (FLATHER_VEL)
                  do i = nfi(n),nli(n)
Knut's avatar
Knut committed
693
                     if ( az(i,j-1) .ne. 0 ) then
Knut's avatar
Knut committed
694 695 696 697
!                    Note (KK): approximate interface depths at vel-time stage
!                               by spatial mean at last sse-time stage
                     depth = _HALF_*(D(i,j-1)+D(i,j))
!                    Note (KK): note approximation of sse at vel-time stage
Knut's avatar
Knut committed
698 699
                     V(i,j-1) = ramp*bdy_data_v(kl)*depth &
                                + _HALF_*sqrt(g*depth)*(z(i,j)-ramp*bdy_data(kl))
Knut's avatar
Knut committed
700
                     end if
Knut's avatar
Knut committed
701
                     k = k+1
Knut's avatar
Knut committed
702
                     kl = kl + 1
Knut's avatar
Knut committed
703 704 705
                  end do
               case (CLAMPED_VEL,CLAMPED)
                  do i = nfi(n),nli(n)
Knut's avatar
Knut committed
706
                     if ( az(i,j-1) .ne. 0 ) then
Knut's avatar
Knut committed
707 708 709
!                    Note (KK): approximate interface depths at vel-time stage
!                               by spatial mean at last sse-time stage
                     depth = _HALF_*(D(i,j-1)+D(i,j))
Knut's avatar
Knut committed
710
                     V(i,j-1) = ramp*bdy_data_v(kl)*depth
Knut's avatar
Knut committed
711
                     end if
Knut's avatar
Knut committed
712
                     k = k+1
Knut's avatar
Knut committed
713
                     kl = kl + 1
Knut's avatar
Knut committed
714 715 716 717 718 719 720
                  end do
            end select
         end do
         l = l + NEB
         do n = 1,NSB
            l = l+1
            k = bdy_index(l)
Knut's avatar
Knut committed
721
            kl = bdy_index_l(l)
Knut's avatar
Knut committed
722
            j = sj(n)
Knut's avatar
Knut committed
723 724 725
            select case (bdy_2d_type(l))
               case (FLATHER_VEL)
                  do i = sfi(n),sli(n)
Knut's avatar
Knut committed
726
                     if ( az(i,j+1) .ne. 0 ) then
Knut's avatar
Knut committed
727 728 729 730
!                    Note (KK): approximate interface depths at vel-time stage
!                               by spatial mean at last sse-time stage
                     depth = _HALF_*(D(i,j)+D(i,j+1))
!                    Note (KK): note approximation of sse at vel-time stage
Knut's avatar
Knut committed
731 732
                     V(i,j) = ramp*bdy_data_v(kl)*depth &
                              - _HALF_*sqrt(g*depth)*(z(i,j)-ramp*bdy_data(kl))
Knut's avatar
Knut committed
733
                     end if
Knut's avatar
Knut committed
734
                     k = k+1
Knut's avatar
Knut committed
735
                     kl = kl + 1
Knut's avatar
Knut committed
736
                  end do
Knut's avatar
Knut committed
737
               case (CLAMPED_VEL,CLAMPED)
Knut's avatar
Knut committed
738
                  do i = sfi(n),sli(n)
Knut's avatar
Knut committed
739
                     if ( az(i,j+1) .ne. 0 ) then
Knut's avatar
Knut committed
740 741 742
!                    Note (KK): approximate interface depths at vel-time stage
!                               by spatial mean at last sse-time stage
                     depth = _HALF_*(D(i,j)+D(i,j+1))
Knut's avatar
Knut committed
743
                     V(i,j) = ramp*bdy_data_v(kl)*depth
Knut's avatar
Knut committed
744
                     end if
Knut's avatar
Knut committed
745
                     k = k+1
Knut's avatar
Knut committed
746
                     kl = kl + 1
Knut's avatar
Knut committed
747 748
                  end do
            end select
Knut's avatar
Knut committed
749 750 751
         end do

   end select
752 753 754 755 756 757 758 759

#ifdef DEBUG
   write(debug,*) 'leaving do_bdy_2d()'
   write(debug,*)
#endif
   return
   end subroutine do_bdy_2d
!EOC
Knut's avatar
Knut committed
760 761 762 763 764 765 766 767 768 769 770 771 772 773 774 775 776 777 778 779 780 781 782
!-----------------------------------------------------------------------
!BOP
!
! !IROUTINE:  LOGICAL function bdy2d_active -
!
! !INTERFACE:
   logical function bdy2d_active(type_2d)
!
! !DESCRIPTION:
!
! !USES:
   IMPLICIT NONE
!
! !INPUT PARAMETERS:
   integer,intent(in)  :: type_2d
!
! !REVISION HISTORY:
!  Original author(s): Knut Klingbeil
!EOP
!-----------------------------------------------------------------------
!BOC

   select case (type_2d)
783 784
      case (CONSTANT)
         bdy2d_active = .false.
Knut's avatar
Knut committed
785 786 787 788 789 790 791 792 793 794 795 796 797 798
      case (CLAMPED)
         bdy2d_active = .true.
      case (ZERO_GRADIENT)
         bdy2d_active = .false.
      case (SOMMERFELD)
         bdy2d_active = .false.
      case (CLAMPED_ELEV)
         bdy2d_active = .true.
      case (FLATHER_ELEV)
         bdy2d_active = .true.
      case (FLATHER_VEL)
         bdy2d_active = .true.
      case (CLAMPED_VEL)
         bdy2d_active = .true.
799 800
      case default
         bdy2d_active = .false.
Knut's avatar
Knut committed
801 802 803 804 805 806 807 808 809 810 811 812 813 814 815 816 817 818 819 820 821 822 823 824 825 826 827 828
   end select

   return
   end function bdy2d_active
!EOC
!-----------------------------------------------------------------------
!BOP
!
! !IROUTINE:  LOGICAL function bdy2d_need_elev -
!
! !INTERFACE:
   logical function bdy2d_need_elev(type_2d)
!
! !DESCRIPTION:
!
! !USES:
   IMPLICIT NONE
!
! !INPUT PARAMETERS:
   integer,intent(in)  :: type_2d
!
! !REVISION HISTORY:
!  Original author(s): Knut Klingbeil
!EOP
!-----------------------------------------------------------------------
!BOC

   select case (type_2d)
829 830
      case (CONSTANT)
         bdy2d_need_elev = .false.
Knut's avatar
Knut committed
831 832 833 834 835 836 837 838 839 840 841 842 843 844
      case (CLAMPED)
         bdy2d_need_elev = .true.
      case (ZERO_GRADIENT)
         bdy2d_need_elev = .false.
      case (SOMMERFELD)
         bdy2d_need_elev = .false.
      case (CLAMPED_ELEV)
         bdy2d_need_elev = .true.
      case (FLATHER_ELEV)
         bdy2d_need_elev = .true.
      case (FLATHER_VEL)
         bdy2d_need_elev = .true.
      case (CLAMPED_VEL)
         bdy2d_need_elev = .false.
845 846
      case default
         bdy2d_need_elev = .false.
Knut's avatar
Knut committed
847 848 849 850 851 852 853 854 855 856 857 858 859 860 861 862 863 864 865 866 867 868 869 870 871 872 873 874
   end select

   return
   end function bdy2d_need_elev
!EOC
!-----------------------------------------------------------------------
!BOP
!
! !IROUTINE:  LOGICAL function bdy2d_need_vel -
!
! !INTERFACE:
   logical function bdy2d_need_vel(type_2d)
!
! !DESCRIPTION:
!
! !USES:
   IMPLICIT NONE
!
! !INPUT PARAMETERS:
   integer,intent(in)  :: type_2d
!
! !REVISION HISTORY:
!  Original author(s): Knut Klingbeil
!EOP
!-----------------------------------------------------------------------
!BOC

   select case (type_2d)
875 876
      case (CONSTANT)
         bdy2d_need_vel = .false.
Knut's avatar
Knut committed
877 878 879 880 881 882 883 884 885 886 887 888 889 890
      case (CLAMPED)
         bdy2d_need_vel = .true.
      case (ZERO_GRADIENT)
         bdy2d_need_vel = .false.
      case (SOMMERFELD)
         bdy2d_need_vel = .false.
      case (CLAMPED_ELEV)
         bdy2d_need_vel = .false.
      case (FLATHER_ELEV)
         bdy2d_need_vel = .true.
      case (FLATHER_VEL)
         bdy2d_need_vel = .true.
      case (CLAMPED_VEL)
         bdy2d_need_vel = .true.
891 892
      case default
         bdy2d_need_vel = .false.
Knut's avatar
Knut committed
893 894 895 896 897
   end select

   return
   end function bdy2d_need_vel
!EOC
898 899 900 901 902 903 904
!-----------------------------------------------------------------------

   end module bdy_2d

!-----------------------------------------------------------------------
! Copyright (C) 2001 - Karsten Bolding and Hans Burchard               !
!-----------------------------------------------------------------------