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(IO, PK, PEK, OWOOD, PBSLAI_NITRO, PLAT, PBIOMASS_LEAF )  
!
!   ###############################################################
!!**  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_isba_options_n, ONLY : isba_options_t
USE modd_isba_n, ONLY : isba_p_t, isba_pe_t
!
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
!
TYPE(isba_options_t), INTENT(INOUT) :: IO
TYPE(isba_p_t), INTENT(INOUT) :: PK
TYPE(isba_pe_t), INTENT(INOUT) :: PEK
!
!*      0.1    declarations of arguments
!
LOGICAL, DIMENSION(:), INTENT(IN) :: OWOOD
!
REAL,   DIMENSION(:), INTENT(IN) :: PBSLAI_NITRO     ! ratio of biomass to LAI
REAL,   DIMENSION(:), INTENT(IN) :: PLAT             ! latitude of each grid point
!
REAL,   DIMENSION(:), INTENT(OUT)   :: PBIOMASS_LEAF ! temporary leaf biomass
!
!*      0.2    declarations of local variables
!
REAL                            :: ZBMCOEF
REAL,    DIMENSION(SIZE(PEK%XLAI,1))  :: ZXSEFOLD        ! e-folding time for senescence corrected (days)
REAL,    DIMENSION(SIZE(PEK%XLAI,1))  :: ZLAIB_NITRO     ! LAI correction parameter used in sefold calculation
REAL,    DIMENSION(SIZE(PEK%XLAI,1))  :: ZASSIM          ! assimilation
REAL,    DIMENSION(SIZE(PEK%XLAI,1))  :: ZBIOMASST       ! leaf + active structural biomass
!
REAL, DIMENSION(SIZE(PEK%XLAI,1),SIZE(PEK%XBIOMASS,2))  :: ZINCREASE
REAL, DIMENSION(SIZE(PEK%XLAI,1),SIZE(PEK%XBIOMASS,2))  :: ZBIOMASS      ! temporary biomass reservoirs
REAL, DIMENSION(SIZE(PEK%XLAI,1),SIZE(PEK%XBIOMASS,2))  :: ZDECLINE      ! biomass decline (storage+mortality) (kgDM m-2 day-1)
REAL, DIMENSION(SIZE(PEK%XLAI,1),SIZE(PEK%XBIOMASS,2))  :: ZSTORAGE      ! storage (part of decline kgDM m-2 day-1)
REAL, DIMENSION(SIZE(PEK%XLAI,1))                   :: ZMORT_LEAF    ! leaf mortality
!
REAL, DIMENSION(SIZE(PEK%XLAI,1))                   :: ZWORK,ZRESP
LOGICAL, DIMENSION(SIZE(PEK%XLAI,1))                :: GMASK_ASSIM
!
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(pek%XLAI,1)
!
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:        
!
pek%XBIOMASS(:,1) = max(pek%XBIOMASS(:,1),0.0)
!
! current leaf biomass value:
!
pbiomass_leaf(:) = pek%XBIOMASS(:,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 (io%LTR_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(pek%XGMES(:)>0.0)
      zwork(:) = 0.321*log(pek%XGMES(:)*1000.)
      zwork(:) = exp(zwork(:))*pek%XLAI(:)/zlaib_nitro(:)
ENDWHERE
! before optimization
!ZXSEFOLD(:)= PEK%XSEFOLD(:) * MAX(((PEK%XGMES(:)*1000.)**0.321)*PEK%XLAI(:)/ZLAIB_NITRO(:), 1.) * ...
zxsefold(:) = pek%XSEFOLD(:) * max(1.0,zwork(:)) * min(1.0,pek%XANFM(:)/pk%XANMAX(:)) / 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(pek%XSEFOLD(:)/xday/10.0,zxsefold(:))
!
! senesence of active biomass
!
zdecline(:,1) = min(pbiomass_leaf(:)-pek%XLAIMIN(:)*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(:) = pek%XANDAY(:)*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/pek%XSEFOLD(:)))
  zstorage(:,1)  = zbiomass(:,2) - pek%XBIOMASS(:,2) + zdecline(:,2) + pek%XRESP_BIOMASS(:,2)
  !
ELSE WHERE
  !
  ! 3.2 - Senescence phase
  !
  ! the active structural biomass dies exponentially at the lowest rate
  !
  zstorage(:,1) = 0.0
  zdecline(:,2) = pek%XBIOMASS(:,2) * (1.0-exp(-1.0*xday/pek%XSEFOLD(:)))
  zbiomass(:,2) = pek%XBIOMASS(:,2) - zdecline(:,2) - pek%XRESP_BIOMASS(:,2)
  !
  !  Avoid negative values of biomass
  !  No test on ZDECLINE(:,2) as it is not used after, or recalculated
  !  No test on PEK%XRESP_BIOMASS(:,2,1) as it should be smaller than PEK%XBIOMASS(:,2,1)
  !  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)    = pek%XBIOMASS(:,2) * (1.0 - exp(-1.0*xday/pek%XSEFOLD(:)))
  zbiomasst(:)     = pek%XBIOMASS(:,1) + pek%XBIOMASS(:,2) - zdecline(:,2) - pek%XRESP_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)    = pek%XBIOMASS(:,1) - pbiomass_leaf(:)
  zstorage(:,1)    = zbiomass(:,2) - pek%XBIOMASS(:,2) + zdecline(:,2) + pek%XRESP_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) = pek%XBIOMASS(:,3)
