nitro_decline.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 nitro_decline(HPHOTO, HRESPSL, OTR_ML, KSPINW ,            &
                pbslai_nitro, psefold, pgmes, panmax, panday,         &
                plat, plaimin, pvegtype, ptau_wood,                   &
                panfm, plai, pbiomass, presp_biomass, pbiomass_leaf,  &
                pincrease ,pturnover                                  )  
!
!   ###############################################################
!!**  NITRO_DECLINE 
!!
!!    PURPOSE
!!    -------
!!
!!**  METHOD
!!    ------
!!     Calvet and Soussana (2001) and Gibelin et al. (2006) for nitrogen dilution.
!!     Gibelin et al. (2008) : New biomass reservoirs, and new method for allocation, 
!!     mortality and respiration.
!!
!!    EXTERNAL
!!    --------
!!    none
!!
!!    IMPLICIT ARGUMENTS
!!    ------------------
!!      
!!    none
!!
!!    REFERENCE
!!    ---------
!!
!! Calvet and Soussana (2001), "Modelling CO2-enrichment effects using an
!! interactive vegetation SVAT scheme", Agricultural and Forest Meteorology, Vol. 108
!! pp. 129-152
!! Gibelin et al. (2008), "Modelling energy and CO2 fluxes with an interactive vegetation 
!! land surface model - Evaluation at high and middle latitudes", 
!! Agricultural and Forest Meteorology, Vol. 148 , pp. 1611-1628
!!      
!!    AUTHOR
!!    ------
!!
!!      A.-L. Gibelin           * Meteo-France *
!!      (following Belair)
!!
!!    MODIFICATIONS
!!    -------------
!!      Original    27/01/03 
!!
!!      P Le Moigne  09/2005 : AGS modifs of L. Jarlan
!!      A.L. Gibelin 04/2009 : BIOMASS and RESP_BIOMASS arrays
!!      A.L. Gibelin 04/2009 : Suppress unused arguments
!!      A.L. Gibelin 04/2009 : Suppress unused modules and add ONLY
!!      A.L. Gibelin 04/2009 : adaptation to SURFEX environment
!!      A.   Barbu   01/2011 : modification of active biomass,leaf reservoir (see nitro_decline.f90)
!!      C.   Delire  04/2012 : spinup wood carbon
!!      R.   Alkama  04/2012 : 19 vegtype rather than 12
!!      B.   Decharme 05/2012: Optimization
!!                              ZCC_NITRO and ZBIOMASST_LIM in modd_co2v_par.F90
!!      C.   Delire   01/2014 : sapwood respiration from IBIS

!
!-------------------------------------------------------------------------------
!
!*       0.     DECLARATIONS
!               ------------
!
USE modd_csts,           ONLY : xpi, xday
USE modd_co2v_par,       ONLY : xpcco2, xcc_nit, xca_nit, xmc, &
                                xmco2, xcc_nitro, xbiomasst_lim 
USE modd_data_cover_par, ONLY : nvt_tebd, nvt_bone, nvt_trbe, nvt_trbd, &
                                nvt_tebe, nvt_tene, nvt_bobd, nvt_bond, &
                                nvt_shrb
!
USE yomhook   ,ONLY : lhook,   dr_hook
USE parkind1  ,ONLY : jprb
!
IMPLICIT NONE
!
!*      0.1    declarations of arguments
!
 CHARACTER(LEN=*),     INTENT(IN) :: hphoto           ! Kind of photosynthesis
!                                                    ! 'NON'
!                                                    ! 'AGS'
!                                                    ! 'LAI'
!                                                    ! 'AST'
!                                                    ! 'LST'
 CHARACTER(LEN=3),     INTENT(IN) :: hrespsl          ! Soil Respiration
