diag_inline_seafluxn.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 diag_inline_seaflux_n (DGS, DGSI, S, &
                                  ptstep, pta, pqa, &
     ppa, pps, prhoa, pzona,                             &
     pmera, pht, phw, pcd, pcdn, pch, pce, pri, phu,     &
     pz0h, pqsat, psfth, psftq, psfzon, psfmer,     &
     pdir_sw, psca_sw, plw, pdir_alb, psca_alb, &
     pemis, ptrad, prain, psnow,                         & 
     pcd_ice, pcdn_ice, pch_ice, pce_ice, pri_ice,       &
     pz0_ice, pz0h_ice, pqsat_ice, psfth_ice, psftq_ice, &
     psfzon_ice, psfmer_ice )
                                          
!     #####################################################################################
!
!!****  *DIAG_INLINE_SEAFLUX_n * - computes diagnostics during SEAFLUX time-step
!!
!!    PURPOSE
!!    -------
!
!!**  METHOD
!!    ------
!!
!!    REFERENCE
!!    ---------
!!      
!!
!!    AUTHOR
!!    ------
!!     V. Masson 
!!
!!    MODIFICATIONS
!!    -------------
!!      Original    01/2004
!!      Modified    01/2006 : sea flux parameterization.
!!      B. Decharme 08/2009 : Diag for Earth System Model Coupling
!!      S. Riette   06/2009 CLS_2M becomes CLS_TQ, CLS_TQ and CLS_WIND have one
!!                          more argument (height of diagnostic)
!!      B. Decharme 04/2013 : Add EVAP and SUBL diag
!!      S. Senesi   01/2014 ! introduce fractional seaice and sea-ice model 
!!------------------------------------------------------------------
!