!
IF (io%CPHOTO=='NIT') THEN
  !
  ! senesence of deep-structural biomass
  !
  zdecline(:,3) = zbiomass(:,3)*(1.0-exp(-1.0*xday/pek%XSEFOLD(:)))          
  !
  ! 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 (io%CPHOTO=='NCB') THEN
  !
  ! 4.2 - Evolution of the other reservoirs
  ! 4.2.1 - senesence, avoiding negative values of biomass
  !
  zdecline(:,3) = min(pek%XBIOMASS(:,3)*(1.0-exp(-1.0*xday/(pek%XSEFOLD(:)/4.))), &
                      pek%XBIOMASS(:,3)-pek%XRESP_BIOMASS(:,3))            
  zdecline(:,4) = min(pek%XBIOMASS(:,4)*(1.0-exp(-1.0*xday/pek%XSEFOLD(:))), &
                      pek%XBIOMASS(:,4)-pek%XRESP_BIOMASS(:,4))
  !
  WHERE (owood(:))
    ! Woody
    zdecline(:,5) = min(pek%XBIOMASS(:,5)*(1.0-exp(-1.0*xday/pk%XTAU_WOOD(:))), &
                        pek%XBIOMASS(:,5)-pek%XRESP_BIOMASS(:,5))
    zdecline(:,6) = pek%XBIOMASS(:,6)*(1.0-exp(-1.0*xday/pk%XTAU_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.owood(:))
      ! 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.owood(:))
      ! Herbaceous
      zstorage(:,4)  = 0.
      !
      zincrease(:,4) = zincrease(:,4) + zstorage(:,2) + zstorage(:,3)
      !
  END WHERE
  !
  WHERE (.NOT.gmask_assim(:).AND.owood(:))
      ! 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.owood(:))
      !  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 (io%CRESPSL=='CNT') THEN
    pk%XTURNOVER(:,1) = zmort_leaf(:)*1000.*xpcco2/xday
    pk%XTURNOVER(:,2) = (zdecline(:,2) - zstorage(:,2))*1000.*xpcco2/xday
    pk%XTURNOVER(:,3) = (zdecline(:,3) - zstorage(:,3))*1000.*xpcco2/xday
    pk%XTURNOVER(:,4) = (zdecline(:,4) - zstorage(:,4))*1000.*xpcco2/xday
    pk%XTURNOVER(:,5) = (zdecline(:,5) - zstorage(:,5))*1000.*xpcco2/xday
    pk%XTURNOVER(:,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) - pek%XRESP_BIOMASS(:,3)
!
! Add net accumulated CO2 assimilation 
pbiomass_leaf(:) = pbiomass_leaf(:) + zassim(:)
!
! re-initialisation of biomass compartments values: X(day) <-- X(day-1)
pek%XBIOMASS(:,1) = pbiomass_leaf(:)
pek%XBIOMASS(:,2) = zbiomass(:,2)
pek%XBIOMASS(:,3) = zbiomass(:,3)
!
! re-initialisation of respiration and assimilation terms
pek%XRESP_BIOMASS(:,2) = 0.0
pek%XRESP_BIOMASS(:,3) = 0.0
pek%XANFM(:) = 0.0
!
!
! 4.2.4 - evolution of reservoirs
!
IF (io%CPHOTO=='NIT') THEN
  !
  pek%XBIOMASS(:,3) = max(pek%XBIOMASS(:,3),0.0)
  !
ELSEIF (io%CPHOTO=='NCB') THEN
  !
  zbiomass(:,4) = pek%XBIOMASS(:,4) + zincrease(:,4) - zdecline(:,4) - pek%XRESP_BIOMASS(:,4)
  !
!
  zbiomass(:,5) = pek%XBIOMASS(:,5)
  zbiomass(:,6) = pek%XBIOMASS(:,6)
  zresp(:) = pek%XRESP_BIOMASS(:,5)
!
  DO jspin = 1, io%NSPINW
    DO ji = 1,ini
       IF(owood(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)/pk%XTAU_WOOD(ji)))
         zdecline(ji,6) = zbiomass(ji,6)*(1.0-exp((-1.0*xday)/pk%XTAU_WOOD(ji)))
         IF (pek%XBIOMASS(ji,5).GT.0.0) &
                 zresp(ji) = pek%XRESP_BIOMASS(ji,5)/pek%XBIOMASS(ji,5) * zbiomass(ji,5)  
       ELSE   
         !Herbaceous
         zbiomass(ji,5) = 0.
         zbiomass(ji,6) = 0.
       ENDIF
    ENDDO
  ENDDO
!
  pek%XBIOMASS(:,4) = zbiomass(:,4)
  pek%XBIOMASS(:,5) = zbiomass(:,5)
  pek%XBIOMASS(:,6) = zbiomass(:,6)
  !
  pek%XRESP_BIOMASS(:,4) = 0.0
  pek%XRESP_BIOMASS(:,5) = 0.0
  !
  pk%XINCREASE(:,:) = zincrease(:,:)
 
ENDIF
!
IF (lhook) CALL dr_hook('NITRO_DECLINE',1,zhook_handle)
!
END SUBROUTINE nitro_decline