isba_snow_agr.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 isba_snow_agr(KK, PK, PEK, DMK, DK, DEK, &
                               OMEB, OMEB_LITTER, PEXNS, PEXNA, PTA, PQA,  &
                               PZREF, PUREF, PDIRCOSZW, PVMOD, PRR, PSR,   &
                               PEMIS, PALB, PUSTAR, PLES3L, PLEL3L,        &
                               PEVAP3L, PQS3L, PALB3L, PGSFCSNOW,          &
                               PZGRNDFLUX, PFLSN_COR, PEMIST, PPALPHAN )
!     ##########################################################################
!
!
!!****  *ISBA_SNOW_AGR* aggregates snow free and snow fluxes
!!
!!    PURPOSE
!!    -------
!     
!!**  METHOD
!!    ------
!
!!    EXTERNAL
!!    --------
!!
!!    IMPLICIT ARGUMENTS
!!    ------------------ 
!!
!!      
!!    REFERENCE
!!    ---------
!!
!!    Noilhan and Planton (1989)
!!      
!!    AUTHOR
!!    ------
!!      V. Masson           * Meteo-France *
!!      (following A. Boone)
!!
!!    MODIFICATIONS
!!    -------------
!!      Original    10/03/95 
!!      B. Decharme 01/2009  Floodplains 
!!      B. Decharme 01/2010  Effective surface temperature (for diag)
!!      B. Decharme 09/2012  Bug total sublimation flux: no DEK%XLESL
!!      B. Decharme 04/2013  Bug wrong radiative temperature
!!                           Sublimation diag flux
!!                           Qs for 3l or crocus (needed for coupling with atm)
!!      A. Boone    11/2014  MEB
!-------------------------------------------------------------------------------
!
!*       0.     DECLARATIONS
!               ------------
USE modd_isba_n, ONLY : isba_k_t, isba_p_t, isba_pe_t
USE modd_diag_n, ONLY : diag_t
USE modd_diag_evap_isba_n, ONLY : diag_evap_isba_t
USE modd_diag_misc_isba_n, ONLY : diag_misc_isba_t
!
USE modd_surf_par,   ONLY : xundef
!
USE yomhook   ,ONLY : lhook,   dr_hook
USE parkind1  ,ONLY : jprb
!
IMPLICIT NONE
!
!*      0.1    declarations of arguments
!              -------------------------
!
!
!* general variables
!  -----------------
!
TYPE(isba_k_t), INTENT(INOUT) :: KK
TYPE(isba_p_t), INTENT(INOUT) :: PK
TYPE(isba_pe_t), INTENT(INOUT) :: PEK
TYPE(diag_t), INTENT(INOUT) :: DK
TYPE(diag_evap_isba_t), INTENT(INOUT) :: DEK
TYPE(diag_misc_isba_t), INTENT(INOUT) :: DMK
!
LOGICAL,              INTENT(IN)  :: OMEB       ! True = patch with multi-energy balance 
!                                               ! False = patch with classical ISBA
LOGICAL, INTENT(IN)               :: OMEB_LITTER !True = litter option activated
!                                                 ! over the ground
!
!* surface and atmospheric parameters
!  ----------------------------------
!
REAL, DIMENSION(:), INTENT(IN)  :: PEXNS     ! Exner function at the surface
REAL, DIMENSION(:), INTENT(IN)  :: PEXNA     ! Exner function
REAL, DIMENSION(:), INTENT(IN)  :: PTA       ! air temperature
REAL, DIMENSION(:), INTENT(IN)  :: PQA       ! air specific humidity
REAL, DIMENSION(:), INTENT(IN)  :: PZREF     ! reference height of the first atmospheric level
REAL, DIMENSION(:), INTENT(IN)  :: PUREF     ! reference height of the wind
REAL, DIMENSION(:), INTENT(IN)  :: PDIRCOSZW ! Cosinus of the angle between the normal to the surface and the vertical
REAL, DIMENSION(:), INTENT(IN)  :: PVMOD     ! module of the horizontal wind
REAL, DIMENSION(:), INTENT(IN)  :: PRR       ! Rain rate (in kg/m2/s)
REAL, DIMENSION(:), INTENT(IN)  :: PSR       ! Snow rate (in kg/m2/s)
!
!* surface parameters
!  ------------------
!
REAL, DIMENSION(:), INTENT(IN)  :: PALB       ! albedo 
REAL, DIMENSION(:), INTENT(IN)  :: PEMIS      ! emissivity
!  'D95'     : represents aggregated (snow + flood + snow-flood-free) albedo and emissivity
!  '3-L'     : represents                    flood + snow-flood-free  albedo and emissivity
!  'MEB+3-L' : represents aggregated (snow + flood + snow-flood-free) albedo and emissivity
!
!
!* snow fractions
!  --------------
!
REAL, DIMENSION(:), INTENT(IN)  :: PPALPHAN   ! fraction of the the explicit veg.
!                                             ! canopy buried by snow
!
!* ISBA-SNOW3L variables/parameters:
!  ---------------------------------
!
! Prognostic variables:
!
REAL, DIMENSION(:),   INTENT(IN) :: PALB3L      ! Snow albedo
REAL, DIMENSION(:),   INTENT(IN) :: PQS3L       ! Surface humidity
! 
! Diagnostics:
!
REAL, DIMENSION(:), INTENT(IN)    :: PZGRNDFLUX ! snow/soil-biomass interface flux (W/m2)
REAL, DIMENSION(:), INTENT(IN)    :: PFLSN_COR  ! snow/soil-biomass correction flux (W/m2)
!
REAL, DIMENSION(:), INTENT(INOUT) :: PGSFCSNOW  ! heat flux from snow sfc to sub sfc layers (W/m2)
REAL, DIMENSION(:), INTENT(IN)    :: PLES3L     ! sublimation from ISBA-ES(3L)
REAL, DIMENSION(:), INTENT(IN)    :: PLEL3L     ! evaporation heat flux of water in the snow (W/m2)
REAL, DIMENSION(:), INTENT(INOUT) :: PEVAP3L    ! evaporation flux over snow from ISBA-ES (kg/m2/s)
!
!* diagnostic variables
!  --------------------
!
REAL, DIMENSION(:), INTENT(INOUT) :: PEMIST   ! total surface emissivity
!
!* surface fluxes
!  --------------
!
REAL, DIMENSION(:), INTENT(INOUT) :: PUSTAR   ! friction velocity
!
!*      0.2    declarations of local variables
!
REAL, DIMENSION(SIZE(PTA))       :: ZWORK
!
REAL(KIND=JPRB) :: ZHOOK_HANDLE
!
!-------------------------------------------------------------------------------
!
IF (lhook) CALL dr_hook('ISBA_SNOW_AGR',0,zhook_handle)
!
zwork(:) = 0.
!
IF(omeb)THEN
  !
  ! Snow free (ground-based snow) diagnostics: canopy and ground blended (W m-2):
  ! NOTE that the effects of snow cover *fraction* are implicitly *included* in these fluxes 
  ! so do NOT multiply by snow fraction.
  !
  dek%XRN_SN_FR   (:) = dek%XSWNET_V(:) + dek%XSWNET_G(:) + dek%XLWNET_V(:) + dek%XLWNET_G(:)
  dek%XH_SN_FR    (:) = dek%XH_CV(:) + dek%XH_GN(:)
  IF (omeb_litter) THEN
    dek%XLEG_SN_FR (:) = dek%XLELITTER (:)
    dek%XLEGI_SN_FR(:) = dek%XLELITTERI(:)
    dek%XLEG (:) = dek%XLELITTER (:)
    dek%XLEGI(:) = dek%XLELITTERI(:)
  ELSE
    dek%XLEG_SN_FR (:) = dek%XLEG (:)
    dek%XLEGI_SN_FR(:) = dek%XLEGI(:)
  ENDIF

  dek%XLEV_SN_FR (:) = dek%XLEV_CV (:)
  dek%XLETR_SN_FR(:) = dek%XLETR_CV(:) 
  ! NOTE for now, this is same as total Ustar (includes snow)   
  dek%XUSTAR_SN_FR(:) = pustar(:)        
  ! LER does not include intercepted snow sublimation
  dek%XLER_SN_FR  (:) = dek%XLEV_CV(:) - dek%XLETR_CV(:) 

  dek%XLEI_SN_FR  (:) = dek%XLEGI(:) + dek%XLEI_FLOOD(:) + dek%XLES(:) + dek%XLES_CV(:)
  ! LE includes intercepted snow sublimation
  dek%XLE_SN_FR   (:) = dek%XLEG_SN_FR(:) + dek%XLEGI_SN_FR(:) + dek%XLEV_SN_FR(:) + &
                 dek%XLES_CV(:) + dek%XLE_FLOOD(:) + dek%XLEI_FLOOD(:)
  dek%XGFLUX_SN_FR(:) = dek%XRN_SN_FR(:) - dek%XH_SN_FR(:) - dek%XLE_SN_FR(:)
  !
  pemist(:) = pemis(:)
  !
  ! Effective surface temperature (for diag): for MEB:

  zwork(:) =  ppalphan(:)*pek%XPSN(:)
  dk%XTS(:) = (1.0 - zwork(:))*pek%XTC(:) + zwork(:)*dmk%XSNOWTEMP(:,1)
  !
  ! Total heat FLUX into snow/soil/vegetation surface:
  !
  dk%XGFLUX(:) = dk%XRN(:) - dk%XH(:) - pek%XLE(:) + dmk%XHPSNOW(:) 
  !
