source: branches/publications/ORCHIDEE_2.2_r7266/ORCHIDEE/src_stomate/stomate_turnover.f90 @ 7541

Last change on this file since 7541 was 7541, checked in by fabienne.maignan, 2 years ago
  1. Zhang publication on coupling factor
File size: 43.2 KB
Line 
1! =================================================================================================================================
2! MODULE       : stomate_turnover.f90
3!
4! CONTACT      : orchidee-help _at_ listes.ipsl.fr
5!
6! LICENCE      : IPSL (2006)
7! This software is governed by the CeCILL licence see ORCHIDEE/ORCHIDEE_CeCILL.LIC
8!
9!>\BRIEF        This module manages the end of the growing season and calculates herbivory and turnover of leaves, fruits, fine roots.
10!!
11!!\n DESCRIPTION: This subroutine calculates leaf senescence due to climatic conditions or as a
12!! function of leaf age and new LAI, and subsequent turnover of the different plant biomass compartments (sections 1 to 6),
13!! herbivory (section 7), fruit turnover for trees (section 8) and sapwood conversion (section 9).
14!!
15!! RECENT CHANGE(S): None
16!!
17!! SVN          :
18!! $HeadURL: svn://forge.ipsl.jussieu.fr/orchidee/branches/ORCHIDEE_2_2/ORCHIDEE/src_stomate/stomate_turnover.f90 $
19!! $Date: 2017-10-18 11:15:06 +0200 (Wed, 18 Oct 2017) $
20!! $Revision: 4693 $
21!! \n
22!_ ================================================================================================================================
23
24MODULE stomate_turnover
25
26  ! modules used:
27  USE xios_orchidee
28  USE ioipsl_para
29  USE stomate_data
30  USE constantes
31  USE pft_parameters
32
33  IMPLICIT NONE
34
35  ! private & public routines
36
37  PRIVATE
38  PUBLIC turn, turn_clear
39
40  LOGICAL, SAVE                          :: firstcall_turnover = .TRUE.           !! first call (true/false)
41!$OMP THREADPRIVATE(firstcall_turnover)
42
43CONTAINS
44
45
46!! ================================================================================================================================
47!! SUBROUTINE   : turn_clear
48!!
49!>\BRIEF        Set flag ::firstcall_turnover to .TRUE., and therefore activate section 1
50!!              of subroutine turn which writes a message to the output.
51!!               
52!_ ================================================================================================================================
53
54  SUBROUTINE turn_clear
55    firstcall_turnover=.TRUE.
56  END SUBROUTINE turn_clear
57
58
59!! ================================================================================================================================
60!! SUBROUTINE    : turn
61!!
62!>\BRIEF         Calculate turnover of leaves, roots, fruits and sapwood due to aging or climatic
63!!               induced senescence. Calculate herbivory.
64!!
65!! DESCRIPTION : This subroutine determines the turnover of leaves and fine roots (and stems for grasses)
66!! and simulates following processes:
67!! 1. Mean leaf age is calculated from leaf ages of separate leaf age classes. Should actually
68!!    be recalculated at the end of the routine, but it does not change too fast. The mean leaf
69!!    age is calculated using the following equation:
70!!    \latexonly
71!!    \input{turnover_lma_update_eqn1.tex}
72!!    \endlatexonly
73!!    \n
74!! 2. Meteorological senescence: the detection of the end of the growing season and shedding
75!!    of leaves, fruits and fine roots due to unfavourable meteorological conditions.
76!!    The model distinguishes three different types of "climatic" leaf senescence, that do not
77!!    change the age structure: sensitivity to cold temperatures, to lack of water, or both.
78!!    If meteorological conditions are fulfilled, a flag ::senescence is set to TRUE. Note
79!!    that evergreen species do not experience climatic senescence.
80!!    Climatic senescence is triggered by sensitivity to cold temperatures where the critical
81!!    temperature for senescence is calculated using the following equation:
82!!    \latexonly
83!!    \input{turnover_temp_crit_eqn2.tex}
84!!    \endlatexonly
85!!    \n
86!!    Climatic senescence is triggered by sensitivity to lack of water availability where the
87!!    moisture availability critical level is calculated using the following equation:
88!!    \latexonly
89!!    \input{turnover_moist_crit_eqn3.tex}
90!!    \endlatexonly
91!!    \n
92!!    Climatic senescence is triggered by sensitivity to temperature or to lack of water where
93!!    critical temperature and moisture availability are calculated as above.\n
94!!    Trees in climatic senescence lose their fine roots at the same rate as they lose their leaves.
95!!    The rate of biomass loss of both fine roots and leaves is presribed through the equation:
96!!    \latexonly
97!!    \input{turnover_clim_senes_biomass_eqn4.tex}
98!!    \endlatexonly
99!!    \n
100!!    with ::leaffall(j) a PFT-dependent time constant which is given in
101!!    ::stomate_constants. In grasses, leaf senescence is extended to the whole plant
102!!    (all carbon pools) except to its carbohydrate reserve.   
103!! 3. Senescence due to aging: the loss of leaves, fruits and  biomass due to aging
104!!    At a certain age, leaves fall off, even if the climate would allow a green plant
105!!    all year round. Even if the meteorological conditions are favorable for leaf maintenance,
106!!    plants, and in particular, evergreen trees, have to renew their leaves simply because the
107!!    old leaves become inefficient. Roots, fruits (and stems for grasses) follow leaves.
108!!    The ??senescence?? rate varies with leaf age. Note that plant is not declared senescent
109!!    in this case (wchich is important for allocation: if the plant loses leaves because of
110!!    their age, it can renew them). The leaf turnover rate due to aging of leaves is calculated
111!!    using the following equation:
112!!    \latexonly
113!!    \input{turnover_age_senes_biomass_eqn5.tex}
114!!    \endlatexonly
115!!    \n
116!!    Drop all leaves if there is a very low leaf mass during senescence. After this, the biomass
117!!    of different carbon pools both for trees and grasses is set to zero and the mean leaf age
118!!    is reset to zero. Finally, the leaf fraction and leaf age of the different leaf age classes
119!!    is set to zero. For deciduous trees: next to leaves, also fruits and fine roots are dropped.
120!!    For grasses: all aboveground carbon pools, except the carbohydrate reserves are affected:
121!! 4. Update the leaf biomass, leaf age class fraction and the LAI
122!!    Older leaves will fall more frequently than younger leaves and therefore the leaf age
123!!    distribution needs to be recalculated after turnover. The fraction of biomass in each
124!!    leaf class is updated using the following equation:
125!!    \latexonly
126!!    \input{turnover_update_LeafAgeDistribution_eqn6.tex}
127!!    \endlatexonly
128!!    \n
129!! 5. Simulate herbivory activity and update leaf and fruits biomass. Herbivore activity
130!!    affects the biomass of leaves and fruits as well as stalks (only for grasses).
131!!    However, herbivores do not modify leaf age structure.
132!! 6. Calculates fruit turnover for trees. Trees simply lose their fruits with a time
133!!    constant ::tau_fruit(j), that is set to 90 days for all PFTs in ::stomate_constants
134!! 7. Convert sapwood to heartwood for trees and update heart and softwood above and
135!!    belowground biomass. Sapwood biomass is converted into heartwood biomass
136!!    with a time constant tau ::tau_sap(j) of 1 year. Note that this biomass conversion
137!!    is not added to "turnover" as the biomass is not lost. For the updated heartwood,
138!!    the sum of new heartwood above and new heartwood below after converting sapwood to
139!!    heartwood, is saved as ::hw_new(:). Creation of new heartwood decreases the age of
140!!    the plant ??carbon?? with a factor that is determined by: old heartwood ::hw_old(:)
141!!    divided by the new heartwood ::hw_new(:)
142!!
143!! RECENT CHANGE(S) : None
144!!
145!! MAIN OUTPUT VARIABLES: ::Biomass of leaves, fruits, fine roots and sapwood above (latter for grasses only),
146!!                        ::Update LAI, ::Update leaf age distribution with new leaf age class fraction
147!!
148!! REFERENCE(S) :
149!! - Krinner, G., N. Viovy, N. de Noblet-Ducoudre, J. Ogee, J. Polcher, P.
150!! Friedlingstein, P. Ciais, S. Sitch and I.C. Prentice (2005), A dynamic global
151!! vegetation model for studies of the coupled atmosphere-biosphere system, Global
152!! Biogeochemical Cycles, 19, doi:10.1029/2003GB002199.
153!! - McNaughton, S. J., M. Oesterheld, D. A. Frank and K. J. Williams (1989),
154!! Ecosystem-level patterns of primary productivity and herbivory in terrestrial habitats,
155!! Nature, 341, 142-144, 1989.
156!! - Sitch, S., C. Huntingford, N. Gedney, P. E. Levy, M. Lomas, S. L. Piao, , Betts, R., Ciais, P., Cox, P.,
157!! Friedlingstein, P., Jones, C. D., Prentice, I. C. and F. I. Woodward : Evaluation of the terrestrial carbon 
158!! cycle, future plant geography and climate-carbon cycle feedbacks using 5 dynamic global vegetation
159!! models (dgvms), Global Change Biology, 14(9), 2015–2039, 2008.
160!!
161!! FLOWCHART    :
162!! \latexonly
163!! \includegraphics[scale=0.5]{turnover_flowchart_1.png}
164!! \includegraphics[scale=0.5]{turnover_flowchart_2.png}
165!! \endlatexonly
166!! \n
167!_ ================================================================================================================================
168
169  SUBROUTINE turn (npts, dt, PFTpresent, &
170       herbivores, &
171       maxmoiavail_lastyear, minmoiavail_lastyear, &
172       moiavail_week, moiavail_month, t2m_longterm, t2m_month, t2m_week, veget_cov_max, &
173       gdd_from_growthinit, leaf_age, leaf_frac, age, lai, biomass, &
174       turnover, senescence,turnover_time)
175
176    !! 0. Variable and parameter declaration
177
178    !! 0.1 Input variables
179
180    INTEGER(i_std), INTENT(in)                                 :: npts                 !! Domain size - number of grid cells
181                                                                                       !! (unitless)
182    REAL(r_std), INTENT(in)                                    :: dt                   !! time step (dt_days)
183    LOGICAL, DIMENSION(npts,nvm), INTENT(in)                   :: PFTpresent           !! PFT exists (true/false)
184    REAL(r_std), DIMENSION(npts,nvm), INTENT(in)               :: herbivores           !! time constant of probability of a leaf to
185                                                                                       !! be eaten by a herbivore (days)
186    REAL(r_std), DIMENSION(npts,nvm), INTENT(in)               :: maxmoiavail_lastyear !! last year's maximum moisture availability
187                                                                                       !! (0-1, unitless)
188    REAL(r_std), DIMENSION(npts,nvm), INTENT(in)               :: minmoiavail_lastyear !! last year's minimum moisture availability
189                                                                                       !! (0-1, unitless)
190    REAL(r_std), DIMENSION(npts,nvm), INTENT(in)               :: moiavail_week        !! "weekly" moisture availability
191                                                                                       !! (0-1, unitless)
192    REAL(r_std), DIMENSION(npts,nvm), INTENT(in)               :: moiavail_month       !! "monthly" moisture availability
193                                                                                       !! (0-1, unitless)
194    REAL(r_std), DIMENSION(npts), INTENT(in)                   :: t2m_longterm         !! "long term" 2 meter reference
195                                                                                       !! temperatures (K)
196    REAL(r_std), DIMENSION(npts), INTENT(in)                   :: t2m_month            !! "monthly" 2-meter temperatures (K)
197    REAL(r_std), DIMENSION(npts), INTENT(in)                   :: t2m_week             !! "weekly" 2 meter temperatures (K)
198    REAL(r_std), DIMENSION(npts,nvm), INTENT(in)               :: veget_cov_max        !! "maximal" coverage fraction of a PFT (LAI
199                                                                                       !! -> infinity) on ground (unitless)
200    REAL(r_std), DIMENSION(npts,nvm), INTENT(in)               :: gdd_from_growthinit  !! gdd senescence for crop
201
202    !! 0.2 Output variables
203
204    REAL(r_std), DIMENSION(npts,nvm,nparts,nelements), INTENT(out) :: turnover         !! Turnover @tex ($gC m^{-2}$) @endtex
205    LOGICAL, DIMENSION(npts,nvm), INTENT(out)                  :: senescence           !! is the plant senescent? (true/false)
206                                                                                       !! (interesting only for deciduous trees:
207                                                                                       !! carbohydrate reserve)
208    !! 0.3 Modified variables
209
210    REAL(r_std), DIMENSION(npts,nvm,nleafages), INTENT(inout)  :: leaf_age             !! age of the leaves (days)
211    REAL(r_std), DIMENSION(npts,nvm,nleafages), INTENT(inout)  :: leaf_frac            !! fraction of leaves in leaf age class
212                                                                                       !! (0-1, unitless)
213    REAL(r_std), DIMENSION(npts,nvm), INTENT(inout)            :: age                  !! age (years)
214    REAL(r_std), DIMENSION(npts,nvm), INTENT(inout)            :: lai                  !! leaf area index @tex ($m^2 m^{-2}$)
215                                                                                       !! @endtex
216    REAL(r_std), DIMENSION(npts,nvm,nparts,nelements), INTENT(inout) :: biomass        !! biomass @tex ($gC m^{-2}$) @endtex
217    REAL(r_std), DIMENSION(npts,nvm), INTENT(inout)            :: turnover_time        !! turnover_time of grasses (days)
218
219    !! 0.4 Local  variables
220
221    REAL(r_std), DIMENSION(npts,nvm)                           :: leaf_meanage         !! mean age of the leaves (days)
222    REAL(r_std), DIMENSION(npts)                               :: dturnover            !! Intermediate variable for turnover ??
223                                                                                       !! @tex ($gC m^{-2}$) @endtex
224    REAL(r_std), DIMENSION(npts)                               :: moiavail_crit        !! critical moisture availability, function
225                                                                                       !! of last year's moisture availability
226                                                                                       !! (0-1, unitless)
227    REAL(r_std), DIMENSION(npts)                               :: tl                   !! long term annual mean temperature, (C)
228    REAL(r_std), DIMENSION(npts)                               :: t_crit               !! critical senescence temperature, function
229                                                                                       !! of long term annual temperature (K)
230    LOGICAL, DIMENSION(npts)                                   :: shed_rest            !! shed the remaining leaves? (true/false)
231    REAL(r_std), DIMENSION(npts)                               :: sapconv              !! Sapwood conversion @tex ($gC m^{-2}$)
232                                                                                       !! @endtex
233    REAL(r_std), DIMENSION(npts)                               :: hw_old               !! old heartwood mass @tex ($gC m^{-2}$)
234                                                                                       !! @endtex
235    REAL(r_std), DIMENSION(npts)                               :: hw_new               !! new heartwood mass @tex ($gC m^{-2}$)
236                                                                                       !! @endtex
237    REAL(r_std), DIMENSION(npts)                               :: lm_old               !! old leaf mass @tex ($gC m^{-2}$) @endtex
238    REAL(r_std), DIMENSION(npts,nleafages)                     :: delta_lm             !! leaf mass change for each age class @tex
239                                                                                       !! ($gC m^{-2}$) @endtex
240    REAL(r_std), DIMENSION(npts)                               :: turnover_rate        !! turnover rate (unitless)
241    REAL(r_std), DIMENSION(npts,nvm)                           :: leaf_age_crit        !! critical leaf age (days)
242    REAL(r_std), DIMENSION(npts,nvm)                           :: new_turnover_time    !! instantaneous turnover time (days)
243    INTEGER(i_std)                                             :: j,m,k                !! Index (unitless)
244
245!_ ================================================================================================================================
246
247    IF (printlev>=3) WRITE(numout,*) 'Entering turnover'
248
249    !! 1. first call - output messages
250
251    IF ( firstcall_turnover ) THEN
252
253       IF (printlev >=2 ) THEN
254          WRITE(numout,*) 'turnover:'
255          WRITE(numout,*) ' > minimum mean leaf age for senescence (days) (::min_leaf_age_for_senescence) : '&
256               ,min_leaf_age_for_senescence
257       END IF
258       firstcall_turnover = .FALSE.
259
260
261    ENDIF
262
263    !! 2. Initializations
264
265    !! 2.1 set output to zero
266    turnover(:,:,:,:) = zero
267    new_turnover_time(:,:) = zero
268    senescence(:,:) = .FALSE.
269
270    !! 2.2 Recalculate mean leaf age
271    !      Mean leaf age is recalculated from leaf ages of separate leaf age classes. Should actually be recalculated at the
272    !      end of this routine, but it does not change too fast.
273    !      The mean leaf age is calculated using the following equation:
274    !      \latexonly
275    !      \input{turnover_lma_update_eqn1.tex}
276    !      \endlatexonly
277    !      \n
278    leaf_meanage(:,:) = zero
279
280    DO m = 1, nleafages
281       leaf_meanage(:,:) = leaf_meanage(:,:) + leaf_age(:,:,m) * leaf_frac(:,:,m)
282    ENDDO
283
284    !! 3. Climatic senescence
285
286    !     Three different types of "climatic" leaf senescence,
287    !     that do not change the age structure.
288    DO j = 2,nvm ! Loop over # PFTs
289
290       !! 3.1 Determine if there is climatic senescence.
291       !      The climatic senescence can be of three types:
292       !      sensitivity to cold temperatures, to lack of water, or both. If meteorological conditions are
293       !      fulfilled, a flag senescence is set to TRUE.
294       !      Evergreen species do not experience climatic senescence.
295
296       SELECT CASE ( senescence_type(j) )
297
298
299       CASE ('crop' )!for crop senescence is based on a GDD criterium as in crop models
300          WHERE ( ( biomass(:,j,ileaf,icarbon) .GT. zero ) .AND. &
301               ( leaf_meanage(:,j) .GT. min_leaf_age_for_senescence(j) ) .AND. &
302               ( gdd_from_growthinit(:,j) .GT.  gdd_senescence(j)))
303             senescence(:,j) = .TRUE.
304          ENDWHERE
305
306       CASE ( 'cold' )
307
308          !! 3.1.1 Summergreen species: Climatic senescence is triggered by sensitivity to cold temperatures
309          !        Climatic senescence is triggered by sensitivity to cold temperatures as follows:
310          !        If biomass is large enough (i.e. when it is greater than zero),
311          !        AND (i.e. when leaf mean age is above a certain PFT-dependent treshold ::min_leaf_age_for_senescence,
312          !        which is given in ::stomate_constants),     
313          !        AND the monthly temperature is low enough (i.e. when monthly temperature ::t2m_month(:) is below a critical
314          !        temperature ::t_crit(:), which is calculated in this module),
315          !        AND the temperature tendency is negative (i.e. when weekly temperatures ::t2m_week(:) are lower than monthly
316          !        temperatures ::t2m_month(:))
317          !        If these conditions are met, senescence is set to TRUE.
318          !
319          !        The critical temperature for senescence is calculated using the following equation:
320          !        \latexonly
321          !        \input{turnover_temp_crit_eqn2.tex}
322          !        \endlatexonly
323          !        \n
324          !
325          ! Critical temperature for senescence may depend on long term annual mean temperature
326          tl(:) = t2m_longterm(:) - ZeroCelsius
327          t_crit(:) = ZeroCelsius + senescence_temp(j,1) + &
328               tl(:) * senescence_temp(j,2) + &
329               tl(:)*tl(:) * senescence_temp(j,3)
330
331          WHERE ( ( biomass(:,j,ileaf,icarbon) .GT. zero ) .AND. &
332               ( leaf_meanage(:,j) .GT. min_leaf_age_for_senescence(j) ) .AND. &
333               ( t2m_month(:) .LT. t_crit(:) ) .AND. ( t2m_week(:) .LT. t2m_month(:) ) )
334
335             senescence(:,j) = .TRUE.
336
337          ENDWHERE
338
339       CASE ( 'dry' )
340
341          !! 3.1.2 Raingreen species: Climatic senescence is triggered by sensitivity to lack of water availability
342          !        Climatic senescence is triggered by sensitivity to lack of water availability as follows: 
343          !        If biomass is large enough (i.e. when it is greater than zero),
344          !        AND (i.e. when leaf mean age is above a certain PFT-dependent treshold ::min_leaf_age_for_senescence,
345          !        which is given in ::stomate_constants),     
346          !        AND the moisture availability drops below a critical level (i.e. when weekly moisture availability
347          !        ::moiavail_week(:,j) is below a critical moisture availability ::moiavail_crit(:),
348          !        which is calculated in this module),
349          !        If these conditions are met, senescence is set to TRUE.
350          !
351          !        The moisture availability critical level is calculated using the following equation:
352          !        \latexonly
353          !        \input{turnover_moist_crit_eqn3.tex}
354          !        \endlatexonly
355          !        \n
356          moiavail_crit(:) = &
357               MIN( MAX( minmoiavail_lastyear(:,j) + hum_frac(j) * &
358               ( maxmoiavail_lastyear(:,j) - minmoiavail_lastyear(:,j) ), &
359               senescence_hum(j) ), &
360               nosenescence_hum(j) )
361
362          WHERE ( ( biomass(:,j,ileaf,icarbon) .GT. zero ) .AND. &
363               ( leaf_meanage(:,j) .GT. min_leaf_age_for_senescence(j) ) .AND. &
364               ( moiavail_week(:,j) .LT. moiavail_crit(:) ) )
365
366             senescence(:,j) = .TRUE.
367
368          ENDWHERE
369
370       CASE ( 'mixed' )
371
372          !! 3.1.3 Mixed criterion: Climatic senescence is triggered by sensitivity to temperature or to lack of water 
373          !        Climatic senescence is triggered by sensitivity to temperature or to lack of water availability as follows:
374          !        If biomass is large enough (i.e. when it is greater than zero),
375          !        AND (i.e. when leaf mean age is above a certain PFT-dependent treshold ::min_leaf_age_for_senescence,
376          !        which is given in ::stomate_constants),     
377          !        AND the moisture availability drops below a critical level (i.e. when weekly moisture availability
378          !        ::moiavail_week(:,j) is below a critical moisture availability ::moiavail_crit(:), calculated in this module),
379          !        OR
380          !        the monthly temperature is low enough (i.e. when monthly temperature ::t2m_month(:) is below a critical
381          !        temperature ::t_crit(:), calculated in this module),
382          !        AND the temperature tendency is negative (i.e. when weekly temperatures ::t2m_week(:) are lower than
383          !        monthly temperatures ::t2m_month(:)).
384          !        If these conditions are met, senescence is set to TRUE.
385          moiavail_crit(:) = &
386               MIN( MAX( minmoiavail_lastyear(:,j) + hum_frac(j) * &
387               (maxmoiavail_lastyear(:,j) - minmoiavail_lastyear(:,j) ), &
388               senescence_hum(j) ), &
389               nosenescence_hum(j) )
390
391          tl(:) = t2m_longterm(:) - ZeroCelsius
392          t_crit(:) = ZeroCelsius + senescence_temp(j,1) + &
393               tl(:) * senescence_temp(j,2) + &
394               tl(:)*tl(:) * senescence_temp(j,3)
395
396          IF ( is_tree(j) ) THEN
397             ! critical temperature for senescence may depend on long term annual mean temperature
398             WHERE ( ( biomass(:,j,ileaf,icarbon) .GT. zero ) .AND. &
399                  ( leaf_meanage(:,j) .GT. min_leaf_age_for_senescence(j) ) .AND. &
400                  ( ( moiavail_week(:,j) .LT. moiavail_crit(:) ) .OR. &
401                  ( ( t2m_month(:) .LT. t_crit(:) ) .AND. ( t2m_week(:) .LT. t2m_month(:) ) ) ) )
402                senescence(:,j) = .TRUE.
403             ENDWHERE
404          ELSE
405
406            turnover_time(:,j) = max_turnover_time(j)
407
408            WHERE ( ( biomass(:,j,ileaf,icarbon) .GT. zero ) .AND. &
409                 ( leaf_meanage(:,j) .GT. min_leaf_age_for_senescence(j) ) .AND. &
410                 ( ( moiavail_week(:,j) .LT. moiavail_crit(:) )))
411                turnover_time(:,j) = max_turnover_time(j) * &
412                     (1.-   (1.- (moiavail_week(:,j)/  moiavail_crit(:)))**2)           
413            ENDWHERE
414
415            WHERE ( turnover_time(:,j) .LT. min_turnover_time(j) )               
416               turnover_time(:,j) = min_turnover_time(j)                         
417            ENDWHERE                                                                     
418
419            WHERE ((( biomass(:,j,ileaf,icarbon) .GT. zero ) .AND. &
420                ( leaf_meanage(:,j) .GT. min_leaf_age_for_senescence(j) ) .AND. &
421                ((t2m_month(:) .LT. t_crit(:)) .AND. (lai(:,j) .GT. lai_max(j)/4.) .OR. &
422                (t2m_month(:) .LT. ZeroCelsius)) .AND. ( t2m_week(:) .LT. t2m_month(:) )))
423               turnover_time(:,j)= leaffall(j) 
424            ENDWHERE
425           
426         ENDIF
427
428
429       !! Evergreen species do not experience climatic senescence
430       CASE ( 'none' )
431
432         
433       !! In case no climatic senescence type is recognized.
434       CASE default
435
436          WRITE(numout,*) '  turnover: don''t know how to treat this PFT.'
437          WRITE(numout,*) '  number (::j) : ',j
438          WRITE(numout,*) '  senescence type (::senescence_type(j)) : ',senescence_type(j)
439
440          CALL ipslerr_p(3,"turn","Dont know how to treat this PFT.","","")
441
442       END SELECT
443
444       !! 3.2 Drop leaves and roots, plus stems and fruits for grasses
445
446       IF ( is_tree(j) ) THEN
447
448          !! 3.2.1 Trees in climatic senescence lose their fine roots at the same rate as they lose their leaves.
449          !        The rate of biomass loss of both fine roots and leaves is presribed through the equation:
450          !        \latexonly
451          !        \input{turnover_clim_senes_biomass_eqn4.tex}
452          !        \endlatexonly
453          !        \n
454          !         with ::leaffall(j) a PFT-dependent time constant which is given in ::stomate_constants),
455          WHERE ( senescence(:,j) )
456
457             turnover(:,j,ileaf,icarbon) = biomass(:,j,ileaf,icarbon) * dt / leaffall(j)
458             turnover(:,j,iroot,icarbon) = biomass(:,j,iroot,icarbon) * dt / leaffall(j)
459
460          ENDWHERE
461
462       ELSE
463
464          !! 3.2.2 In grasses, leaf senescence is extended to the whole plant
465          !        In grasses, leaf senescence is extended to the whole plant (all carbon pools) except to its
466          !        carbohydrate reserve.     
467
468          IF (senescence_type(j) .EQ. 'crop') THEN
469             ! 3.2.2.1 crops with 'crop' phenological model
470             WHERE ( senescence(:,j) )
471                turnover(:,j,ileaf,icarbon) = biomass(:,j,ileaf,icarbon) * dt / leaffall(j)
472                turnover(:,j,iroot,icarbon) = biomass(:,j,iroot,icarbon) * dt / leaffall(j)
473                turnover(:,j,isapabove,icarbon) = biomass(:,j,isapabove,icarbon) * dt / leaffall(j)
474                turnover(:,j,ifruit,icarbon) = biomass(:,j,ifruit,icarbon) * dt /leaffall(j)
475             ENDWHERE
476          ELSE
477          ! 3.2.2.2 grass or crops based on 'mixed' phenological model
478             WHERE (turnover_time(:,j) .LT. max_turnover_time(j)) 
479                turnover(:,j,ileaf,icarbon) = biomass(:,j,ileaf,icarbon) * dt / turnover_time(:,j)
480                turnover(:,j,isapabove,icarbon) = biomass(:,j,isapabove,icarbon) * dt / turnover_time(:,j)
481                turnover(:,j,iroot,icarbon) = biomass(:,j,iroot,icarbon) * dt / turnover_time(:,j) 
482                turnover(:,j,ifruit,icarbon) = biomass(:,j,ifruit,icarbon) * dt / turnover_time(:,j)
483             ENDWHERE
484          ENDIF
485       ENDIF      ! tree/grass
486       biomass(:,j,ileaf,icarbon) = biomass(:,j,ileaf,icarbon) - turnover(:,j,ileaf,icarbon)
487       biomass(:,j,isapabove,icarbon) = biomass(:,j,isapabove,icarbon) - turnover(:,j,isapabove,icarbon)
488       biomass(:,j,iroot,icarbon) = biomass(:,j,iroot,icarbon) - turnover(:,j,iroot,icarbon)
489       biomass(:,j,ifruit,icarbon) = biomass(:,j,ifruit,icarbon) - turnover(:,j,ifruit,icarbon)
490    ENDDO        ! loop over PFTs
491
492    !! 4. Leaf fall
493    !     At a certain age, leaves fall off, even if the climate would allow a green plant
494    !     all year round. Even if the meteorological conditions are favorable for leaf maintenance,
495    !     plants, and in particular, evergreen trees, have to renew their leaves simply because the
496    !     old leaves become inefficient.   
497    !     Roots, fruits (and stems) follow leaves. The decay rate varies with leaf age.
498    !     Note that plant is not declared senescent in this case (wchich is important for allocation:
499    !     if the plant loses leaves because of their age, it can renew them).
500    !
501    !     The leaf turnover rate due to aging of leaves is calculated using the following equation:
502    !     \latexonly
503    !     \input{turnover_age_senes_biomass_eqn5.tex}
504    !     \endlatexonly
505    !     \n
506    DO j = 2,nvm ! Loop over # PFTs
507
508       !! save old leaf mass
509       lm_old(:) = biomass(:,j,ileaf,icarbon)
510
511       !! initialize leaf mass change in age class
512       delta_lm(:,:) = zero
513
514       IF ( is_tree(j) .OR. (.NOT. natural(j)) ) THEN
515
516          !! 4.1 Trees: leaves, roots, fruits roots and fruits follow leaves.
517
518          !! 4.1.1 Critical age: prescribed for trees
519          leaf_age_crit(:,j) = leafagecrit(j)
520
521       ELSE
522
523          !! 4.2 Grasses: leaves, roots, fruits, sap follow leaves.
524
525          !! 4.2.1 Critical leaf age depends on long-term temperature
526          !        Critical leaf age depends on long-term temperature
527          !        generally, lower turnover in cooler climates.
528          leaf_age_crit(:,j) = &
529               MIN( leafagecrit(j) * leaf_age_crit_coeff(1) , &
530               MAX( leafagecrit(j) * leaf_age_crit_coeff(2) , &
531               leafagecrit(j) - leaf_age_crit_coeff(3) * &
532               ( t2m_longterm(:)-ZeroCelsius - leaf_age_crit_tref ) ) )
533
534       END IF
535         
536       ! 4.2.2 Loop over leaf age classes
537       DO m = 1, nleafages
538
539          turnover_rate(:) = zero
540
541          WHERE ( leaf_age(:,j,m) .GT. leaf_age_crit(:,j)/2. )
542
543             turnover_rate(:) =  &
544                  MIN( 0.99_r_std, dt / ( leaf_age_crit(:,j) * &
545                  ( leaf_age_crit(:,j) / leaf_age(:,j,m) )**quatre ) )
546
547          ENDWHERE
548         
549          dturnover(:) = biomass(:,j,ileaf,icarbon) * leaf_frac(:,j,m) * turnover_rate(:)
550          turnover(:,j,ileaf,icarbon) = turnover(:,j,ileaf,icarbon) + dturnover(:)
551          biomass(:,j,ileaf,icarbon) = biomass(:,j,ileaf,icarbon) - dturnover(:)
552
553          ! save leaf mass change
554          delta_lm(:,m) = - dturnover(:)
555         
556          dturnover(:) = biomass(:,j,iroot,icarbon) * leaf_frac(:,j,m) * turnover_rate(:)
557          turnover(:,j,iroot,icarbon) = turnover(:,j,iroot,icarbon) + dturnover(:)
558          biomass(:,j,iroot,icarbon) = biomass(:,j,iroot,icarbon) - dturnover(:)
559         
560          dturnover(:) = biomass(:,j,ifruit,icarbon) * leaf_frac(:,j,m) * turnover_rate(:)
561          turnover(:,j,ifruit,icarbon) = turnover(:,j,ifruit,icarbon) + dturnover(:)
562          biomass(:,j,ifruit,icarbon) = biomass(:,j,ifruit,icarbon) - dturnover(:)
563         
564          IF (.NOT. is_tree(j)) THEN
565             dturnover(:) = biomass(:,j,isapabove,icarbon) * leaf_frac(:,j,m) * turnover_rate(:)
566             turnover(:,j,isapabove,icarbon) = turnover(:,j,isapabove,icarbon) + dturnover(:)
567             biomass(:,j,isapabove,icarbon) = biomass(:,j,isapabove,icarbon) - dturnover(:)
568          ENDIF
569         
570       ENDDO
571
572       !! 4.3 Recalculate the fraction of leaf biomass in each leaf age class.
573       !      Older leaves will fall more fast than younger leaves and therefore
574       !      the leaf age distribution needs to be recalculated after turnover.
575       !      The fraction of biomass in each leaf class is updated using the following equation:
576       !      \latexonly
577       !      \input{turnover_update_LeafAgeDistribution_eqn6.tex}
578       !      \endlatexonly
579       !      \n
580       !
581       !      new fraction = new leaf mass of that fraction / new total leaf mass
582       !                   = (old fraction*old total leaf mass ::lm_old(:) + biomass change of that fraction ::delta_lm(:,m)  ) /
583       !                     new total leaf mass ::biomass(:,j,ileaf
584       DO m = 1, nleafages
585         
586          WHERE ( biomass(:,j,ileaf,icarbon) .GT. min_sechiba )
587             leaf_frac(:,j,m) = ( leaf_frac(:,j,m)*lm_old(:) + delta_lm(:,m) ) / biomass(:,j,ileaf,icarbon)
588          ELSEWHERE
589             leaf_frac(:,j,m) = zero
590          ENDWHERE
591       
592       ENDDO
593       
594    ENDDO         ! loop over PFTs
595
596    !! 5. New (provisional) LAI
597    !     ::lai(:,j) is determined from the leaf biomass ::biomass(:,j,ileaf,icarbon) and the
598    !     specific leaf surface :: sla(j) (m^2 gC^{-1})
599    !     The leaf area index is updated using the following equation:
600    !     \latexonly
601    !     \input{turnover_update_LAI_eqn7.tex}
602    !     \endlatexonly
603    !     \n
604
605    !    lai(:,ibare_sechiba) = zero
606    !    DO j = 2, nvm ! Loop over # PFTs
607    !        lai(:,j) = biomass(:,j,ileaf,icarbon) * sla(j)
608    !    ENDDO
609
610    !! 6. Definitely drop all leaves if there is a very low leaf mass during senescence.
611
612    !     Both for deciduous trees and grasses same conditions are checked:
613    !     If biomass is large enough (i.e. when it is greater than zero),
614    !     AND when senescence is set to true
615    !     AND the leaf biomass drops below a critical minimum biomass level (i.e. when it is lower than half
616    !     the minimum initial LAI ::lai_initmin(j) divided by the specific leaf area ::sla(j),
617    !     ::lai_initmin(j) is set to 0.3 in stomate_data.f90 and sla is a constant that is set to 0.015366 m2/gC),
618    !     If these conditions are met, the flag ::shed_rest(:) is set to TRUE.
619    !
620    !     After this, the biomass of different carbon pools both for trees and grasses is set to zero
621    !     and the mean leaf age is reset to zero.
622    !     Finally, the leaf fraction and leaf age of the different leaf age classes is set to zero.
623    DO j = 2,nvm ! Loop over # PFTs
624
625       shed_rest(:) = .FALSE.
626
627       !! 6.1 For deciduous trees: next to leaves, also fruits and fine roots are dropped
628       !      For deciduous trees: next to leaves, also fruits and fine roots are dropped: fruit ::biomass(:,j,ifruit)
629       !      and fine root ::biomass(:,j,iroot) carbon pools are set to zero.
630       IF ( is_tree(j) .AND. ( senescence_type(j) .NE. 'none' ) ) THEN
631
632          ! check whether we shed the remaining leaves
633          WHERE ( ( biomass(:,j,ileaf,icarbon) .GT. zero ) .AND. senescence(:,j) .AND. &
634               ( biomass(:,j,ileaf,icarbon) .LT. (lai_initmin(j) / 2.)/sla(j) )             )
635
636             shed_rest(:) = .TRUE.
637
638             turnover(:,j,ileaf,icarbon)  = turnover(:,j,ileaf,icarbon) + biomass(:,j,ileaf,icarbon)
639             turnover(:,j,iroot,icarbon)  = turnover(:,j,iroot,icarbon) + biomass(:,j,iroot,icarbon)
640             turnover(:,j,ifruit,icarbon) = turnover(:,j,ifruit,icarbon) + biomass(:,j,ifruit,icarbon)
641
642             biomass(:,j,ileaf,icarbon)  = zero
643             biomass(:,j,iroot,icarbon)  = zero
644             biomass(:,j,ifruit,icarbon) = zero
645
646             ! reset leaf age and lai
647             leaf_meanage(:,j) = zero
648             lai(:,j) = zero
649          ENDWHERE
650
651       ENDIF
652
653       !! 6.2 For grasses: all aboveground carbon pools, except the carbohydrate reserves are affected:
654       !      For grasses: all aboveground carbon pools, except the carbohydrate reserves are affected:
655       !      fruit ::biomass(:,j,ifruit,icarbon), fine root ::biomass(:,j,iroot,icarbon) and sapwood above
656       !      ::biomass(:,j,isapabove,icarbon) carbon pools are set to zero.
657       IF ( .NOT. is_tree(j) ) THEN
658
659          ! Shed the remaining leaves if LAI very low.
660          WHERE ( ( biomass(:,j,ileaf,icarbon) .GT. zero ) .AND. senescence(:,j) .AND. &
661               (  biomass(:,j,ileaf,icarbon) .LT. (lai_initmin(j) / 2.)/sla(j) ))
662
663             shed_rest(:) = .TRUE.
664
665             turnover(:,j,ileaf,icarbon) = turnover(:,j,ileaf,icarbon) + biomass(:,j,ileaf,icarbon)
666             turnover(:,j,isapabove,icarbon) = turnover(:,j,isapabove,icarbon) + biomass(:,j,isapabove,icarbon)
667             turnover(:,j,iroot,icarbon) = turnover(:,j,iroot,icarbon) + biomass(:,j,iroot,icarbon)
668             turnover(:,j,ifruit,icarbon) = turnover(:,j,ifruit,icarbon) + biomass(:,j,ifruit,icarbon)
669
670             biomass(:,j,ileaf,icarbon) = zero
671             biomass(:,j,isapabove,icarbon) = zero
672             biomass(:,j,iroot,icarbon) = zero
673             biomass(:,j,ifruit,icarbon) = zero
674
675             ! reset leaf age and lai
676             leaf_meanage(:,j) = zero
677             lai(:,j) = zero
678          ENDWHERE
679
680       ENDIF
681
682       !! 6.3 Reset the leaf age structure: the leaf fraction and leaf age of the different leaf age classes is set to zero.
683     
684       DO m = 1, nleafages
685
686          WHERE ( shed_rest(:) )
687
688             leaf_age(:,j,m) = zero
689             leaf_frac(:,j,m) = zero
690
691          ENDWHERE
692
693       ENDDO
694
695    ENDDO          ! loop over PFTs
696   
697    !! 7. Herbivore activity: elephants, cows, gazelles but no lions.
698 
699    !     Herbivore activity affects the biomass of leaves and fruits as well
700    !     as stalks (only for grasses). Herbivore activity does not modify leaf
701    !     age structure. Herbivores ::herbivores(:,j) is the time constant of
702    !     probability of a leaf to be eaten by a herbivore, and is calculated in
703    !     ::stomate_season. following Mc Naughton et al. [1989].
704
705    IF ( ok_herbivores ) THEN
706
707       ! If the herbivore activity is allowed (if ::ok_herbivores is true, which is set in run.def),
708       ! remove the amount of biomass consumed by herbivory from the leaf biomass ::biomass(:,j,ileaf,icarbon) and
709       ! the fruit biomass ::biomass(:,j,ifruit,icarbon).
710       ! The daily amount consumed equals the biomass multiplied by 1 day divided by the time constant ::herbivores(:,j).
711       DO j = 2,nvm ! Loop over # PFTs
712
713          IF ( is_tree(j) ) THEN
714
715             !! For trees: only the leaves and fruit carbon pools are affected
716
717             WHERE (biomass(:,j,ileaf,icarbon) .GT. zero)
718                ! added by shilong
719                WHERE (herbivores(:,j).GT. min_sechiba)
720                   dturnover(:) = biomass(:,j,ileaf,icarbon) * dt / herbivores(:,j)
721                   turnover(:,j,ileaf,icarbon) = turnover(:,j,ileaf,icarbon) + dturnover(:)
722                   biomass(:,j,ileaf,icarbon) = biomass(:,j,ileaf,icarbon) - dturnover(:)
723
724                   dturnover(:) = biomass(:,j,ifruit,icarbon) * dt / herbivores(:,j)
725                   turnover(:,j,ifruit,icarbon) = turnover(:,j,ifruit,icarbon) + dturnover(:)
726                   biomass(:,j,ifruit,icarbon) = biomass(:,j,ifruit,icarbon) - dturnover(:)
727                ENDWHERE
728             ENDWHERE
729
730          ELSE
731
732             ! For grasses: all aboveground carbon pools are affected: leaves, fruits and sapwood above
733             WHERE ( biomass(:,j,ileaf,icarbon) .GT. zero )
734                ! added by shilong
735                WHERE (herbivores(:,j) .GT. min_sechiba)
736                   dturnover(:) = biomass(:,j,ileaf,icarbon) * dt / herbivores(:,j)
737                   turnover(:,j,ileaf,icarbon) = turnover(:,j,ileaf,icarbon) + dturnover(:)
738                   biomass(:,j,ileaf,icarbon) = biomass(:,j,ileaf,icarbon) - dturnover(:)
739
740                   dturnover(:) = biomass(:,j,isapabove,icarbon) * dt / herbivores(:,j)
741                   turnover(:,j,isapabove,icarbon) = turnover(:,j,isapabove,icarbon) + dturnover(:)
742                   biomass(:,j,isapabove,icarbon) = biomass(:,j,isapabove,icarbon) - dturnover(:)
743
744                   dturnover(:) = biomass(:,j,ifruit,icarbon) * dt / herbivores(:,j)
745                   turnover(:,j,ifruit,icarbon) = turnover(:,j,ifruit,icarbon) + dturnover(:)
746                   biomass(:,j,ifruit,icarbon) = biomass(:,j,ifruit,icarbon) - dturnover(:)
747                ENDWHERE
748
749             ENDWHERE
750
751          ENDIF  ! tree/grass?
752
753       ENDDO    ! loop over PFTs
754
755    ENDIF ! end herbivores
756
757    !! 8. Fruit turnover for trees
758
759    !     Fruit turnover for trees: trees simply lose their fruits with a time constant ::tau_fruit(j),
760    !     that is set to 90 days for all PFTs in ::stomate_constants
761
762    DO k = 1,nelements 
763       DO j = 2,nvm ! Loop over # PFTs
764          IF ( is_tree(j) ) THEN
765
766             dturnover(:) = biomass(:,j,ifruit,k) * dt / tau_fruit(j)
767             turnover(:,j,ifruit,k) = turnover(:,j,ifruit,k) + dturnover(:)
768             biomass(:,j,ifruit,k) = biomass(:,j,ifruit,k) - dturnover(:)
769             
770          ENDIF
771       ENDDO       ! loop over PFTs
772    END DO
773
774    !! 9 Conversion of sapwood to heartwood both for aboveground and belowground sapwood and heartwood.
775
776    !   Following LPJ (Sitch et al., 2003), sapwood biomass is converted into heartwood biomass
777    !   with a time constant tau ::tau_sap(j) of 1 year.
778    !   Note that this biomass conversion is not added to "turnover" as the biomass is not lost!
779    DO j = 2,nvm ! Loop over # PFTs
780
781       IF ( is_tree(j) ) THEN
782
783          !! For the recalculation of age in 9.2 (in case the vegetation is not dynamic ie. ::ok_dgvm is FALSE),
784          !! the heartwood above and below is stored in ::hw_old(:).
785          IF ( .NOT. ok_dgvm ) THEN
786             hw_old(:) = biomass(:,j,iheartabove,icarbon) + biomass(:,j,iheartbelow,icarbon)
787          ENDIF
788
789          !! 9.1 Calculate the rate of sapwood to heartwood conversion
790          !      Calculate the rate of sapwood to heartwood conversion with the time constant ::tau_sap(j)
791          !      and update aboveground and belowground sapwood ::biomass(:,j,isapabove) and ::biomass(:,j,isapbelow)
792          !      and heartwood ::biomass(:,j,iheartabove) and ::biomass(:,j,iheartbelow).
793
794          DO k = 1,nelements
795
796             ! Above the ground
797             sapconv(:) = biomass(:,j,isapabove,k) * dt / tau_sap(j)
798             biomass(:,j,isapabove,k) = biomass(:,j,isapabove,k) - sapconv(:)
799             biomass(:,j,iheartabove,k) =  biomass(:,j,iheartabove,k) + sapconv(:)
800             
801             ! Below the ground
802             sapconv(:) = biomass(:,j,isapbelow,k) * dt / tau_sap(j)
803             biomass(:,j,isapbelow,k) = biomass(:,j,isapbelow,k) - sapconv(:)
804             biomass(:,j,iheartbelow,k) =  biomass(:,j,iheartbelow,k) + sapconv(:)
805
806          END DO
807
808          !! 9.2 If the vegetation is not dynamic, the age of the plant is decreased.
809          !      The updated heartwood, the sum of new heartwood above and new heartwood below after
810          !      converting sapwood to heartwood, is saved as ::hw_new(:) .
811          !      Creation of new heartwood decreases the age of the plant with a factor that is determined by:
812          !      old heartwood ::hw_old(:) divided by the new heartwood ::hw_new(:)
813          IF ( .NOT. ok_dgvm ) THEN
814
815             hw_new(:) = biomass(:,j,iheartabove,icarbon) + biomass(:,j,iheartbelow,icarbon)
816
817             WHERE ( hw_new(:) .GT. min_sechiba )
818
819                age(:,j) = age(:,j) * hw_old(:)/hw_new(:)
820
821             ENDWHERE
822
823          ENDIF
824
825       ENDIF
826
827    ENDDO       ! loop over PFTs
828
829
830    CALL xios_orchidee_send_field("HERBIVORES",herbivores)
831    CALL xios_orchidee_send_field("LEAF_AGE",leaf_meanage)
832   
833
834    ! Write mean leaf age and time constant of probability of a leaf to be eaten by a herbivore
835    ! to the stomate output file.
836    CALL histwrite_p (hist_id_stomate, 'LEAF_AGE', itime, &
837         leaf_meanage, npts*nvm, horipft_index)
838    CALL histwrite_p (hist_id_stomate, 'HERBIVORES', itime, &
839         herbivores, npts*nvm, horipft_index)
840
841    IF (printlev>=4) WRITE(numout,*) 'Leaving turnover'
842
843  END SUBROUTINE turn
844
845END MODULE stomate_turnover
Note: See TracBrowser for help on using the repository browser.