sso_z0_frictionn.F90

Go to the documentation of this file.

!SFX_LIC Copyright 1994-2014 CNRS, Meteo-France and Universite Paul Sabatier
!SFX_LIC This is part of the SURFEX software governed by the CeCILL-C licence
!SFX_LIC version 1. See LICENSE, CeCILL-C_V1-en.txt and CeCILL-C_V1-fr.txt  
!SFX_LIC for details. version 1.
!     ################################################################################
SUBROUTINE sso_z0_friction_n (USS, &
                              psea,puref,prhoa,pu,pv,ppew_a_coef,ppew_b_coef,psfu,psfv)
!     ################################################################################
!
!!****  *SSO_Z0_FRICTION_n * - Computes subgrid-scale orography friction
!                                  according to several options:
!                                CROUGH='Z01D' : orographic roughness length
!                                CROUGH='Z04D' : orographic roughness length
!                                                variable with wind direction
!!
!!    PURPOSE
!!    -------
!
!!**  METHOD
!!    ------
!!
!!    REFERENCE
!!    ---------
!!      
!!
!!    AUTHOR
!!    ------
!!     V. Masson 
!!
!!    MODIFICATIONS
!!    -------------
!!      Original    05/2010
!!      E. Martin   01/2012 Correction masque (compatibilité XUNDEF)
!!      B. Decharme 09/2012 new wind implicitation and sea fraction
!!      B. Decharme 06/2013 CIMPLICIT_WIND in MODD_REPROD_OPER
!!      J. Escobar  05/2014 for bug with ifort/10, replace WHERE by IF
!!      J. Escobar  06/2015 bug with gfortran ZZ0EFF to small, change with > XSURF_EPSILON
!----------------------------------------------------------------
!
!
USE modd_surf_atm_sso_n, ONLY : surf_atm_sso_t
!
USE modd_reprod_oper, ONLY : cimplicit_wind
!
USE modd_surf_par,         ONLY : xundef, xsurf_epsilon
USE modd_csts,             ONLY : xkarman, xpi
!
USE yomhook   ,ONLY : lhook,   dr_hook
USE parkind1  ,ONLY : jprb
!
IMPLICIT NONE
!
!*      0.1    declarations of arguments
!
!
TYPE(surf_atm_sso_t), INTENT(INOUT) :: uss
!
REAL, DIMENSION(:), INTENT(IN)    :: psea      ! Sea fraction                          (-)
REAL, DIMENSION(:), INTENT(IN)    :: puref     ! Wind forcing height                   (m)
REAL, DIMENSION(:), INTENT(IN)    :: prhoa     ! air density                           (kg/m3)
REAL, DIMENSION(:), INTENT(IN)    :: pu        ! zonal wind                            (m/s)
REAL, DIMENSION(:), INTENT(IN)    :: pv        ! meridian wind                         (m/s)
REAL, DIMENSION(:), INTENT(IN)    :: ppew_a_coef! implicit coefficients                (m2s/kg)
REAL, DIMENSION(:), INTENT(IN)    :: ppew_b_coef! needed if HCOUPLING='I'              (m/s)
REAL, DIMENSION(:), INTENT(INOUT) :: psfu      ! zonal momentum flux                   (Pa)
REAL, DIMENSION(:), INTENT(INOUT) :: psfv      ! meridian momentum flux                (Pa)
!
!*      0.2    declarations of local variables
!
REAL, DIMENSION(SIZE(PU))    :: zwind   ! wind strength (m/s)
REAL, DIMENSION(SIZE(PU))    :: zwork   ! work array
REAL, DIMENSION(SIZE(PU))    :: zdir    ! wind direction (rad., clockwise)
REAL, DIMENSION(SIZE(PU))    :: zalfa   ! angle between z0eff J axis and wind direction (rad., clockwise) 
REAL, DIMENSION(SIZE(PU))    :: zcos2, zsin2
REAL, DIMENSION(SIZE(PU))    :: zz0eff  ! Momentum Roughness length
REAL, DIMENSION(SIZE(PU))    :: zcd     ! drag coefficient
REAL, DIMENSION(SIZE(PU))    :: zustar2 ! square of friction velocity
REAL, DIMENSION(SIZE(PU))    :: zsso_sfu! zonal orographic momentum flux
REAL, DIMENSION(SIZE(PU))    :: zsso_sfv! meridian orographic momentum flux
LOGICAL, DIMENSION(SIZE(PU)) :: gmask   ! mask where SSO exists
INTEGER                      :: ii
REAL(KIND=JPRB) :: zhook_handle

