! Light competition ! ! If canopy is almost closed (fpc > fpc_crit), then trees outcompete grasses. ! fpc_crit is normally fpc_crit. ! Here, fpc ("foilage protected cover") also takes into account the minimum fraction ! of space covered by trees through branches etc. This is done to prevent strong relative ! changes of fpc from one day to another for deciduous trees at the beginning of their ! growing season, which would yield to strong cutbacks (see 3.2.1.1.2) ! No competition between woody pfts (height of individuals is not considered) ! ! Exception: when one woody pft is overwhelming (i.e. fpc > fpc_crit). In that ! case, eliminate all other woody pfts and reduce dominant pft to fpc_crit. ! Age of individuals is not considered. In reality, light competition would more ! easily kill young individuals, thus increasing the mean age of the stand. ! Exclude agricultural pfts from competition ! ! SZ: added light competition for the static case if the mortality is not ! assumed to be constant. ! other modifs: ! -1 FPC is now always calculated from lm_lastyearmax*sla, since the aim of this DGVM is ! to represent community ecology effects; seasonal variations in establishment related to phenology ! may be relevant, but beyond the scope of a 1st generation DGVM ! -2 problem, if agriculture is present, fpc can never reach 1.0 since natural veget_max < 1.0. To ! correct for this, ind must be recalculated to correspond to the natural density... ! since ind is 1/m2 grid cell, this can be achived by dividing ind by the agricultural fraction ! ! $Header: /home/ssipsl/CVSREP/ORCHIDEE/src_stomate/lpj_light.f90,v 1.8 2009/01/06 15:01:25 ssipsl Exp $ ! IPSL (2006) ! This software is governed by the CeCILL licence see ORCHIDEE/ORCHIDEE_CeCILL.LIC ! MODULE lpj_light ! modules used: USE ioipsl USE stomate_constants USE constantes_veg IMPLICIT NONE ! private & public routines PRIVATE PUBLIC light, light_clear ! first call LOGICAL, SAVE :: firstcall = .TRUE. CONTAINS SUBROUTINE light_clear firstcall=.TRUE. END SUBROUTINE light_clear SUBROUTINE light (npts, dt, & veget_max, fpc_max, PFTpresent, cn_ind, lai, maxfpc_lastyear, & lm_lastyearmax, ind, biomass, veget_lastlight, bm_to_litter, mortality) ! ! 0 declarations ! ! 0.1 input ! Domain size INTEGER(i_std), INTENT(in) :: npts ! Time step (days) REAL(r_std), INTENT(in) :: dt ! Is pft there LOGICAL, DIMENSION(npts,nvm), INTENT(in) :: PFTpresent ! crown area of individuals (m**2) REAL(r_std), DIMENSION(npts,nvm), INTENT(in) :: cn_ind ! leaf area index OF AN INDIVIDUAL PLANT REAL(r_std), DIMENSION(npts,nvm), INTENT(in) :: lai ! last year's maximum fpc for each natural PFT, on ground REAL(r_std), DIMENSION(npts,nvm), INTENT(in) :: maxfpc_lastyear ! last year's maximum leafmass for each natural PFT, on ground REAL(r_std), DIMENSION(npts,nvm), INTENT(in) :: lm_lastyearmax ! last year's maximum fpc for each natural PFT, on ground REAL(r_std), DIMENSION(npts,nvm), INTENT(in) :: veget_max ! last year's maximum fpc for each natural PFT, on ground REAL(r_std), DIMENSION(npts,nvm), INTENT(in) :: fpc_max ! 0.2 modified fields ! Number of individuals / m2 REAL(r_std), DIMENSION(npts,nvm), INTENT(inout) :: ind ! biomass (gC/(m**2 of ground)) REAL(r_std), DIMENSION(npts,nvm,nparts), INTENT(inout) :: biomass ! Vegetation cover after last light competition REAL(r_std), DIMENSION(npts,nvm), INTENT(inout) :: veget_lastlight ! biomass taken away (gC/m**2) REAL(r_std), DIMENSION(npts,nvm,nparts), INTENT(inout) :: bm_to_litter ! fraction of individuals that died this time step REAL(r_std), DIMENSION(npts,nvm), INTENT(inout) :: mortality ! 0.3 local ! maximum total number of grass individuals in a closed canopy REAL(r_std), PARAMETER :: grass_mercy = 0.01 ! minimum fraction of trees that survive even in a closed canopy REAL(r_std), PARAMETER :: tree_mercy = 0.01 ! for diagnosis of fpc increase, compare today's fpc to last year's maximum (T) or ! to fpc of last time step (F)? LOGICAL, PARAMETER :: annual_increase = .TRUE. ! index INTEGER(i_std) :: i,j,k,m ! total natural fpc REAL(r_std), DIMENSION(npts) :: sumfpc ! fraction of natural vegetation at grid cell level REAL(r_std), DIMENSION(npts) :: fracnat ! total natural woody fpc REAL(r_std) :: sumfpc_wood ! change in total woody fpc REAL(r_std) :: sumdelta_fpc_wood ! maximum wood fpc REAL(r_std) :: maxfpc_wood ! which woody pft is maximum INTEGER(i_std) :: optpft_wood ! total natural grass fpc REAL(r_std) :: sumfpc_grass ! this year's foliage protected cover on natural part of the grid cell REAL(r_std), DIMENSION(npts,nvm) :: fpc_nat ! fpc change within last year REAL(r_std), DIMENSION(nvm) :: deltafpc ! Relative change of number of individuals for trees REAL(r_std) :: reduct ! Fraction of plants that survive REAL(r_std), DIMENSION(nvm) :: survive ! FPC for static mode REAL(r_std), DIMENSION(npts) :: fpc_real ! FPC mortality for static mode REAL(r_std), DIMENSION(npts) :: lai_ind ! number of grass PFTs present in the grid box ! INTEGER(i_std) :: num_grass ! New total grass fpc REAL(r_std) :: sumfpc_grass2 ! fraction of plants that dies each day (1/day) REAL(r_std), DIMENSION(npts,nvm) :: light_death ! Relative change of number of individuals for trees REAL(r_std) :: fpc_dec ! ========================================================================= IF (bavard.GE.3) WRITE(numout,*) 'Entering light' ! ! 1 first call ! IF ( firstcall ) THEN WRITE(numout,*) 'light:' WRITE(numout,*) ' > Maximum total number of grass individuals in' WRITE(numout,*) ' a closed canopy: ', grass_mercy WRITE(numout,*) ' > Minimum fraction of trees that survive even in' WRITE(numout,*) ' a closed canopy: ', tree_mercy WRITE(numout,*) ' > For trees, minimum fraction of crown area covered' WRITE(numout,*) ' (due to its branches etc.)', min_cover WRITE(numout,*) ' > for diagnosis of fpc increase, compare today''s fpc' IF ( annual_increase ) THEN WRITE(numout,*) ' to last year''s maximum.' ELSE WRITE(numout,*) ' to fpc of the last time step.' ENDIF firstcall = .FALSE. ENDIF IF (control%ok_dgvm) THEN ! ! 2 fpc characteristics ! ! 2.0 Only natural part of the grid cell: ! calculate fraction of natural and agricultural (1-fracnat) surface fracnat(:) = 1. DO j = 2,nvm IF ( .NOT. natural(j) ) THEN fracnat(:) = fracnat(:) - veget_max(:,j) ENDIF ENDDO ! ! 2.1 calculate fpc on natural part of grid cell. ! fpc_nat(:,:)=zero fpc_nat(:,ibare_sechiba)=un DO j = 2, nvm IF ( natural(j) ) THEN ! 2.1.1 natural PFTs IF ( tree(j) ) THEN ! 2.1.1.1 trees: minimum cover due to stems, branches etc. ! DO i = 1, npts ! IF (lai(i,j) == val_exp) THEN ! fpc_nat(i,j) = cn_ind(i,j) * ind(i,j) ! ELSE ! fpc_nat(i,j) = cn_ind(i,j) * ind(i,j) * & ! MAX( ( 1._r_std - exp( -lai(i,j) * ext_coeff(j) ) ), min_cover ) ! ENDIF ! ENDDO !NV : modif from SZ version : fpc is based on veget_max, not veget. WHERE(fracnat(:).GE.min_stomate) ! WHERE(LAI(:,j) == val_exp) ! fpc_nat(:,j) = cn_ind(:,j) * ind(:,j) / fracnat(:) ! ELSEWHERE ! fpc_nat(:,j) = cn_ind(:,j) * ind(:,j) / fracnat(:) * & ! MAX( ( 1._r_std - exp( - lm_lastyearmax(:,j) * sla(j) * ext_coeff(j) ) ), min_cover ) ! ENDWHERE fpc_nat(:,j) = cn_ind(:,j) * ind(:,j) / fracnat(:) ENDWHERE ELSE !NV : modif from SZ version : fpc is based on veget_max, not veget. WHERE(fracnat(:).GE.min_stomate) ! WHERE(LAI(:,j) == val_exp) ! fpc_nat(:,j) = cn_ind(:,j) * ind(:,j) / fracnat(:) ! ELSEWHERE ! fpc_nat(:,j) = cn_ind(:,j) * ind(:,j) / fracnat(:) * & ! ( 1._r_std - exp( - lm_lastyearmax(:,j) * sla(j) * ext_coeff(j) ) ) ! ENDWHERE fpc_nat(:,j) = cn_ind(:,j) * ind(:,j) / fracnat(:) ENDWHERE !!$ ! 2.1.1.2 bare ground !!$ IF (j == ibare_sechiba) THEN !!$ fpc_nat(:,j) = cn_ind(:,j) * ind(:,j) !!$ !!$ ! 2.1.1.3 grasses !!$ ELSE !!$ DO i = 1, npts !!$ IF (lai(i,j) == val_exp) THEN !!$ fpc_nat(i,j) = cn_ind(i,j) * ind(i,j) !!$ ELSE !!$ fpc_nat(i,j) = cn_ind(i,j) * ind(i,j) * & !!$ ( 1._r_std - exp( -lai(i,j) * ext_coeff(j) ) ) !!$ ENDIF !!$ ENDDO !!$ ENDIF ENDIF ! tree/grass ELSE ! 2.1.2 agricultural PFTs: not present on natural part fpc_nat(:,j) = zero ENDIF ! natural/agricultural ENDDO ! ! 2.2 sum natural fpc for every grid point ! sumfpc(:) = zero DO j = 2,nvm !SZ bug correction MERGE: need to subtract agricultural area! sumfpc(:) = sumfpc(:) + fpc_nat(:,j) ENDDO ! ! 3 Light competition ! light_death(:,:) = zero DO i = 1, npts ! SZ why this loop and not a vector statement ? ! Only if vegetation cover is dense IF ( sumfpc(i) .GT. fpc_crit ) THEN ! fpc change for each pft ! There are two possibilities: either we compare today's fpc with the fpc after the last ! time step, or we compare it to last year's maximum fpc of that PFT. In the first case, ! the fpc increase will be strong for seasonal PFTs at the beginning of the growing season. ! As for trees, the cutback is proportional to this increase, this means that seasonal trees ! will be disadvantaged compared to evergreen trees. In the original LPJ model, with its ! annual time step, the second method was used (this corresponds to annual_increase=.TRUE.) IF ( annual_increase ) THEN deltafpc(:) = MAX( (fpc_nat(i,:)-maxfpc_lastyear(i,:)), 0._r_std ) ELSE deltafpc(:) = MAX( (fpc_nat(i,:)-veget_lastlight(i,:)), 0._r_std ) ENDIF ! default: survive survive(:) = 1.0 ! ! 3.1 determine some characteristics of the fpc distribution ! sumfpc_wood = zero sumdelta_fpc_wood = zero maxfpc_wood = zero optpft_wood = 0 sumfpc_grass = zero ! num_grass = 0 DO j = 2,nvm ! only natural pfts IF ( natural(j) ) THEN IF ( tree(j) ) THEN ! trees ! total woody fpc sumfpc_wood = sumfpc_wood + fpc_nat(i,j) ! how much did the woody fpc increase sumdelta_fpc_wood = sumdelta_fpc_wood + deltafpc(j) ! which woody pft is preponderant IF ( fpc_nat(i,j) .GT. maxfpc_wood ) THEN optpft_wood = j maxfpc_wood = fpc_nat(i,j) ENDIF ELSE ! grasses ! total (natural) grass fpc sumfpc_grass = sumfpc_grass + fpc_nat(i,j) ! number of grass PFTs present in the grid box ! IF ( PFTpresent(i,j) ) THEN ! num_grass = num_grass + 1 ! ENDIF ENDIF ! tree or grass ENDIF ! natural ENDDO ! loop over pfts ! ! 3.2 light competition: assume wood outcompetes grass ! !SZ !!$ IF (sumfpc_wood .GE. fpc_crit ) THEN ! ! 3.2.1 all allowed natural space is covered by wood: ! cut back trees to fpc_crit. ! Original DGVM: kill grasses. Modified: we let a very ! small fraction of grasses survive. ! DO j = 2,nvm ! only present and natural pfts compete IF ( PFTpresent(i,j) .AND. natural(j) ) THEN IF ( tree(j) ) THEN ! ! 3.2.1.1 tree ! ! no single woody pft is overwhelming ! (original DGVM: tree_mercy = 0.0 ) ! The reduction rate is proportional to the ratio deltafpc/fpc. IF (sumfpc_wood .GE. fpc_crit .AND. fpc_nat(i,j) .GT. min_stomate .AND. & sumdelta_fpc_wood .GT. min_stomate) THEN ! reduct = MIN( ( ( deltafpc(j)/sumdelta_fpc_wood * & ! (sumfpc_wood-fpc_crit) ) / fpc_nat(i,j) ), & ! ( 1._r_std - tree_mercy ) ) reduct = un - MIN((fpc_nat(i,j)-(sumfpc_wood-fpc_crit) & * deltafpc(j)/sumdelta_fpc_wood)/fpc_nat(i,j), un ) ELSE ! tree fpc didn't icrease or it started from nothing reduct = zero ENDIF survive(j) = un - reduct ELSE ! ! 3.2.1.2 grass: let a very small fraction survive (the sum of all ! grass individuals may make up a maximum cover of ! grass_mercy [for lai -> infinity]). ! In the original DGVM, grasses were killed in that case, ! corresponding to grass_mercy = 0. ! ! survive(j) = ( grass_mercy / REAL( num_grass,r_std ) ) / ind(i,j) ! survive(j) = MIN( 1._r_std, survive(j) ) IF(sumfpc_grass .GE. 1.0-MIN(fpc_crit,sumfpc_wood).AND. & sumfpc_grass.GE.min_stomate) THEN fpc_dec=(sumfpc_grass-1.+MIN(fpc_crit,sumfpc_wood))*fpc_nat(i,j)/sumfpc_grass reduct=fpc_dec ELSE reduct = zero ENDIF survive(j) = ( un - reduct ) ENDIF ! tree or grass ENDIF ! pft there and natural ENDDO ! loop over pfts !SZ !!$ ELSE !!$ !!$ ! !!$ ! 3.2.2 not too much wood so that grasses can subsist !!$ ! !!$ !!$ ! new total grass fpc !!$ sumfpc_grass2 = fpc_crit - sumfpc_wood !!$ !!$ DO j = 2,nvm !!$ !!$ ! only present and natural PFTs compete !!$ !!$ IF ( PFTpresent(i,j) .AND. natural(j) ) THEN !!$ !!$ IF ( tree(j) ) THEN !!$ !!$ ! no change for trees !!$ !!$ survive(j) = 1.0 !!$ !!$ ELSE !!$ !!$ ! grass: fractional loss is the same for all grasses !!$ !!$ IF ( sumfpc_grass .GT. min_stomate ) THEN !!$ survive(j) = sumfpc_grass2 / sumfpc_grass !!$ ELSE !!$ survive(j)= zero !!$ ENDIF !!$ !!$ ENDIF !!$ !!$ ENDIF ! pft there and natural !!$ !!$ ENDDO ! loop over pfts !!$ !!$ ENDIF ! sumfpc_wood > fpc_crit ! ! 3.3 update output variables ! DO j = 2,nvm IF ( PFTpresent(i,j) .AND. natural(j) ) THEN bm_to_litter(i,j,:) = bm_to_litter(i,j,:) + & biomass(i,j,:) * ( un - survive(j) ) biomass(i,j,:) = biomass(i,j,:) * survive(j) IF ( control%ok_dgvm ) THEN ind(i,j) = ind(i,j) * survive(j) ENDIF ! fraction of plants that dies each day. ! exact formulation: light_death(i,j) = un - survive(j) / dt light_death(i,j) = ( un - survive(j) ) / dt ENDIF ! pft there and natural ENDDO ! loop over pfts ENDIF ! sumfpc > fpc_crit ENDDO ! loop over grid points ! ! 4 recalculate fpc on natural part of grid cell (for next light competition) ! DO j = 2,nvm IF ( natural(j) ) THEN ! ! 4.1 natural PFTs ! IF ( tree(j) ) THEN ! 4.1.1 trees: minimum cover due to stems, branches etc. DO i = 1, npts !NVMODIF ! IF (lai(i,j) == val_exp) THEN ! veget_lastlight(i,j) = cn_ind(i,j) * ind(i,j) ! ELSE ! veget_lastlight(i,j) = & ! cn_ind(i,j) * ind(i,j) * & ! MAX( ( un - exp( -lai(i,j) * ext_coeff(j) ) ), min_cover ) ! ENDIF !! veget_lastlight(i,j) = cn_ind(i,j) * ind(i,j) IF (lai(i,j) == val_exp) THEN veget_lastlight(i,j) = cn_ind(i,j) * ind(i,j) ELSE veget_lastlight(i,j) = & cn_ind(i,j) * ind(i,j) * & MAX( ( un - EXP( - lm_lastyearmax(i,j) * sla(j) * ext_coeff(j) ) ), min_cover ) ENDIF ENDDO ELSE ! 4.1.2 grasses DO i = 1, npts !NVMODIF ! IF (lai(i,j) == val_exp) THEN ! veget_lastlight(i,j) = cn_ind(i,j) * ind(i,j) ! ELSE ! veget_lastlight(i,j) = cn_ind(i,j) * ind(i,j) * & ! ( un - exp( -lai(i,j) * ext_coeff(j) ) ) ! ENDIF !!veget_lastlight(i,j) = cn_ind(i,j) * ind(i,j) IF (lai(i,j) == val_exp) THEN veget_lastlight(i,j) = cn_ind(i,j) * ind(i,j) ELSE veget_lastlight(i,j) = cn_ind(i,j) * ind(i,j) * & ( un - exp( - lm_lastyearmax(i,j) * sla(j) * ext_coeff(j) ) ) ENDIF ENDDO ENDIF ! tree/grass ELSE ! ! 4.2 agricultural PFTs: not present on natural part ! veget_lastlight(:,j) = zero ENDIF ! natural/agricultural ENDDO ELSE ! static light_death(:,:)=0.0 DO j = 2, nvm IF ( natural(j) ) THEN ! 2.1.1 natural PFTs, in the one PFT only case there needs to be no special case for grasses, ! neither a redistribution of mortality (delta fpc) WHERE( ind(:,j)*cn_ind(:,j) .GT. min_stomate ) lai_ind(:)=sla(j) * lm_lastyearmax(:,j) / ( ind(:,j) * cn_ind(:,j) ) ELSEWHERE lai_ind(:)=zero ENDWHERE fpc_nat(:,j) = cn_ind(:,j) * ind(:,j) * & MAX( ( 1._r_std - exp( - ext_coeff(j) * lai_ind(:) ) ), min_cover ) WHERE(fpc_nat(:,j).GT.fpc_max(:,j)) light_death(:,j)=MIN(1.0,1.0-fpc_max(:,j)/fpc_nat(:,j)) ENDWHERE DO k=1,nparts bm_to_litter(:,j,k)=bm_to_litter(:,j,k)+light_death(:,j)*biomass(:,j,k) biomass(:,j,k)=biomass(:,j,k)-light_death(:,j)*biomass(:,j,k) ENDDO ind(:,j)=ind(:,j)-light_death(:,j)*ind(:,j) ! if (j==10) print *,'ind10bis=',ind(:,j),light_death(:,j)*ind(:,j) ENDIF ENDDO light_death(:,:)=light_death(:,:)/dt ENDIF ! ! 5 history ! CALL histwrite (hist_id_stomate, 'LIGHT_DEATH', itime, & light_death, npts*nvm, horipft_index) IF (bavard.GE.4) WRITE(numout,*) 'Leaving light' END SUBROUTINE light END MODULE lpj_light