ELSE
!
! * 2. Using an explicit snow scheme option with composite soil/veg ISBA:
!      ------------------------------------------------------------------
!
   IF(pek%TSNOW%SCHEME == '3-L' .OR. pek%TSNOW%SCHEME == 'CRO')THEN
!
!     Save fluxes from Force-Restore snow/explicit snow-free
!     portion of grid box (vegetation/soil):
!
      dek%XLEG_SN_FR  (:) = dek%XLEG (:)
      dek%XLEGI_SN_FR (:) = dek%XLEGI(:)
      dek%XLEV_SN_FR  (:) = dek%XLEV (:)
      dek%XLETR_SN_FR (:) = dek%XLETR(:)
      dek%XLER_SN_FR  (:) = dek%XLER (:)
      dek%XRN_SN_FR   (:) = dk%XRN   (:)
      dek%XH_SN_FR    (:) = dk%XH    (:)
      dek%XUSTAR_SN_FR(:) = pustar(:)      

      dek%XLE_SN_FR   (:) = pek%XLE(:)
      dek%XGFLUX_SN_FR(:) = dk%XGFLUX(:)
!  
      dek%XLEI_SN_FR  (:)= dek%XLEGI(:) + dek%XLEI_FLOOD(:) + dek%XLES(:)
!
!     Effective surface temperature (for diag):
!
      dk%XTS(:) = (1.-pek%XPSN(:))*pek%XTG(:,1)+pek%XPSN(:)*dmk%XSNOWTEMP(:,1)
