init_teb_greenroof_pgdn.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 init_teb_greenroof_pgd_n (DTCO, U, OCH_BIO_FLUX, G, PGREENROOF, TOP, IO, S, K, P, PEK, DTV, GB, &
                                     HPROGRAM, HINIT, OPATCH1, KI, KVERSION, PCO2, PRHOA)
!#############################################################
!
!!****  *INIT_TEB_GREENROOF_PGD_n* - routine to initialize ISBA
!!
!!    PURPOSE
!!    -------
!!
!!**  METHOD
!!    ------
!!
!!    EXTERNAL
!!    --------
!!
!!
!!    IMPLICIT ARGUMENTS
!!    ------------------
!!
!!    REFERENCE
!!    ---------
!!
!!
!!    AUTHOR
!!    ------
!!      A. Lemonsu  *Meteo France*
!!
!!    MODIFICATIONS
!!    -------------
!!      Original    09/2009
!!                  11/2013 (B. Decharme) No exp profile with DIF
!-------------------------------------------------------------------------------
!
!*       0.    DECLARATIONS
!              ------------
!
USE modd_data_cover_n, ONLY : data_cover_t
USE modd_surf_atm_n, ONLY : surf_atm_t
USE modd_sso_n, ONLY : sso_t, sso_init
USE modd_sfx_grid_n, ONLY : grid_t
USE modd_teb_option_n, ONLY : teb_options_t
!
USE modd_isba_options_n, ONLY : isba_options_t
USE modd_isba_n, ONLY : isba_s_t, isba_k_t, isba_p_t, isba_pe_t
!
USE modd_data_isba_n, ONLY : data_isba_t
USE modd_gr_biog_n, ONLY : gr_biog_t
!
USE modd_agri_n, ONLY : agri_t
!
USE modd_type_date_surf
USE modd_type_snow
!
USE modd_data_cover_par,       ONLY: nvegtype
USE modd_surf_par,             ONLY: xundef, nundef

USE modd_sgh_par,              ONLY: xf_decay
!
USE modi_read_prep_greenroof_snow
USE modi_get_luout
USE modi_allocate_teb_veg_pgd
USE modi_read_pgd_teb_greenroof_n
USE modi_convert_patch_isba
USE modi_init_from_data_teb_veg_n
USE modi_init_veg_pgd_n
USE modi_exp_decay_soil_fr
USE modi_abor1_sfx
USE modi_av_pgd
!
USE mode_teb_veg
!
USE yomhook   ,ONLY : lhook,   dr_hook
USE parkind1  ,ONLY : jprb
!
IMPLICIT NONE
!
!*       0.1   Declarations of arguments
!              -------------------------
!
!
TYPE(data_cover_t), INTENT(INOUT) :: DTCO
TYPE(surf_atm_t), INTENT(INOUT) :: U
LOGICAL, INTENT(IN) :: OCH_BIO_FLUX
TYPE(grid_t), INTENT(INOUT) :: G
REAL, DIMENSION(:), INTENT(IN) :: PGREENROOF
TYPE(teb_options_t), INTENT(INOUT) :: TOP
!
TYPE(isba_options_t), INTENT(INOUT) :: IO
TYPE(isba_s_t), INTENT(INOUT) :: S
TYPE(isba_k_t), INTENT(INOUT) :: K
TYPE(isba_p_t), INTENT(INOUT) :: P
TYPE(isba_pe_t), INTENT(INOUT) :: PEK
!
TYPE(data_isba_t), INTENT(INOUT) :: DTV
TYPE(gr_biog_t), INTENT(INOUT) :: GB
!
 CHARACTER(LEN=6),                   INTENT(IN)  :: HPROGRAM  ! program calling surf. schemes
 CHARACTER(LEN=3),                   INTENT(IN)  :: HINIT     ! choice of fields to initialize
