bottom_friction.F90 6.34 KB
Newer Older
kbk's avatar
kbk committed
1
!$Id: bottom_friction.F90,v 1.7 2006-03-01 15:54:07 kbk Exp $
gotm's avatar
gotm committed
2 3 4 5
#include "cppdefs.h"
!-----------------------------------------------------------------------
!BOP
!
6
! !IROUTINE: bottom_friction - calculates the 2D bottom friction.
gotm's avatar
gotm committed
7 8 9 10 11 12
!
! !INTERFACE:
   subroutine bottom_friction(runtype)
!
! !DESCRIPTION:
!
hb's avatar
hb committed
13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34
! In this routine the bottom friction for the external (vertically integrated)
! mode is calculated. This is done separately for the $U$-equation in the
! U-points and for the $V$-equation in the V-points. 
! The drag coefficient $R$ for the external mode is given in eq.\
! (\ref{bottom_vert}) on page \pageref{bottom_vert}.
! For {\tt runtype=1} (only vertically integrated calculations), the
! bottom roughness length is depending on the bed friction 
! velocity $u_*^b$ and the molecular viscosity $\nu$:
!
! \begin{equation}\label{Defz0b}
! z_0^b = 0.1 \frac{\nu}{u_*^b} + \left(z^b_0\right)_{\min},
! \end{equation}
! 
! see e.g.\ \cite{KAGAN95}, i.e.\ the given roughness may be increased
! by viscous effects.
! After this, the drag coefficient is multiplied by the absolute value of the
! local velocity, which is alculated by dividing the local transports by the
! local water depths and by properly interpolating these velocities 
! to the U- and V-points. The resulting fields are {\tt ru}, representing
! $R\sqrt{u^2+v^2}$ on the U-points and {\tt rv}, representing
! this quantity on the V-points.
!
gotm's avatar
gotm committed
35 36
! !USES:
   use parameters, only: kappa,avmmol
kbk's avatar
kbk committed
37
   use domain, only: imin,imax,jmin,jmax,au,av,min_depth
gotm's avatar
gotm committed
38 39 40 41
   use variables_2d
   IMPLICIT NONE
!
! !INPUT PARAMETERS:
kbk's avatar
kbk committed
42
   integer, intent(in)                 :: runtype
gotm's avatar
gotm committed
43 44 45 46 47
!
! !INPUT/OUTPUT PARAMETERS:
!
! !OUTPUT PARAMETERS:
!
kbk's avatar
kbk committed
48 49 50
! !REVISION HISTORY:
!  Original author(s): Karsten Bolding & Hans Burchard
!
gotm's avatar
gotm committed
51
!  !LOCAL VARIABLES:
kbk's avatar
kbk committed
52 53 54
   integer                   :: i,j
   REALTYPE                  :: uloc(E2DFIELD),vloc(E2DFIELD)
   REALTYPE                  :: HH(E2DFIELD),fricvel(E2DFIELD)
gotm's avatar
gotm committed
55 56 57 58 59 60 61 62 63
!EOP
!-----------------------------------------------------------------------
!BOC
#ifdef DEBUG
   integer, save :: Ncall = 0
   Ncall = Ncall+1
   write(debug,*) 'bottom_friction() # ',Ncall
#endif

kbk's avatar
kbk committed
64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105
#ifdef DEBUG
   if(Ncall .eq. 1) then
      STDERR 'bottom_friction(): checking for 0 depth values'
      do j=jmin,jmax
         do i=imin,imax
            if (av(i,j) .eq. 1) then
               if(DU(i,j) .eq. _ZERO_) then
                  STDERR 'DU=0 uv_advect ',i,j,i,j,av(i,j)
               end if
               if(DU(i-1,j) .eq. _ZERO_) then
                  STDERR 'DU=0 uv_advect ',i-1,j,i,j,av(i,j)
               end if
               if(DU(i,j+1) .eq. _ZERO_) then
                  STDERR 'DU=0 uv_advect ',i,j+1,i,j,av(i,j)
               end if
               if(DU(i-1,j+1) .eq. _ZERO_) then
                  STDERR 'DU=0 uv_advect ',i-1,j+1,i,j,av(i,j)
               end if
            end if
         end do
      end do
      do j=jmin,jmax
         do i=imin,imax
            if (au(i,j) .eq. 1) then
               if(DV(i,j) .eq. _ZERO_) then
                  STDERR 'DV=0 uv_advect ',i,j,i,j,au(i,j)
               end if
               if(DV(i+1,j) .eq. _ZERO_) then
                  STDERR 'DV=0 uv_advect ',i+1,j,i,j,au(i,j)
               end if
               if(DV(i,j-1) .eq. _ZERO_) then
                  STDERR 'DV=0 uv_advect ',i,j-1,i,j,au(i,j)
               end if
               if(DV(i+1,j-1) .eq. _ZERO_) then
                  STDERR 'DV=0 uv_advect ',i+1,j-1,i,j,au(i,j)
               end if
            end if
         end do
      end do
   end if