!                                                    ! 'DEF' = Norman 1992
!                                                    ! 'PRM' = Rivalland PhD Thesis (2003)
!                                                    ! 'CNT' = CENTURY model (Gibelin 2008)
LOGICAL,              INTENT(IN) :: otr_ml           ! new TR
INTEGER, INTENT(IN)              :: kspinw           ! wood spinup
!
REAL,   DIMENSION(:), INTENT(IN) :: pbslai_nitro     ! ratio of biomass to LAI
REAL,   DIMENSION(:), INTENT(IN) :: psefold          ! e-folding time for senescence (s)
REAL,   DIMENSION(:), INTENT(IN) :: pgmes            ! mesophyll conductance (m s-1)
REAL,   DIMENSION(:), INTENT(IN) :: panmax           ! maximum photosynthesis rate
REAL,   DIMENSION(:), INTENT(IN) :: panday           ! daily net CO2 accumulation
REAL,   DIMENSION(:), INTENT(IN) :: plat             ! latitude of each grid point
REAL,   DIMENSION(:), INTENT(IN) :: plaimin          ! minimum LAI
REAL, DIMENSION(:,:), INTENT(IN) :: pvegtype         ! fraction of each vegetation
REAL,   DIMENSION(:), INTENT(IN) :: ptau_wood        ! residence time in wood (s)
REAL,   DIMENSION(:), INTENT(IN) :: plai             ! leaf area index (LAI) 
!
REAL,   DIMENSION(:), INTENT(INOUT) :: panfm         ! maximum leaf assimilation
REAL, DIMENSION(:,:), INTENT(INOUT) :: pbiomass      ! biomass reservoirs
REAL, DIMENSION(:,:), INTENT(INOUT) :: presp_biomass ! cumulated daily biomass respiration (kgDM m-2 day-1)
!
REAL,   DIMENSION(:), INTENT(OUT)   :: pbiomass_leaf ! temporary leaf biomass
REAL, DIMENSION(:,:), INTENT(OUT)   :: pincrease     ! increment of biomass
REAL, DIMENSION(:,:), INTENT(OUT)   :: pturnover     ! biomass turnover going into litter (gC m-2 s-1)
!
!*      0.2    declarations of local variables
!
REAL                            :: zbmcoef
REAL,    DIMENSION(SIZE(PLAI))  :: zxsefold        ! e-folding time for senescence corrected (days)
REAL,    DIMENSION(SIZE(PLAI))  :: zlaib_nitro     ! LAI correction parameter used in sefold calculation
REAL,    DIMENSION(SIZE(PLAI))  :: zassim          ! assimilation
REAL,    DIMENSION(SIZE(PLAI))  :: zbiomasst       ! leaf + active structural biomass
!
REAL, DIMENSION(SIZE(PLAI),SIZE(PBIOMASS,2))  :: zincrease
REAL, DIMENSION(SIZE(PLAI),SIZE(PBIOMASS,2))  :: zbiomass      ! temporary biomass reservoirs
REAL, DIMENSION(SIZE(PLAI),SIZE(PBIOMASS,2))  :: zdecline      ! biomass decline (storage+mortality) (kgDM m-2 day-1)
REAL, DIMENSION(SIZE(PLAI),SIZE(PBIOMASS,2))  :: zstorage      ! storage (part of decline kgDM m-2 day-1)
REAL, DIMENSION(SIZE(PLAI))                   :: zmort_leaf    ! leaf mortality
!
REAL, DIMENSION(SIZE(PLAI))                   :: zwork,zresp
LOGICAL, DIMENSION(SIZE(PLAI))                :: gmask_assim
LOGICAL, DIMENSION(SIZE(PLAI))                :: gwoody
!
REAL(KIND=JPRB) :: zhook_handle
!
INTEGER :: jspin, ji, ini
!
! correspondence between array indices and biomass compartments
! LEAF = 1
! STRUCT_ACT = 2
! STRUCT_PAS = 3
! STRUCT_BELOW = 4
! WOOD_ABOVE = 5
! WOOD_BELOW = 6
!
!-------------------------------------------------------------------------------
!
! 1 - Initialisations
! -------------------
!
IF (lhook) CALL dr_hook('NITRO_DECLINE',0,zhook_handle)
!
ini = SIZE(plai)
!
zxsefold(:)         = 0.0
zlaib_nitro(:)      = 0.0
zbiomasst(:)        = 0.0
zassim(:)           = 0.0
zbiomass(:,:)       = 0.0
zdecline(:,:)       = 0.0
zincrease(:,:)      = 0.0
zstorage(:,:)       = 0.0
zmort_leaf(:)       = 0.0
!---------------------------------------------------
!
zbmcoef     = xmc/(xmco2*xpcco2)
!
!-----------------------------------------------------------------
!avoid possible but unlikely negative values for biomass:        
!
pbiomass(:,1) = max(pbiomass(:,1),0.0)
!
! current leaf biomass value:
!
pbiomass_leaf(:) = pbiomass(:,1)
!
!-------------------------------------------------------------------------------
!
! Once a day (at midnight),repartition of net assimilation and mortality 
! into different biomass compartments.
!
! 2 - Evolution of leaf biomass and senescence calculations
! ---------------------------------------------------------
!
! coef c for biomass in kg/m2 now in modd_co2v_par.F90 (XCC_NITRO)
!
! LAI correction for shadow effect
IF (otr_ml) THEN
  zlaib_nitro(:) = 5.30