!
!     Effective surface radiating temperature:
!
      dk%XALBT (:) = palb(:)*(1.-pek%XPSN(:)) + pek%XPSN(:)*palb3l(:)
      pemist(:) = pemis(:)*(1.-pek%XPSN(:)) + pek%XPSN(:)*pek%TSNOW%EMIS(:)
!  
      dk%XTSRAD(:) = ( ((1.-pek%XPSN(:))*pemis(:)*pek%XTG(:,1)**4 +   &
                               pek%XPSN(:) *pek%TSNOW%EMIS(:)*dmk%XSNOWTEMP(:,1)**4     &
                            )/pemist(:) )**(0.25)  
!
!     Calculate actual fluxes from snow-free natural
!     portion of surface: NET flux from surface is the sum of
!     fluxes from snow free and snow covered portions
!     of natural portion of grid box when *ISBA-ES* in force.
!     when NOT in use, then these fluxes equal those above.
!
      dk%XRN   (:) = (1.-pek%XPSN(:)) * dk%XRN(:) + pek%XPSN(:) * dmk%XRNSNOW(:)
      dk%XH    (:) = (1.-pek%XPSN(:)) * dk%XH (:) + pek%XPSN(:) * dmk%XHSNOW(:)
!  
      dek%XLEG (:) = (1.-pek%XPSNG(:)-kk%XFFG(:)) * dek%XLEG(:)  
      dek%XLEGI(:) = (1.-pek%XPSNG(:)-kk%XFFG(:)) * dek%XLEGI(:)  
      dek%XLEV (:) = (1.-pek%XPSNV(:)-kk%XFFV(:)) * dek%XLEV(:)   
      dek%XLETR(:) = (1.-pek%XPSNV(:)-kk%XFFV(:)) * dek%XLETR(:)  
      dek%XLER (:) = (1.-pek%XPSNV(:)-kk%XFFV(:)) * dek%XLER(:)  