LOGICAL,                            INTENT(IN)  :: OPATCH1 ! flag to read PGD fields in the file
INTEGER,                            INTENT(IN)  :: KI        ! number of points
INTEGER,                            INTENT(IN)  :: KVERSION  ! version number of the file being read
REAL,             DIMENSION(KI),    INTENT(IN)  :: PCO2        ! CO2 concentration (kg/m3)
REAL,             DIMENSION(KI),    INTENT(IN)  :: PRHOA       ! air density
!
!
!
!*       0.2   Declarations of local variables
!              -------------------------------
!
TYPE(sso_t) :: YSS
TYPE(agri_t) :: YAG
!
INTEGER           :: JILU     ! loop increment
INTEGER           :: ILUOUT   ! unit of output listing file
!
INTEGER           :: IDECADE  ! decade of simulation
!
INTEGER :: JVEG, JL, JI  ! loop counter on layers
!
REAL, DIMENSION(KI)               :: ZF
REAL, DIMENSION(KI)               :: ZWORK
!
!*       0.3   Soil parameter values for organic matter - from Lawrence and Slater (2008):
!              ----------------------------------------------------------------------------------
!
REAL, PARAMETER   :: ZWSAT_OM      = 0.9       ! Porosity of OM (m3/m3)
REAL, PARAMETER   :: ZCONDSAT_OM   = 2.8e-4    ! Saturated hydraulic conductivity for OM (m/s)
REAL, PARAMETER   :: ZMPOTSAT_OM   = -10.3e-3  ! Saturated matric potential for OM (m)
REAL, PARAMETER   :: ZBCOEF_OM     = 2.7       ! CH78 b-parameter for OM (-)
!
REAL, PARAMETER   :: ZCONDDRY_OM   = 0.05      ! Dry thermal conductivity for OM (W/m/K)
REAL, PARAMETER   :: ZCONDSLD_OM   = 0.25      ! Soil solids thermal conductivity for OM (W/m/K)
REAL, PARAMETER   :: ZHCAPSOIL_OM  = 2.5e+6    ! Soil heat capacity for OM
!
REAL, PARAMETER   :: ZMPOT_WWILT   = -150.     ! Matric potential at wilting point (m)
REAL, PARAMETER   :: ZHYDCOND_WFC  = 1.157e-9  ! Hydraulic conductivity at field capacity (m/s)
!
REAL, DIMENSION(0) :: ZTDEEP_CLI, ZGAMMAT_CLI, ZTHRESHOLD
!
REAL(KIND=JPRB) :: ZHOOK_HANDLE
!
!-------------------------------------------------------------------------------
!
!               Initialisation for IO
!
IF (lhook) CALL dr_hook('INIT_TEB_GREENROOF_PGD_n',0,zhook_handle)
!
 CALL get_luout(hprogram,iluout)
!
 CALL sso_init(yss)
!-------------------------------------------------------------------------------
!
!*       1.     Reading of snow configuration:
!               ------------------------------
!
!* initialization of snow scheme (TSNOW defined in MODD_TEB_GREENROOF_n)
!
IF (hinit=='PRE') THEN
   CALL read_prep_greenroof_snow(hprogram,pek%TSNOW%SCHEME,pek%TSNOW%NLAYER)
!
   IF (pek%TSNOW%SCHEME.NE.'3-L' .AND. pek%TSNOW%SCHEME.NE.'CRO' .AND. io%CISBA=='DIF') THEN
    CALL abor1_sfx("INIT_TEB_GREENROOF_n: WITH CISBA_GR = DIF, CSNOW MUST BE 3-L OR CRO")
  ENDIF
ENDIF
!
!-------------------------------------------------------------------------------
!
!*       2.     Physiographic fields
!               --------------------
!
!
!*       2.1    Cover, soil and orographic fields:
!               ---------------------------------
!
 CALL allocate_teb_veg_pgd(pek, s, k, p, opatch1, ki, nvegtype, io%NGROUND_LAYER )
!
IF (top%TTIME%TDATE%MONTH /= nundef) THEN
  idecade = 3 * ( top%TTIME%TDATE%MONTH - 1 ) + min(top%TTIME%TDATE%DAY-1,29) / 10 + 1
ELSE
  idecade = 1