ELSE
  zlaib_nitro(:) = max( 5.76-0.64*atan(abs(plat(:))*xpi/180.),3.8 )
ENDIF
!
!
! leaf life expectancy
!
zwork(:) = 0.0
WHERE(pgmes(:)>0.0)
      zwork(:) = 0.321*log(pgmes(:)*1000.)
      zwork(:) = exp(zwork(:))*plai(:)/zlaib_nitro(:)
ENDWHERE
! before optimization
!ZXSEFOLD(:)= PSEFOLD(:) * MAX(((PGMES(:)*1000.)**0.321)*PLAI(:)/ZLAIB_NITRO(:), 1.) * ...
zxsefold(:) = psefold(:) * max(1.0,zwork(:)) * min(1.0,panfm(:)/panmax(:)) / xday
!
! avoid possible but unlikely division by zero
!
zxsefold(:) = max(1.0e-8,zxsefold(:))
!
! limitation of leaf life expectancy
!
! OLD   ZXSEFOLD(:) = MAX(5.,ZXSEFOLD(:))
! Following Marita's work limitation of the senesence
zxsefold(:) = max(psefold(:)/xday/10.0,zxsefold(:))
!
! senesence of active biomass
!
zdecline(:,1) = min(pbiomass_leaf(:)-plaimin(:)*pbslai_nitro(:), &
                    pbiomass_leaf(:)*(1.0-exp(-1.0/zxsefold(:))))
!
! avoid negative values due to computation precision
!
zdecline(:,1) = max(zdecline(:,1),0.0)
!
! current leaf biomass with assimilation and senescence
!
pbiomass_leaf(:) = pbiomass_leaf(:) - zdecline(:,1)
!
! daily active biomass assimilation
!
zassim(:) = panday(:)*zbmcoef
!
!-------------------------------------------------------------------------------
!
! 3 - Evolution of active structural biomass
! ------------------------------------------
!
zwork(:) = 0.0
WHERE(pbiomass_leaf(:)>0.0)
      zwork(:) = (1.0/(1.0-xca_nit))*log(pbiomass_leaf(:)/xcc_nitro)
      zwork(:) = exp(zwork(:))
ENDWHERE
!
WHERE (zassim(:) >= zdecline(:,1))
  !
  ! 3.1 - Growing phase : plant nitrogen decline theory
  !
  ! the growth allometric law is applied
  ! repartition of total biomass    
  !
  !before optimization
  !ZBIOMASST(:)= MAX(PBIOMASS_LEAF(:), (PBIOMASS_LEAF(:)/XCC_NITRO)**(1.0/(1.0-XCA_NIT)))  
  zbiomasst(:) = max(pbiomass_leaf(:), zwork(:))
  !
  ! active structural biomass increment and storage
  !
  zbiomass(:,2)  = zbiomasst(:)  - pbiomass_leaf(:)
  zdecline(:,2)  = zbiomass(:,2) * (1.0-exp(-1.0*xday/psefold(:)))
  zstorage(:,1)  = zbiomass(:,2) - pbiomass(:,2) + zdecline(:,2) + presp_biomass(:,2)
  !