!
!     Total evapotranspiration flux (kg/m2/s):
!
      dk%XEVAP(:) = (dek%XLEV(:) + dek%XLEG(:))/pk%XLVTT(:) + dek%XLEGI(:)/pk%XLSTT(:) + &
                                   dek%XLE_FLOOD(:)/pk%XLVTT(:) + dek%XLEI_FLOOD(:)/pk%XLSTT(:) + &
                                   pek%XPSN(:) * pevap3l(:)
!
!     ISBA-ES/SNOW3L fluxes:
!
      dek%XLES     (:) = pek%XPSN(:) * ples3l(:)
      dek%XLESL    (:) = pek%XPSN(:) * plel3l(:)
      dmk%XRNSNOW   (:) = pek%XPSN(:) * dmk%XRNSNOW   (:)
      dmk%XHSNOW    (:) = pek%XPSN(:) * dmk%XHSNOW    (:)
      dmk%XGFLUXSNOW(:) = pek%XPSN(:) * dmk%XGFLUXSNOW(:)
      dmk%XSNOWHMASS(:) = pek%XPSN(:) * dmk%XSNOWHMASS(:)  ! (J m-2)
      dmk%XHPSNOW   (:) = pek%XPSN(:) * dmk%XHPSNOW   (:)
      pgsfcsnow(:) = pek%XPSN(:) * pgsfcsnow(:)
      pevap3l(:) = pek%XPSN(:) * pevap3l(:)
!
!     Total heat flux between snow and soil
!
      dmk%XGRNDFLUX(:) = pek%XPSN(:) * (pzgrndflux(:)+pflsn_cor(:))
      dek%XMELTADV(:) = pek%XPSN(:) * dek%XMELTADV(:)
!
!     Total evaporative flux (W/m2) :
!
      pek%XLE(:) = dek%XLEG(:) + dek%XLEV(:) + dek%XLES(:) + dek%XLESL(:) + &
                   dek%XLEGI(:) + dek%XLE_FLOOD(:) + dek%XLEI_FLOOD(:)
!
!     Total sublimation flux (W/m2) :
!
      dk%XLEI  (:) = dek%XLES(:) + dek%XLEGI(:) + dek%XLEI_FLOOD(:)
!
!     Total sublimation flux (kg/m2/s):
!
      dk%XSUBL (:) = dk%XLEI(:)/pk%XLSTT(:)
!
!     Total FLUX into snow/soil/vegetation surface:
!
      dk%XGFLUX(:) = dk%XRN(:) - dk%XH(:) - pek%XLE(:) + dmk%XHPSNOW(:)  
!
!     surface humidity:
!
      dk%XQS   (:) = (1.-pek%XPSN(:)) * dk%XQS(:) + pek%XPSN(:) * pqs3l(:)