!
!
!
USE modd_diag_seaflux_n, ONLY : diag_seaflux_t
USE modd_diag_seaice_n, ONLY : diag_seaice_t
USE modd_seaflux_n, ONLY : seaflux_t
!
USE modd_csts,           ONLY : xtts
USE modd_surf_par,       ONLY : xundef
USE modd_sfx_oasis,      ONLY : lcpl_sea
!
USE modd_types_glt,     ONLY : t_glt
USE modd_glt_param ,    ONLY : gelato_dim=>nx
USE mode_glt_stats ,    ONLY : glt_avhicem, glt_avhsnwm
USE modi_param_cls
USE modi_cls_tq
USE modi_cls_wind
USE modi_diag_surf_budget_sea
USE modi_diag_surf_budgetc_sea
USE modi_diag_cpl_esm_sea
!
USE modi_seaflux_albedo
!
USE yomhook   ,ONLY : lhook,   dr_hook
USE parkind1  ,ONLY : jprb
!
IMPLICIT NONE
!
!*      0.1    declarations of arguments
!
!
TYPE(diag_seaflux_t), INTENT(INOUT) :: dgs
TYPE(diag_seaice_t), INTENT(INOUT) :: dgsi
TYPE(seaflux_t), INTENT(INOUT) :: s
!
REAL,               INTENT(IN) :: ptstep ! atmospheric time-step                 (s)
REAL, DIMENSION(:), INTENT(IN) :: pta    ! atmospheric temperature
REAL, DIMENSION(:), INTENT(IN) :: pqa    ! atmospheric specific humidity
REAL, DIMENSION(:), INTENT(IN) :: ppa    ! atmospheric level pressure
REAL, DIMENSION(:), INTENT(IN) :: pps    ! surface pressure
REAL, DIMENSION(:), INTENT(IN) :: prhoa  ! air density
REAL, DIMENSION(:), INTENT(IN) :: pzona  ! zonal wind
REAL, DIMENSION(:), INTENT(IN) :: pmera  ! meridian wind
REAL, DIMENSION(:), INTENT(IN) :: pht    ! atmospheric level height
REAL, DIMENSION(:), INTENT(IN) :: phw    ! atmospheric level height for wind
REAL, DIMENSION(:), INTENT(IN) :: pcd    ! drag coefficient for momentum
REAL, DIMENSION(:), INTENT(IN) :: pcdn   ! neutral drag coefficient
REAL, DIMENSION(:), INTENT(IN) :: pch    ! drag coefficient for heat
REAL, DIMENSION(:), INTENT(IN) :: pce    ! drag coefficient for vapor
REAL, DIMENSION(:), INTENT(IN) :: pri    ! Richardson number
REAL, DIMENSION(:), INTENT(IN) :: phu    ! near-surface humidity
REAL, DIMENSION(:), INTENT(IN) :: pz0h   ! roughness length for heat
REAL, DIMENSION(:), INTENT(IN) :: pqsat  ! humidity at saturation
REAL, DIMENSION(:), INTENT(IN) :: psfzon ! zonal friction
REAL, DIMENSION(:), INTENT(IN) :: psfmer ! meridian friction
REAL, DIMENSION(:), INTENT(IN) :: psfth  ! heat flux  (W/m2)
REAL, DIMENSION(:), INTENT(IN) :: psftq  ! water flux (kg/m2/s)
REAL, DIMENSION(:,:),INTENT(IN):: pdir_sw   ! direct  solar radiation (on horizontal surf.)
!                                           !                                       (W/m2)
REAL, DIMENSION(:,:),INTENT(IN):: psca_sw   ! diffuse solar radiation (on horizontal surf.)
!                                           !                                       (W/m2)
REAL, DIMENSION(:), INTENT(IN) :: plw       ! longwave radiation (on horizontal surf.)
REAL, DIMENSION(:), INTENT(IN) :: ptrad     ! radiative temperature                 (K)
REAL, DIMENSION(:,:),INTENT(IN):: pdir_alb  ! direct albedo for each spectral band  (-)
REAL, DIMENSION(:,:),INTENT(IN):: psca_alb  ! diffuse albedo for each spectral band (-)
REAL, DIMENSION(:), INTENT(IN) :: pemis     ! emissivity                            (-)
!
REAL, DIMENSION(:), INTENT(IN) :: prain     ! Rainfall (kg/m2/s)
REAL, DIMENSION(:), INTENT(IN) :: psnow     ! Snowfall (kg/m2/s)
!
REAL, DIMENSION(:), INTENT(IN)    :: pcd_ice    ! drag coefficient for momentum
REAL, DIMENSION(:), INTENT(IN)    :: pcdn_ice   ! neutral drag coefficient
REAL, DIMENSION(:), INTENT(IN)    :: pch_ice    ! drag coefficient for heat
REAL, DIMENSION(:), INTENT(IN)    :: pce_ice    ! drag coefficient for vapor
REAL, DIMENSION(:), INTENT(IN)    :: pri_ice    ! Richardson number
REAL, DIMENSION(:), INTENT(IN)    :: pz0_ice    ! roughness length for momentum
REAL, DIMENSION(:), INTENT(IN)    :: pz0h_ice   ! roughness length for heat
REAL, DIMENSION(:), INTENT(IN)    :: pqsat_ice  ! humidity at saturation
REAL, DIMENSION(:), INTENT(IN)    :: psfth_ice  ! heat flux  (W/m2)
REAL, DIMENSION(:), INTENT(IN)    :: psftq_ice  ! water flux (kg/m2/s)
REAL, DIMENSION(:), INTENT(IN)    :: psfzon_ice ! zonal friction
REAL, DIMENSION(:), INTENT(IN)    :: psfmer_ice ! meridian friction
!
!*      0.2    declarations of local variables
!
LOGICAL                         :: gsic
REAL, DIMENSION(SIZE(PTA))      :: zz0w
REAL, DIMENSION(SIZE(PTA))      :: zh
!
REAL(KIND=JPRB) :: zhook_handle
!-------------------------------------------------------------------------------------
!
IF (lhook) CALL dr_hook('DIAG_INLINE_SEAFLUX_N',0,zhook_handle)
!
! * Mean surface temperature need to couple with AGCM
!
IF (s%LHANDLE_SIC) THEN
   dgs%XTS   (:) = (1 - s%XSIC(:)) * s%XSST(:) + s%XSIC(:) * s%XTICE(:)
   dgs%XTSRAD(:) = ptrad(:)
ELSE
   dgs%XTS   (:) = s%XSST (:)
   dgs%XTSRAD(:) = ptrad(:)
ENDIF
!
IF (.NOT. s%LSBL) THEN
!
  IF (dgs%N2M==1) THEN        
    CALL param_cls(pta, s%XSST, pqa, ppa, prhoa, pzona, pmera, pht, phw, &
         psfth, psftq, psfzon, psfmer,                                 &
         dgs%XT2M, dgs%XQ2M, dgs%XHU2M, dgs%XZON10M, dgs%XMER10M )  
    IF (s%LHANDLE_SIC) THEN
       CALL param_cls(pta, s%XTICE, pqa, ppa, prhoa, pzona, pmera, pht, phw, &
            psfth_ice, psftq_ice, psfzon_ice, psfmer_ice,                  &
            dgs%XT2M_ICE, dgs%XQ2M_ICE, dgs%XHU2M_ICE, dgs%XZON10M_ICE, dgs%XMER10M_ICE  )  
    ENDIF
  ELSE IF (dgs%N2M==2) THEN
    zh(:)=2.          
    CALL cls_tq(pta, pqa, ppa, pps, pht,          &
                  pcd, pch, pri,                  &
                  s%XSST, phu, pz0h, zh,            &
                  dgs%XT2M, dgs%XQ2M, dgs%XHU2M)
    zh(:)=10.                
    CALL cls_wind(pzona, pmera, phw,              &
                    pcd, pcdn, pri, zh,           &
                    dgs%XZON10M, dgs%XMER10M)  
    IF (s%LHANDLE_SIC) THEN
       zh(:)=2.          
       CALL cls_tq(pta, pqa, ppa, pps, pht,  &
            pcd_ice, pch_ice, pri_ice,       &
            s%XTICE, phu, pz0h_ice, zh,        &
            dgs%XT2M_ICE, dgs%XQ2M_ICE, dgs%XHU2M_ICE)  
       zh(:)=10.                
       CALL cls_wind(pzona, pmera, phw,      &
            pcd_ice, pcdn_ice, pri_ice, zh,  &
            dgs%XZON10M_ICE, dgs%XMER10M_ICE  )  
    ENDIF 
  END IF