ELSE WHERE
  !
  ! 3.2 - Senescence phase
  !
  ! the active structural biomass dies exponentially at the lowest rate
  !
  zstorage(:,1) = 0.0
  zdecline(:,2) = pbiomass(:,2) * (1.0-exp(-1.0*xday/psefold(:)))
  zbiomass(:,2) = pbiomass(:,2) - zdecline(:,2) - presp_biomass(:,2)
  !
  !  Avoid negative values of biomass
  !  No test on ZDECLINE(:,2) as it is not used after, or recalculated
  !  No test on PRESP_BIOMASS(:,2) as it should be smaller than PBIOMASS(:,2)
  !  otherwise there are irrealistic values of temperature     
  !
  zbiomass(:,2) = max(zbiomass(:,2),0.0)
  !
  zbiomasst(:) = pbiomass_leaf(:) + zbiomass(:,2)
  !
END WHERE
!
! 3.3 - Flow to the passive structural biomass: cut or growth after senescence
! Biomass is taken from active structural biomass, not from senescence of leaves
! 
zincrease(:,1) = zassim(:)
zincrease(:,2) = zstorage(:,1)
zincrease(:,3) = -min(zstorage(:,1),0.0)
!
zstorage(:,1) = max(0.0,zstorage(:,1))
!
! 3.4 - Mass conservation : leaf biomass sensecence must be >= structural storage
!
WHERE( zstorage(:,1) > zdecline(:,1))
  zdecline(:,2)    = pbiomass(:,2) * (1.0 - exp(-1.0*xday/psefold(:)))
  zbiomasst(:)     = pbiomass(:,1) + pbiomass(:,2) - zdecline(:,2) - presp_biomass(:,2)  
END WHERE
!
zwork(:) = 0.0
WHERE( zbiomasst(:) > 0.0)
      zwork(:) = (1.0-xca_nit)*log(zbiomasst(:))
      zwork(:) = exp(zwork(:))
ENDWHERE
!
WHERE( zstorage(:,1) > zdecline(:,1))
  !   
  !before optimization
  !PBIOMASS_LEAF(:)= ZCC_NITRO * (ZBIOMASST(:)**(1.0-XCA_NIT))
  pbiomass_leaf(:) = xcc_nitro * zwork(:)
  zbiomass(:,2)    = zbiomasst(:)  - pbiomass_leaf(:)
  zdecline(:,1)    = pbiomass(:,1) - pbiomass_leaf(:)
  zstorage(:,1)    = zbiomass(:,2) - pbiomass(:,2) + zdecline(:,2) + presp_biomass(:,2)  
  !
  zincrease(:,2) = zstorage(:,1)
  !
END WHERE
!
!-------------------------------------------------------------------------------
!
! 4 - Evolution of other biomass pools and final calculations
! -----------------------------------------------------------
!
! 4.1 - Mortality of leaf biomass
!
zmort_leaf(:) = max(0.0, zdecline(:,1) - zstorage(:,1))
!
zbiomass(:,3) = pbiomass(:,3)
!
IF (hphoto=='NIT') THEN
  !
  ! senesence of deep-structural biomass
  !
  zdecline(:,3) = zbiomass(:,3)*(1.0-exp(-1.0*xday/psefold(:)))          
  !
  ! threshold value for leaf biomass and total above ground biomass in nitrogen
  ! dilution theory now in modd_co2v_par.F90 (XBIOMASST_LIM)
  !
  ! emergency deep structural biomass
  WHERE((zbiomasst(:) <= xbiomasst_lim) .AND. (zxsefold(:) > 1.0))
    zbiomass(:,3) = zbiomass(:,3) + zmort_leaf(:)
  END WHERE
  !