!
!     near-surface humidity :
!  
      dk%XHU   (:) = (1.-pek%XPSN(:)) * dk%XHU(:) + pek%XPSN(:)
!
!     Momentum fluxes:
!
      pustar(:) = sqrt( (1.-pek%XPSN(:))  * pustar(:)**2  + &
                                 pek%XPSN(:) * dmk%XUSTARSNOW(:)**2 )
!
!     Richardson number and Drag coeff:
!
      CALL comput_ri_drag
!
   ELSE
!
      dk%XTS    (:)  = pek%XTG(:,1)
      dk%XTSRAD (:)  = pek%XTG(:,1)
      dk%XALBT  (:)  = palb(:)
      pemist(:)  = pemis(:)
!  
!     Total sublimation flux (W/m2) :
      dk%XLEI   (:)  = dek%XLES(:) + dek%XLEGI(:) + dek%XLEI_FLOOD(:)
!
!     Total sublimation flux (kg/m2/s):
      dk%XSUBL  (:)  = dk%XLEI(:)/pk%XLSTT(:)
!
   ENDIF
!
ENDIF
!
IF (lhook) CALL dr_hook('ISBA_SNOW_AGR',1,zhook_handle)
!
!-------------------------------------------------------------------------------
CONTAINS
!-------------------------------------------------------------------------------
!
SUBROUTINE comput_ri_drag
!
USE modd_surf_atm, ONLY : ldrag_coef_arp, lrrgust_arp,   &
                          xrrscale, xrrgamma, xutilgust  
!
USE modi_surface_ri
USE modi_surface_aero_cond
USE modi_surface_cd
USE modi_surface_cdch_1darp
USE modi_wind_threshold
!
!*      0.2    declarations of local variables
!
REAL, DIMENSION(SIZE(PTA)) :: ZFP, ZRRCOR, ZVMOD, ZAC, ZRA
!
REAL(KIND=JPRB) :: ZHOOK_HANDLE
!
!-------------------------------------------------------------------------------
!
IF (lhook) CALL dr_hook('ISBA_SNOW_AGR:COMPUT_RI_DRAG',0,zhook_handle)
!
! * Richardson number
!
 CALL surface_ri(dk%XTS, dk%XQS, pexns, pexna, pta, pqa,  &
                pzref, puref, pdircoszw, pvmod, dk%XRI)  
!
! * Wind check
!
zvmod = wind_threshold(pvmod,puref)
!
! * Drag coefficient for heat and momentum
!
IF (ldrag_coef_arp) THEN
   CALL surface_cdch_1darp(pzref, dk%XZ0EFF, dk%XZ0H, zvmod, pta, pek%XTG(:,1), &
                             pqa, dk%XQS, dk%XCD, dk%XCDN, dk%XCH              )
ELSE
   CALL surface_aero_cond(dk%XRI, pzref, puref, zvmod, dk%XZ0, dk%XZ0H, zac, zra, dk%XCH)
   CALL surface_cd(dk%XRI, pzref, puref, dk%XZ0EFF, dk%XZ0H, dk%XCD, dk%XCDN)
ENDIF
!
IF (lrrgust_arp) THEN
   zfp(:)=max(0.0,prr(:)+psr(:))
   zrrcor(:)=sqrt(1.0+((((zfp(:)/(zfp(:)+xrrscale))**xrrgamma)*xutilgust)**2) &
       /(dk%XCD(:)*zvmod(:)**2))  
   dk%XCD  = dk%XCD  * zrrcor
   dk%XCH  = dk%XCH  * zrrcor
   dk%XCDN = dk%XCDN * zrrcor
ENDIF
!
IF (lhook) CALL dr_hook('ISBA_SNOW_AGR:COMPUT_RI_DRAG',1,zhook_handle)
!
END SUBROUTINE comput_ri_drag
!
!-------------------------------------------------------------------------------
!
END SUBROUTINE isba_snow_agr