!
  IF (dgs%N2M>=1) THEN
     IF (s%LHANDLE_SIC) THEN
        !
        dgs%XT2M    = dgs%XT2M    * (1 - s%XSIC) + dgs%XT2M_ICE    * s%XSIC
        dgs%XQ2M    = dgs%XQ2M    * (1 - s%XSIC) + dgs%XQ2M_ICE    * s%XSIC
        dgs%XHU2M   = dgs%XHU2M   * (1 - s%XSIC) + dgs%XHU2M_ICE   * s%XSIC
        !
        dgs%XZON10M(:) = dgs%XZON10M(:) * (1 - s%XSIC(:)) + dgs%XZON10M_ICE(:) * s%XSIC(:)
        dgs%XMER10M(:) = dgs%XMER10M(:) * (1 - s%XSIC(:)) + dgs%XMER10M_ICE(:) * s%XSIC(:)
        dgs%XWIND10M_ICE(:) = sqrt(dgs%XZON10M_ICE(:)**2+dgs%XMER10M_ICE(:)**2)
        !
        dgs%XRI    = pri     * (1 - s%XSIC) + pri_ice     * s%XSIC
        dgs%XRI_ICE=pri_ice
     ELSE
        dgs%XRI    =pri
     ENDIF
    !
    dgs%XT2M_MIN(:) = min(dgs%XT2M_MIN(:),dgs%XT2M(:))
    dgs%XT2M_MAX(:) = max(dgs%XT2M_MAX(:),dgs%XT2M(:))
    !
    dgs%XHU2M_MIN(:) = min(dgs%XHU2M_MIN(:),dgs%XHU2M(:))
    dgs%XHU2M_MAX(:) = max(dgs%XHU2M_MAX(:),dgs%XHU2M(:))
    !
    dgs%XWIND10M(:) = sqrt(dgs%XZON10M(:)**2+dgs%XMER10M(:)**2)
    dgs%XWIND10M_MAX(:) = max(dgs%XWIND10M_MAX(:),dgs%XWIND10M(:))
    !
  ENDIF
!
ELSE
  IF (dgs%N2M>=1) THEN
    dgs%XT2M    = xundef
    dgs%XQ2M    = xundef
    dgs%XHU2M   = xundef
    dgs%XZON10M = xundef
    dgs%XMER10M = xundef
    dgs%XRI     = pri
  ENDIF
ENDIF
!
IF (dgs%LSURF_BUDGET.OR.dgs%LSURF_BUDGETC) THEN
!
  CALL seaflux_albedo(pdir_sw,psca_sw,pdir_alb,psca_alb,dgs%XALBT)
!
  CALL diag_surf_budget_sea(xtts, s%XSST, prhoa, psfth, psfth_ice,    &
                                 psftq, psftq_ice,                     &
                                 pdir_sw, psca_sw, plw, pdir_alb,      &
                                 psca_alb,s%XICE_ALB, pemis, ptrad,      &
                                 psfzon, psfzon_ice, psfmer,           &
                                 psfmer_ice, s%LHANDLE_SIC, s%XSIC, s%XTICE, &
                                 dgs%XRN, dgs%XH, dgs%XLE, dgs%XLE_ICE, dgs%XGFLUX,        &
                                 dgs%XSWD, dgs%XSWU, dgs%XSWBD, dgs%XSWBU, dgs%XLWD, dgs%XLWU, &
                                 dgs%XFMU, dgs%XFMV, dgs%XEVAP, dgs%XSUBL,             &
                                 dgs%XRN_ICE, dgs%XH_ICE, dgs%XGFLUX_ICE,          &
                                 dgs%XSWU_ICE, dgs%XSWBU_ICE, dgs%XLWU_ICE,        &
                                 dgs%XFMU_ICE, dgs%XFMV_ICE                    ) 