!-------------------------------------------------------------------------------------
!
IF (lhook) CALL dr_hook('SSO_Z0_FRICTION_N',0,zhook_handle)
!
zwork(:) = xundef
!
!*      1.     roughness length formalism
!              --------------------------
!
!* wind strength
!
zwind(:) = sqrt(pu(:)**2+pv(:)**2)
!
!* wind direction
!
zdir(:) = 0.
WHERE (zwind(:)>0.)  zdir(:)=atan2(pu(:),pv(:))
!
!* default value
!
gmask(:)=(psea(:)/=1..AND. uss%XZ0REL(:)/=0.)
zz0eff(:) = xundef
!
!*      2.     Constant orographic roughness length
!              ------------------------------------
!
IF (uss%CROUGH=="Z01D") zz0eff(:) = uss%XZ0REL(:)
!
!*      3.     Directionnal roughness length
!              -----------------------------
!
IF (uss%CROUGH=="Z04D") THEN
   DO ii=1,SIZE(gmask)
      IF (gmask(ii)) THEN
         !
         zalfa(ii) = zdir(ii) - uss%XZ0EFFJPDIR(ii) * xpi/180.
         !
         IF    (zalfa(ii)<=-xpi) THEN
            zalfa(ii) = zalfa(ii) + 2.*xpi
         ELSEIF(zalfa(ii)>  xpi) THEN
            zalfa(ii) = zalfa(ii) - 2.*xpi
         END IF
         !
         IF (zalfa(ii)>=-xpi.AND.zalfa(ii)<=xpi) THEN
            !
            zsin2(ii) = sin(zalfa(ii))**2
            zcos2(ii) = cos(zalfa(ii))**2
            !
            IF (zalfa(ii)<0.) THEN
               zz0eff(ii)=uss%XZ0EFFIM(ii)*zsin2(ii)
            ELSE
               zz0eff(ii)=uss%XZ0EFFIP(ii)*zsin2(ii)
            END IF
            !
            IF (zalfa(ii)>=-xpi/2. .AND. zalfa(ii)<xpi/2.) THEN
               zz0eff(ii) = zz0eff(ii) + uss%XZ0EFFJP(ii)*zcos2(ii)
            ELSE
               zz0eff(ii) = zz0eff(ii) + uss%XZ0EFFJM(ii)*zcos2(ii)
            END IF
            !
         END IF
         !    
      END IF
   END DO
ENDIF
!
!*      4.     Friction coefficient
!              --------------------
!
zcd(:) = 0.
zustar2(:) = 0.
!
gmask(:)=(gmask(:).AND.zz0eff(:)>xsurf_epsilon)
!
DO ii=1,SIZE(gmask)
  !
  IF (gmask(ii)) THEN
    !
    !* sets a limit to roughness length
    zz0eff(ii) = min(zz0eff(ii),puref(ii)/uss%XFRACZ0)
    !
    ! neutral case
    zcd(ii) = (xkarman/log(puref(ii)/zz0eff(ii)))**2
  END IF
  !
END DO
!
!*      5.     Friction due to orography
!              -------------------------
!
! Modify flux-form implicit coupling coefficients:
!
IF(cimplicit_wind=='OLD')THEN
! old implicitation
  zustar2(:) =  zcd(:)*zwind(:)*ppew_b_coef(:)   &
             / (1.0-prhoa(:)*zcd(:)*zwind(:)*ppew_a_coef(:))
ELSE
! new implicitation
  zustar2(:) = (zcd(:)*zwind(:)*(2.*ppew_b_coef(:)-zwind(:))   )   &
             / (1.0-2.0*prhoa(:)*zcd(:)*zwind(:)*ppew_a_coef(:))
ENDIF
!
WHERE (gmask(:))
!
  zwork(:) = prhoa(:)*ppew_a_coef(:)*zustar2(:) + ppew_b_coef(:)
  zwork(:) = max(zwork(:),0.)
!
  WHERE(ppew_a_coef(:)/= 0.)
    zustar2(:) = max( ( zwork(:) - ppew_b_coef(:) ) / (prhoa(:)*ppew_a_coef(:)), 0.)
  ENDWHERE
!
END WHERE
!
!*      6.     Projection of friction on wind components
!              -----------------------------------------
!
zsso_sfu(:) = 0.
zsso_sfv(:) = 0.
WHERE (zwind(:)>0.)
  zsso_sfu(:) = - pu(:)/zwind(:) * zustar2(:) * prhoa(:)
  zsso_sfv(:) = - pv(:)/zwind(:) * zustar2(:) * prhoa(:)
END WHERE
!
!*      7.     Adds orographic friction to other sources of friction
!              -----------------------------------------------------
!
psfu(:) = psfu(:) + zsso_sfu(:) * (1.0-psea(:))
psfv(:) = psfv(:) + zsso_sfv(:) * (1.0-psea(:))
!
IF (lhook) CALL dr_hook('SSO_Z0_FRICTION_N',1,zhook_handle)
!
!-------------------------------------------------------------------------------------
!
END SUBROUTINE sso_z0_friction_n