END IF
!
IF (opatch1) THEN

  CALL read_pgd_teb_greenroof_n(och_bio_flux, dtco, dtv, gb, u, &
                                io, s, k, g%NDIM, hprogram,kversion)
  !
  ALLOCATE(s%XVEGTYPE(ki,nvegtype))
  IF (io%LPAR) THEN
    s%XVEGTYPE = dtv%XPAR_VEGTYPE
  ELSE
    !classical ecoclimap case
    DO jveg=1,nvegtype
      CALL av_pgd(dtco, s%XVEGTYPE(:,jveg),top%XCOVER ,dtco%XDATA_VEGTYPE(:,jveg),'GRD','ARI',top%LCOVER)
    END DO
  ENDIF
  DO jveg=1,nvegtype
    WHERE (pgreenroof==0)
      s%XVEGTYPE(:,jveg) = 0.
      s%XVEGTYPE(:,1) = 1.
    END WHERE
  ENDDO  
  !
  ALLOCATE(s%XPATCH(ki,1),p%XPATCH(ki))
  ALLOCATE(s%XVEGTYPE_PATCH(ki,nvegtype,1),p%XVEGTYPE_PATCH(ki,nvegtype))
  s%XPATCH(:,1) = 1.
  p%XPATCH(:) = s%XPATCH(:,1)
  s%XVEGTYPE_PATCH(:,:,1) = s%XVEGTYPE
  p%XVEGTYPE_PATCH(:,:) = s%XVEGTYPE_PATCH(:,:,1)
  p%NSIZE_P = ki
  ALLOCATE(p%NR_P(ki))
  DO ji = 1,SIZE(p%NR_P)
    p%NR_P(ji) = ji
  ENDDO  
  !
  IF (.NOT. io%LPAR) THEN
    CALL convert_patch_isba(dtco, dtv, io, idecade, idecade, top%XCOVER, top%LCOVER,&
                        .false.,'GRD', 1, k, p, pek, &
                        .true., .false., .false., .false., .false., .false., &
                        psoilgrid=io%XSOILGRID  )   
  ELSE
    CALL init_from_data_teb_veg_n(dtv, k, p, pek, idecade, .false., .true., .false.,.false.)
  ENDIF
  ! 
  ALLOCATE(s%XWSN_WR(0,0,1))
  ALLOCATE(s%XRHO_WR(0,0,1))
  ALLOCATE(s%XALB_WR(0,1))
  ALLOCATE(s%XHEA_WR(0,0,1))
  ALLOCATE(s%XAGE_WR(0,0,1))
  ALLOCATE(s%XSG1_WR(0,0,1))
  ALLOCATE(s%XSG2_WR(0,0,1)) 
  ALLOCATE(s%XHIS_WR(0,0,1))
  !
END IF
!
!*       2.2    Physiographic data fields from land cover:
!               -----------------------------------------
!
IF (.NOT. io%LPAR) THEN
  CALL convert_patch_isba(dtco, dtv, io, idecade, idecade, top%XCOVER, top%LCOVER,&
                        .false.,'GRD', 1, k, p, pek, &
                        .false., .true., .false., .false., .false., .false.  )   
ELSE

  CALL init_from_data_teb_veg_n(dtv, k, p, pek, idecade, .false., .false., .true.,.false.)

  IF (io%CISBA=='DIF') CALL init_if_dif(io%NGROUND_LAYER, pgreenroof, p)

END IF
!
 CALL init_if_noveg(pgreenroof, io, s, p, pek)
!
ALLOCATE(k%XVEGTYPE(ki,nvegtype))
k%XVEGTYPE = s%XVEGTYPE
!
ALLOCATE(yss%XAOSIP(0))
!
 CALL init_veg_pgd_n(yss, dtv, io, s, k, k, p, pek, yag, ki,                     &
                      hprogram, 'TOWN  ',iluout, ki, top%TTIME%TDATE%MONTH, &
                      .false., .false., ztdeep_cli, zgammat_cli,            &
                      .false., zthreshold, hinit, pco2, prhoa  )
!
!-------------------------------------------------------------------------------
!
IF (opatch1) THEN
  !
  !*       5.1     Soil thermal characteristics for greenroofs:
  !               ----------------------------------------------
  !
  ! WARNING: must be done before soil hydraulic characteristics (because of WSAT)
  ! Estimation of WSAT_MI for use in HEATCAPZ and THRMCONDZ for mineral fraction
  ! and allow weighted combination with regard to OM & no-OM fractions:
  !
  IF (io%CSCOND=='PL98' .OR. io%CISBA=='DIF') THEN
    DO jl=1,io%NGROUND_LAYER
      k%XHCAPSOIL(:,jl) = s%XSOC(:,jl) * zhcapsoil_om + (1-s%XSOC(:,jl)) * k%XHCAPSOIL(:,jl)  
      k%XCONDDRY (:,jl) = (zconddry_om * k%XCONDDRY(:,jl)) / &
                        ( s%XSOC(:,jl) * k%XCONDDRY(:,jl) + (1-s%XSOC(:,jl)) * zconddry_om )
      k%XCONDSLD (:,jl) = (zcondsld_om * k%XCONDSLD(:,jl)) / &
                        ( s%XSOC(:,jl) * k%XCONDSLD(:,jl) + (1-s%XSOC(:,jl)) * zcondsld_om )
    ENDDO
  END IF
  !
  ! +++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
  ! Validation case : experimental values for Nancy 2011 case
  ! +++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
  ! Substrate layer
  DO jl=1,4
    k%XCONDDRY (:,jl) = 0.15
    k%XHCAPSOIL(:,jl) = 1342000.
  ENDDO
  ! Drainage layer
  DO jl=5,6
    k%XCONDDRY (:,jl) = 0.09
    k%XHCAPSOIL(:,jl) = 331500.
  ENDDO
  ! +++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
  !