!
END IF
!
IF(dgs%LSURF_BUDGETC)THEN
  CALL diag_surf_budgetc_sea(dgs, &
                             ptstep, dgs%XRN, dgs%XH, dgs%XLE, dgs%XLE_ICE, dgs%XGFLUX,  &
                               dgs%XSWD, dgs%XSWU, dgs%XLWD, dgs%XLWU, dgs%XFMU, dgs%XFMV,   &
                               dgs%XEVAP, dgs%XSUBL, s%LHANDLE_SIC,            &
                               dgs%XRN_ICE, dgs%XH_ICE, dgs%XGFLUX_ICE,          &
                               dgs%XSWU_ICE, dgs%XLWU_ICE, dgs%XFMU_ICE, dgs%XFMV_ICE)
ENDIF
!
IF (dgs%LCOEF) THEN
   IF (s%LHANDLE_SIC) THEN 
      !
      !* Transfer coefficients
      !
      dgs%XCD = (1 - s%XSIC) * pcd + s%XSIC * pcd_ice
      dgs%XCH = (1 - s%XSIC) * pch + s%XSIC * pch_ice
      dgs%XCE = (1 - s%XSIC) * pce + s%XSIC * pce_ice
      !
      !* Roughness lengths
      !
      zz0w = ( 1 - s%XSIC ) * 1.0/(log(phw/s%XZ0)    **2)  +  &
                   s%XSIC   * 1.0/(log(phw/pz0_ice)**2)  
      dgs%XZ0  = phw  * exp( - sqrt( 1./  zz0w ))
      zz0w = ( 1 - s%XSIC ) * 1.0/(log(phw/pz0h)    **2)  +  &
                   s%XSIC   * 1.0/(log(phw/pz0h_ice)**2)  
      dgs%XZ0H = phw  * exp( - sqrt( 1./  zz0w ))

      dgs%XCD_ICE  = pcd_ice
      dgs%XCH_ICE  = pch_ice
      dgs%XZ0_ICE  = pz0_ice
      dgs%XZ0H_ICE = pz0h_ice
      !
   ELSE
      !
      !* Transfer coefficients
      !
      dgs%XCD = pcd
      dgs%XCH = pch
      dgs%XCE = pce
      !
      !* Roughness lengths
      !
      dgs%XZ0  = s%XZ0
      dgs%XZ0H = pz0h
   ENDIF
   !
ENDIF
!
IF (dgs%LSURF_VARS) THEN
  !
  !* Humidity at saturation
  !
   IF (s%LHANDLE_SIC) THEN 
      dgs%XQS     = (1 - s%XSIC) * pqsat + s%XSIC * pqsat_ice
      dgs%XQS_ICE = pqsat_ice
   ELSE 
      dgs%XQS = pqsat
   ENDIF
ENDIF
!
! Diags from embedded Seaice model
! CALL DIAG_INLINE_SEAICE() : simply  : 
!
IF (dgsi%LDIAG_SEAICE) THEN
   IF (trim(s%CSEAICE_SCHEME) == 'GELATO') THEN 
      gelato_dim=SIZE(pta)
      dgsi%XSIT  = reshape(glt_avhicem(s%TGLT%dom,s%TGLT%sit),(/gelato_dim/))
      dgsi%XSND  = reshape(glt_avhsnwm(s%TGLT%dom,s%TGLT%sit),(/gelato_dim/))
      dgsi%XMLT  = s%TGLT%oce_all(:,1)%tml
   ELSE
      ! Placeholder for an alternate seaice scheme
   ENDIF
ENDIF
!
! Diags for Earth System Model coupling or for embedded Seaice model
! (we are actually using XCPL_.. variables for feeding the seaice model)
!
gsic=(s%LHANDLE_SIC.AND.(s%CSEAICE_SCHEME /= 'NONE  '))
!
IF (lcpl_sea.OR.gsic) THEN
!
  CALL diag_cpl_esm_sea(s, &
                        ptstep,dgs%XZON10M,dgs%XMER10M,dgs%XFMU,dgs%XFMV,  &
                          dgs%XSWD,dgs%XSWU,dgs%XGFLUX,psftq,prain,    &
                          psnow,plw,s%XTICE,psfth_ice,       &
                          psftq_ice,pdir_sw,psca_sw,       &
                          dgs%XSWU_ICE,dgs%XLWU_ICE,gsic           )
! 
ENDIF
IF (lhook) CALL dr_hook('DIAG_INLINE_SEAFLUX_N',1,zhook_handle)
!
!-------------------------------------------------------------------------------------
!
END SUBROUTINE diag_inline_seaflux_n