brough.f90

Go to the documentation of this file.

!/===========================================================================/
! Copyright (c) 2007, The University of Massachusetts Dartmouth 
! Produced at the School of Marine Science & Technology 
! Marine Ecosystem Dynamics Modeling group
! All rights reserved.
!
! FVCOM has been developed by the joint UMASSD-WHOI research team. For 
! details of authorship and attribution of credit please see the FVCOM
! technical manual or contact the MEDM group.
!
! 
! This file is part of FVCOM. For details, see http://fvcom.smast.umassd.edu 
! The full copyright notice is contained in the file COPYRIGHT located in the 
! root directory of the FVCOM code. This original header must be maintained
! in all distributed versions.
!
! THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS" 
! AND ANY EXPRESS OR  IMPLIED WARRANTIES, INCLUDING,  BUT NOT  LIMITED TO,
! THE IMPLIED WARRANTIES OF MERCHANTABILITY AND  FITNESS FOR A PARTICULAR
! PURPOSE ARE DISCLAIMED.  
!
!/---------------------------------------------------------------------------/
! CVS VERSION INFORMATION
! $Id$
! $Name$
! $Revision$
!/===========================================================================/

!==============================================================================|
!  Calculate Bottom Drag Coefficient based on Bottom Roughness                 !
!   note:                                                                      !
!   when the log function derived from the constant stress log-viscous         !
!   layer is applied to an estuary, if the value of z0 is close to             !
!   (zz(kbm1)-z(kb)*dt1, drag coefficient "cbc" could become a huge            !
!   number due to near-zero value of alog function. In our application         !
!   we simply cutoff at cbc=0.005. One could adjust this cutoff value          !
!   based on observations or his or her experiences.                           !   
!   CALCULATES:   WUBOT(N), WVBOT(N) : BOTTOM SHEAR STRESSES                   !
!==============================================================================|

SUBROUTINE bottom_roughness

  !==============================================================================!
  USE all_vars
  USE mod_utils
  USE mod_wd
  USE mod_par

  IMPLICIT NONE
  INTEGER :: I,II
  REAL(SP), PARAMETER  :: KAPPA = .40_sp   !!VON KARMAN LENGTH SCALE
  REAL(SP), PARAMETER  :: VK2   = .160_sp  !!KAPPA SQUARED
  REAL(SP)             :: ZTEMP,BTPS,RR,U_TAUB,Z0B_GOTM,CFF

! USED IN 2D MODEL ONLY
!   REAL(SP), PARAMETER  :: CONST_CD=.0015_SP      !! CD SET CONSTANT TO THIS VALUE
   REAL(SP), PARAMETER  :: ALFA =  .166667_sp, &   !! POWER OF WATER DEPTH
                           nn   = 0.02_sp          !! FACTOR TO DIVIDE
!   REAL(SP), PARAMETER  :: CFMIN   = .0025_SP, &  !! DEEP WATER VALUE
!                           H_BREAK = 1.0_SP,    & !! 
!                           THETA   = 10._SP,   &  !!
!                           LAMB    = 0.3333333333_SP
  !==============================================================================!

  if(dbg_set(dbg_sbr)) write(ipt,*) "Start: BOTTOM_ROUGHNESS"

  !
  !  SET CONSTANTS
  !

  SELECT CASE(bottom_roughness_type) 
  !==============================================================================|
  CASE(br_orig) !USE ORIGINAL FVCOM FORM FOR BOTTOM FRICTION  |
  !==============================================================================|

     ! SET A EFFECTIVE MAXIMUM FOR CBC USING THE DEPTH
!!     WHERE (DT1 > 3.0_SP)
!!        CBC = VK2/(LOG((ZZ1(:,KBM1)-Z1(:,KB))*DT1/CC_Z0B))**2
!!     ELSEWHERE
!!        CBC = VK2/(LOG((ZZ1(:,KBM1)-Z1(:,KB))*3.0/CC_Z0B))**2
!!     END WHERE
     WHERE (dt1(1:nt) > 3.0_sp)
        cbc(1:nt) = vk2/(log((zz1(1:nt,kbm1)-z1(1:nt,kb))*dt1(1:nt)/cc_z0b(1:nt)))**2
     ELSEWHERE
        cbc(1:nt) = vk2/(log((zz1(1:nt,kbm1)-z1(1:nt,kb))*3.0/cc_z0b(1:nt)))**2
     END WHERE
     
     ! SET A MINIMUM VALUE FOR CBC
     WHERE (cbc < cbcmin)
        cbc=cbcmin
     END WHERE


  !==============================================================================|
  CASE(br_gotm) !GOTM FORMULATION FOR BOTTOM FRICTION    |
  !==============================================================================|

     !----Convert Input Z0B to GOTMS H0B
     !     H0B = Z0B/.03  
     ! DAS fixed bug to match gotm's friction.f90
     DO i=1,n
        u_taub = 0.0_sp
        DO ii=1,40       
           IF (umol <= 0.0_sp) THEN
              z0b_gotm=cc_z0b(i)   !0.03*H0B 
           ELSE
              z0b_gotm=0.1_sp*umol/max(umol,u_taub)+cc_z0b(i) !0.03*H0B
           END IF
           ztemp=(zz1(i,kbm1)-z1(i,kb))*dt1(i)
           rr=kappa/(log((z0b_gotm+ztemp)/z0b_gotm))
           u_taub = rr*sqrt( u(i,kbm1)*u(i,kbm1) + v(i,kbm1)*v(i,kbm1) )
        END DO
        cbc(i) =   rr*rr
     END DO


  CASE DEFAULT
     CALL fatal_error ("BROUGH: UNKNOWN BOTTOM_ROUGHNESS_TYPE:"&
          & ,trim(bottom_roughness_type) )
  END SELECT


  !==============================================================================|
  !  CALCULATE SHEAR STRESS ON BOTTOM  --> WUBOT/WVBOT                           |
  !==============================================================================|
  DO  i = 1, n
     btps= cbc(i)*sqrt(u(i,kbm1)**2+v(i,kbm1)**2)
     wubot(i) = -btps * u(i,kbm1)
     wvbot(i) = -btps * v(i,kbm1)
     cff=0.75_sp*dz1(i,kbm1)*d1(i)
     wubot(i)=sign(1.0_sp,wubot(i))*min(abs(wubot(i)),abs(u(i,kbm1)*cff/dti))
     wvbot(i)=sign(1.0_sp,wvbot(i))*min(abs(wvbot(i)),abs(v(i,kbm1)*cff/dti))
     !for plb case only
  
  END DO


  !==============================================================================|
  !  Calculate shear stress on nodes (x-component, y-component, magnitude) 
  !==============================================================================|
  taubm = sqrt(wubot**2 + wvbot**2) 

     CALL e2n2d(wubot,wubot_n)
     CALL e2n2d(wvbot,wvbot_n)
  taubm_n = sqrt(wubot_n**2 + wvbot_n**2) 

  if(dbg_set(dbg_sbr)) write(ipt,*) "End: BOTTOM_ROUGHNESS"

  RETURN
END SUBROUTINE bottom_roughness
!==============================================================================|