#endif

gotm's avatar
gotm committed
106 107
   do j=jmin,jmax
      do i=imin,imax
108
         if (au(i,j) .gt. 0) then
kbk's avatar
kbk committed
109 110 111
            vloc(i,j)=0.25* ( V(i  ,j  )/DV(i  ,j  )   &
                             +V(i+1,j  )/DV(i+1,j  )   &
                             +V(i  ,j-1)/DV(i  ,j-1)   &
112 113 114 115
                             +V(i+1,j-1)/DV(i+1,j-1) )
         else
            vloc(i,j) = _ZERO_
         end if
gotm's avatar
gotm committed
116 117 118 119 120 121 122
      end do
   end do

!  The x-direction

#ifndef DEBUG
   where (au .gt. 0)
123
      uloc=U/DU
gotm's avatar
gotm committed
124 125 126 127
      HH=max(min_depth,DU)
      ruu=(kappa/log((zub+0.5*HH)/zub))**2
   end where
#else
128
   uloc=U/DU
gotm's avatar
gotm committed
129 130 131 132 133 134 135 136 137 138 139 140 141 142 143 144 145 146 147 148 149 150 151 152 153 154 155 156 157 158 159 160 161 162
   HH=max(min_depth,DU)
   ruu=(zub+0.5*HH)/zub

   do j=jmin,jmax
      do i=imin,imax
         if (ruu(i,j) .le. _ONE_) then
            STDERR 'bottom_friction friction coefficient infinite.'
            STDERR 'i = ',i,' j = ',j
            STDERR 'min_depth = ',min_depth,' DU = ',DU(i,j)
            stop
         end if
      end do
   end do
   where (au .gt. 0)
      ruu=(kappa/log(ruu))**2
   end where
#endif

   if (runtype .eq. 1) then
      where (au .gt. 0)
         fricvel=sqrt(ruu*(uloc**2+vloc**2))
         zub=min(HH,zub0+0.1*avmmol/max(avmmol,fricvel))
         ruu=(zub+0.5*HH)/zub
         ruu=(kappa/log(ruu))**2
      end where
   end if

   where (au .gt. 0)
      ru=ruu*sqrt(uloc**2+vloc**2)
   end where

!  The y-direction
   do j=jmin,jmax
      do i=imin,imax
163
         if (av(i,j) .gt. 0) then
kbk's avatar
kbk committed
164 165 166
            uloc(i,j)=0.25* ( U(i  ,j  )/DU(i  ,j  )   &
                             +U(i-1,j  )/DU(i-1,j  )   &
                             +U(i  ,j+1)/DU(i  ,j+1)   &
167 168 169 170
                             +U(i-1,j+1)/DU(i-1,j+1) )
         else
            uloc(i,j) = _ZERO_
         end if
gotm's avatar
gotm committed
171 172 173 174 175
      end do
   end do

#ifndef DEBUG
   where (av .gt. 0)
176
      vloc=V/DV
gotm's avatar
gotm committed
177 178 179 180
      HH=max(min_depth,DV)
      rvv=(kappa/log((zvb+0.5*HH)/zvb))**2
   end where
#else
181
   vloc=V/DV
gotm's avatar
gotm committed
182 183 184 185 186 187 188 189 190 191 192 193 194 195 196 197 198 199 200
   HH=max(min_depth,DV)
   rvv=(zvb+0.5*HH)/zvb

   do j=jmin,jmax
      do i=imin,imax
         if (rvv(i,j) .le. _ONE_) then
            STDERR 'bottom_friction friction coefficient infinite.'
            STDERR 'i = ',i,' j = ',j
            STDERR 'min_depth = ',min_depth,' DV = ',DU(i,j)
            stop
         end if
      end do
   end do

   where (av .gt. 0)
      rvv=(kappa/log(rvv))**2
   end where
#endif

201
   if (runtype .eq. 1) then
gotm's avatar
gotm committed
202 203
      where (av .gt. 0)
         fricvel=sqrt(rvv*(uloc**2+vloc**2))
204
         zvb=min(HH,zvb0+0.1*avmmol/max(avmmol,fricvel))
gotm's avatar
gotm committed
205 206 207
         rvv=(zvb+0.5*HH)/zvb
         rvv=(kappa/log(rvv))**2
      end where
208
   end if
gotm's avatar
gotm committed
209 210 211 212 213 214 215 216 217 218 219 220 221 222 223 224

   where (av .gt. 0)
      rv=rvv*sqrt(uloc**2+vloc**2)
   end where

#ifdef DEBUG
   write(debug,*) 'Leaving bottom_friction()'
   write(debug,*)
#endif
   return
   end subroutine bottom_friction
!EOC

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