ELSEIF (hphoto=='NCB') THEN
  !
  gwoody = (pvegtype(:,nvt_tebd)+pvegtype(:,nvt_bone)+pvegtype(:,nvt_trbe)+ &
            pvegtype(:,nvt_trbd)+pvegtype(:,nvt_tebe)+pvegtype(:,nvt_tene)+ &
            pvegtype(:,nvt_bobd)+pvegtype(:,nvt_bond)+pvegtype(:,nvt_shrb) >= 0.5)
  !
  ! 4.2 - Evolution of the other reservoirs
  ! 4.2.1 - senesence, avoiding negative values of biomass
  !
  zdecline(:,3) = min(pbiomass(:,3)*(1.0-exp(-1.0*xday/(psefold(:)/4.))), &
                      pbiomass(:,3)-presp_biomass(:,3))            
  zdecline(:,4) = min(pbiomass(:,4)*(1.0-exp(-1.0*xday/psefold(:))), &
                      pbiomass(:,4)-presp_biomass(:,4))
  !
  WHERE (gwoody(:))
    ! Woody
    zdecline(:,5) = min(pbiomass(:,5)*(1.0-exp(-1.0*xday/ptau_wood(:))), &
                      pbiomass(:,5)-presp_biomass(:,5))
    zdecline(:,6) = pbiomass(:,6)*(1.0-exp(-1.0*xday/ptau_wood(:)))
  ELSEWHERE
    ! Herbaceous
    zdecline(:,5) = 0.
    zdecline(:,6) = 0.
  END WHERE
  !
  ! 4.2.2 - storage (part of decline used as input for other reservoirs)
  !
  gmask_assim(:)=(zassim(:) >= zdecline(:,1))
  !
  WHERE (gmask_assim(:))
    !
    ! Remaining mortality is stored in roots.
    zincrease(:,4)   = zmort_leaf(:)
    !      
    ! Growing phase, all leaf decline is used as storage.
    zstorage(:,1)    = zstorage(:,1) + zincrease(:,4)
    zmort_leaf(:)    = zmort_leaf(:) - zincrease(:,4)
    !      
    zstorage(:,2)    = zdecline(:,2)
    zstorage(:,3)    = zdecline(:,3)
    !   
  ELSEWHERE
    !
    ! Senescence, a part of mortality is stored in roots, limited by assimilation rate.
    zincrease(:,4)   = min(max(0.5*zassim(:),0.) , 0.5*zmort_leaf(:))
    !
    zstorage(:,1)    = zstorage(:,1) + zincrease(:,4)
    zmort_leaf(:)    = zmort_leaf(:) - zincrease(:,4)
    !   
  END WHERE
  !
  WHERE(gmask_assim(:).AND.gwoody(:))
      ! Woody
      zstorage(:,4)  = zdecline(:,4)
      !
      zincrease(:,4) = zincrease(:,4) + 0.3* (zstorage(:,2) + zstorage(:,3))
      zincrease(:,5) =                  0.7* (zstorage(:,2) + zstorage(:,3))
      zincrease(:,6) = zstorage(:,4)
      !
  ELSEWHERE(gmask_assim(:).AND..NOT.gwoody(:))
      ! Herbaceous
      zstorage(:,4)  = 0.
      !
      zincrease(:,4) = zincrease(:,4) + zstorage(:,2) + zstorage(:,3)
      !
  END WHERE
  !
  WHERE (.NOT.gmask_assim(:).AND.gwoody(:))
      ! Woody
      ! Senescence, only a part of decline is used as storage
      zstorage(:,2)  = 0.5*zdecline(:,2)
      zstorage(:,3)  = 0.5*zdecline(:,3)
      zstorage(:,4)  = 0.5*zdecline(:,4)
      !
      zincrease(:,5) = zstorage(:,2) + zstorage(:,3)
      zincrease(:,6) = zstorage(:,4)
      !
  ELSEWHERE(.NOT.gmask_assim(:).AND..NOT.gwoody(:))
      !  Herbaceous
      ! Senescence, no storage
      zstorage(:,2)  = 0.
      zstorage(:,3)  = 0.
      zstorage(:,4)  = 0.
      !
  END WHERE
  !
  zstorage(:,5) = 0.
  zstorage(:,6) = 0.
  !
  ! 4.2.3 - mortality (senescence - storage) and turnover
  !
  IF (hrespsl=='CNT') THEN
    pturnover(:,1) = zmort_leaf(:)*1000.*xpcco2/xday
    pturnover(:,2) = (zdecline(:,2) - zstorage(:,2))*1000.*xpcco2/xday
    pturnover(:,3) = (zdecline(:,3) - zstorage(:,3))*1000.*xpcco2/xday
    pturnover(:,4) = (zdecline(:,4) - zstorage(:,4))*1000.*xpcco2/xday
    pturnover(:,5) = (zdecline(:,5) - zstorage(:,5))*1000.*xpcco2/xday
    pturnover(:,6) = (zdecline(:,6) - zstorage(:,6))*1000.*xpcco2/xday
  ENDIF
  !