ENDIF
!
!*       5.2     Soil thermal characteristics:
!               --------------------------------
!
DO jl=1,io%NGROUND_LAYER
  p%XCONDSAT(:,jl) = s%XSOC(:,jl)* zcondsat_om + (1-s%XSOC(:,jl)) * p%XCONDSAT(:,jl)
END DO
!
!
IF (opatch1) THEN
  !
  ! Note that if ISBA/=DIF, always CDIF = 'BC' and CPEDOTF = 'CH78'
  DO jl=1,io%NGROUND_LAYER
    k%XBCOEF  (:,jl) = s%XSOC(:,jl) * zbcoef_om   + (1-s%XSOC(:,jl)) * k%XBCOEF(:,jl)
    k%XMPOTSAT(:,jl) = s%XSOC(:,jl) * zmpotsat_om + (1-s%XSOC(:,jl)) * k%XMPOTSAT(:,jl)
  END DO
  !        
  DO jl=1,io%NGROUND_LAYER
    k%XWSAT (:,jl) =  s%XSOC(:,jl)* zwsat_om +(1-s%XSOC(:,jl))* k%XWSAT(:,jl)
    k%XWWILT(:,jl) = exp(((log(-1*zmpot_wwilt)-log(-1*k%XMPOTSAT(:,jl)))   &
                    / (-1*k%XBCOEF(:,jl)))+log(k%XWSAT(:,jl)))
    k%XWFC  (:,jl) = exp(((log(zhydcond_wfc)-log(p%XCONDSAT(:,jl)))      &
                    / (2*k%XBCOEF(:,jl)+3))+log(k%XWSAT(:,jl)))
  END DO
  !
  ! +++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
  ! Validation case : experimental values for Nancy 2011 case
  ! +++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
  ! Substrate layer
  DO jl=1,4
    k%XWSAT   (:,jl) = 0.674     ! Value tested
    k%XMPOTSAT(:,jl) = -0.932    ! Value tested
    k%XBCOEF  (:,jl) = 3.9       ! Value tested
    k%XWWILT  (:,jl) = 0.15      ! from OBS-NANCY
    k%XWFC    (:,jl) = 0.37      ! from OBS-NANCY
  ENDDO
  ! Drainage layer
  DO jl=5,6
    k%XWSAT   (:,jl) = 0.9       ! Value tested
    k%XMPOTSAT(:,jl) = -0.121    ! Value tested
    k%XBCOEF  (:,jl) = 2.7       ! Value tested
    k%XWWILT  (:,jl) = 0.15      ! sert à initialiser le WG ds la couche
    k%XWFC    (:,jl) = 0.37      ! sert à initialiser le WG ds la couche
  ENDDO
  !
ENDIF
!
! Substrate layer
DO jl=1,4
  p%XCONDSAT(:,jl) = 2.162e-3  ! Value tested
ENDDO
! Drainage layer
DO jl=5,6
  p%XCONDSAT(:,jl) = 3.32e-3   ! Value tested
ENDDO
!
!-------------------------------------------------------------------------------
!
!*       6.1    Initialize of the SGH scheme:'
!               ------------------------------
!
IF(io%CKSAT=='SGH' .AND. io%CISBA/='DIF' .AND. hinit/='PRE')THEN 
  zf(:)=min(4.0/p%XDG(:,2),xf_decay)
  CALL exp_decay_soil_fr(io%CISBA,  zf, p)
ENDIF
!
!-------------------------------------------------------------------------------
!
IF (lhook) CALL dr_hook('INIT_TEB_GREENROOF_PGD_n',1,zhook_handle)
!
!-------------------------------------------------------------------------------
!
!
END SUBROUTINE init_teb_greenroof_pgd_n