hydro_glacier.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 hydro_glacier (I, &
                                 ptstep,psr,psnowrho,psnowswe,pglasto,piceflux)
!     ########################################################################
!
!!****  *HYDRO_GLACIER*  
!!
!!    PURPOSE
!!    -------
!
!     Calculate the ice runoff fluxes over permanent snow area with LGLACIER
!     option
!     
!!**  METHOD
!!    ------
!
!     Direct calculation
!
!!    EXTERNAL
!!    --------
!
!     None
!!
!!    IMPLICIT ARGUMENTS
!!    ------------------
!!
!!      
!!    REFERENCE
!!    ---------
!!      
!!    AUTHOR
!!    ------
!!      B. Decharme     
!!
!!    MODIFICATIONS
!!    -------------
!!      Original    09/09 
!-------------------------------------------------------------------------------
!
!*       0.     DECLARATIONS
!               ------------
!
!
USE modd_isba_n, ONLY : isba_t
!
USE modd_csts,     ONLY : xday
USE modd_snow_par, ONLY : xrhosmax, xhgla, xsnowdmin, xrhosmax_es
!
!
USE yomhook   ,ONLY : lhook,   dr_hook
USE parkind1  ,ONLY : jprb
!
IMPLICIT NONE
!
!*      0.1    declarations of arguments
!
!
TYPE(isba_t), INTENT(INOUT) :: i
!
REAL, INTENT(IN)                     :: ptstep
!                                       KTSTEP = timestep [s]
!
REAL, DIMENSION(:), INTENT(IN)       :: psr
!                                       PSR      = Snowfall    [kg/m²s]
!
REAL, DIMENSION(:,:), INTENT(INOUT)  :: psnowrho
!                                       PSNOWRHO = Snow density [kg/m3]
!
REAL, DIMENSION(:,:), INTENT(INOUT)  :: psnowswe
!                                       PSNOWSWE = Snow water equivalent [kg/m²]
!
REAL, DIMENSION(:), INTENT(INOUT)    :: pglasto
!                                       PGLASTO  = Glacier storage      [kg/m²]
!
REAL, DIMENSION(:), INTENT(OUT)      :: piceflux
!                                       PICEFLUX = Ice flux from the Snowfall reservoir [kg/m²s]
!
!
!*      0.2    declarations of local variables
!
REAL, PARAMETER :: ztau=365.25 !days
!
REAL, DIMENSION(SIZE(PSR)) :: zglasto,zstomax,zflux,zsr,zswe
REAL, DIMENSION(SIZE(PSR)) :: zsnowd
!
REAL              ::zrhosmax
!
INTEGER :: jwrk
REAL(KIND=JPRB) :: zhook_handle
!
!-------------------------------------------------------------------------------
!
!-------------------------------------------------------------------------------
!
IF (lhook) CALL dr_hook('HYDRO_GLACIER',0,zhook_handle)
!
zglasto(:) = pglasto(:)
zstomax(:) = 0.0
zflux(:) = 0.0
zsr(:) = 0.0
zswe(:) = 0.0
!
piceflux(:) = 0.0
!
!-------------------------------------------------------------------------------
!Ice accumulation only if snow amount is > to 33.3 meters of aged snow
!
IF(i%TSNOW%SCHEME/='3-L' .AND. i%TSNOW%SCHEME/='CRO')THEN
  zrhosmax=xrhosmax
  zswe(:)=psnowswe(:,1)
ELSE
  zrhosmax=xrhosmax_es
  DO jwrk=1,SIZE(psnowswe,2)
     zswe(:) = zswe(:) + psnowswe(:,jwrk)
  END DO
ENDIF
!
WHERE(zswe(:)>xhgla*zrhosmax)
     zsr(:) = psr(:)
ELSEWHERE
     zsr(:) = 0.0
ENDWHERE
!
!Snow storage calculation
!
pglasto(:)=(zglasto(:)+ptstep*zsr(:))/(1.0+ptstep/(ztau*xday))
!
!supress numerical artifacs
!
zstomax(:)=zsr(:)*ptstep+zglasto(:)
!
pglasto(:)=min(zstomax(:),pglasto(:))
!
!Ice flux calculation                
!
zflux(:)=(zglasto(:)-pglasto(:))/ptstep+zsr(:)
!      
!supress numerical artifacs
!
piceflux(:) = max(0.0,zflux(:))
pglasto(:) = pglasto(:) + piceflux(:)-zflux(:)             
!
WHERE(pglasto(:)<=1.e-10)pglasto(:)=0.0
!
!-------------------------------------------------------------------------------
!Snow pack update
!
IF(i%TSNOW%SCHEME/='3-L' .AND. i%TSNOW%SCHEME/='CRO')THEN
!
  WHERE(psnowswe(:,1)<=xhgla*zrhosmax)piceflux(:)=0.0
  psnowswe(:,1)=psnowswe(:,1)-piceflux(:)*ptstep
!
ELSE
!
  WHERE(zswe(:)<=xhgla*zrhosmax)piceflux(:)=0.0
!
! Snow total depth
  zsnowd(:) = 0.
  DO jwrk=1,SIZE(psnowswe,2)
     zsnowd(:) = zsnowd(:) + psnowswe(:,jwrk)/psnowrho(:,jwrk)
  END DO
!
! Flux
  DO jwrk=1,SIZE(psnowswe,2)
     zflux(:) = piceflux(:)*(psnowswe(:,jwrk)/psnowrho(:,jwrk)) &
                /max(zsnowd(:),0.0001)
     psnowswe(:,jwrk)=psnowswe(:,jwrk)-zflux(:)*ptstep
  END DO
!
ENDIF
!-------------------------------------------------------------------------------
IF (lhook) CALL dr_hook('HYDRO_GLACIER',1,zhook_handle)
!-------------------------------------------------------------------------------
END SUBROUTINE hydro_glacier