ENDIF
!
!
! 4.3 - Re-initialisations for next time step
!
zbiomass(:,3) = zbiomass(:,3) + zincrease(:,3) - zdecline(:,3) - presp_biomass(:,3)
!
! Add net accumulated CO2 assimilation 
pbiomass_leaf(:) = pbiomass_leaf(:) + zassim(:)
!
! re-initialisation of biomass compartments values: X(day) <-- X(day-1)
pbiomass(:,1) = pbiomass_leaf(:)
pbiomass(:,2) = zbiomass(:,2)
pbiomass(:,3) = zbiomass(:,3)
!
! re-initialisation of respiration and assimilation terms
presp_biomass(:,2) = 0.0
presp_biomass(:,3) = 0.0
panfm(:) = 0.0
!
!
! 4.2.4 - evolution of reservoirs
!
IF (hphoto=='NIT') THEN
  !
  pbiomass(:,3) = max(pbiomass(:,3),0.0)
  !
ELSEIF (hphoto=='NCB') THEN
  !
  zbiomass(:,4) = pbiomass(:,4) + zincrease(:,4) - zdecline(:,4) - presp_biomass(:,4)
  !
!
  zbiomass(:,5) = pbiomass(:,5)
  zbiomass(:,6) = pbiomass(:,6)
  zresp(:) = presp_biomass(:,5)
!
  DO jspin = 1, kspinw
    DO ji = 1,ini
       IF(gwoody(ji))THEN
         !Woody
         zbiomass(ji,5) = zbiomass(ji,5) + zincrease(ji,5) - zdecline(ji,5) - zresp(ji)
         zbiomass(ji,6) = zbiomass(ji,6) + zincrease(ji,6) - zdecline(ji,6)
         zdecline(ji,5) = zbiomass(ji,5)*(1.0-exp((-1.0*xday)/ptau_wood(ji)))
         zdecline(ji,6) = zbiomass(ji,6)*(1.0-exp((-1.0*xday)/ptau_wood(ji)))
         IF (pbiomass(ji,5) .gt. 0.0) zresp(ji) = presp_biomass(ji,5)/pbiomass(ji,5) * zbiomass(ji,5)  
       ELSE   
         !Herbaceous
         zbiomass(ji,5) = 0.
         zbiomass(ji,6) = 0.
       ENDIF
    ENDDO
  ENDDO
!
  pbiomass(:,4) = zbiomass(:,4)
  pbiomass(:,5) = zbiomass(:,5)
  pbiomass(:,6) = zbiomass(:,6)
  !
  presp_biomass(:,4) = 0.0
  presp_biomass(:,5) = 0.0
  !
  pincrease(:,:) = zincrease(:,:)
!  
ENDIF
!
IF (lhook) CALL dr_hook('NITRO_DECLINE',1,zhook_handle)
!
END SUBROUTINE nitro_decline