Documentation/ORCHIDEE_DOC: stomate_alloc.f90

File stomate_alloc.f90, 31.5 KB (added by maignan, 9 years ago)
Line 
1!! IPSL (2006)
2!! This software is governed by the CeCILL licence see ORCHIDEE/ORCHIDEE_CeCILL.LIC
3
4MODULE stomate_alloc
5
6  ! Modules used:
7
8  USE ioipsl
9  USE stomate_constants
10  USE constantes_veg
11
12  IMPLICIT NONE
13
14  ! Private & public routines
15
16  PRIVATE
17  PUBLIC alloc,alloc_clear
18
19 ! Variables shared by all subroutines in this module
20
21 ! Is this the first call? (true/false)
22  LOGICAL, SAVE                                             :: firstcall = .TRUE.
23CONTAINS
24
25
26!! ==============================================================================\n
27!! SUBROUTINE   : alloc_clear
28!! AUTHOR       : P. Friedlingstein
29!! CREATION DATE:
30!!
31!> BRIEF        : Set the flag ::firstcall to .TRUE. and as such activate section
32!! 1.1 of the subroutine alloc (see below).\n
33!! ==============================================================================
34  SUBROUTINE alloc_clear
35    firstcall = .TRUE.
36  END SUBROUTINE alloc_clear
37
38
39
40!! ==============================================================================\n
41!! SUBROUTINE   : alloc
42!! AUTHOR       : P. Friedlingstein
43!! CREATION DATE:
44!!
45!> BRIEF        : Allocate net primairy production (= photosynthesis
46!! minus autothrophic respiration) to: carbon reserves, aboveground sapwood,
47!! belowground sapwood, root, fruits and leaves following Friedlingstein et al. (1999).\n
48!!
49!! DESCRIPTION (definitions, functional, design, flags):\n
50!! The philosophy underluying the scheme is that allocation patterns result from
51!! evolved responses that adjust carbon investments to facilitate capture of most
52!! limiting resources i.e. light, water and mineral nitrogen. The implemented scheme
53!! calculates the limitation of light, water and nitrogen. However, nitrogen is not a
54!! prognostic variable of the model and therefore soil temperature and soil moisture
55!! are used as a proxy for soil nitrogen availability.\n
56!! Sharpe & Rykiel (1991) proposed a generic relationship between the allocation of
57!! carbon to a given plant compartment and the availability of a particular resource:\n
58!! \latexonly
59!! \input{alloc1.tex}
60!! \endlatexonly
61!! \n
62!! where A is the allocation of biomass production (NPP) to a given compartment (either
63!! leaves, stem, or roots). Xi and Yj are resource availabilities (e.g. light, water,
64!! nutrient). For a given plant compartment, a resource can be of type X or Y. An increase
65!! in a X-type resource will increase the allocation to compartment A. An increase in a
66!! Y-type resource will, however, lead to a decrease in carbon allocation to that compartment.
67!! In other words, Y-type resources are those for which uptake increases with increased
68!! investment in the compartment in question. X-type resources, as a consequence of
69!! trade-offs, are the opposite. For example, water is a Y-type resource for root allocation.
70!! Water-limited conditions should promote carbon allocation to roots, which enhance water
71!! uptake and hence minimize plant water stress. Negative relationships between investment
72!! and uptake arise when increased investment in one compartment leads, as required for
73!! conservation of mass, to decreased investment in a component involved in uptake of
74!! that resource.\n
75!!
76!! The implemented scheme allocates carbon to the following components:\n
77!! - Carbon reserves;\n
78!! - Aboveground sapwood;\n
79!! - Belowground sapwood;\n
80!! - Roots;\n
81!! - Fruits/seeds and\n
82!! - Leaves.
83!! \n
84!!
85!! The allocation to fruits and seeds is simply a 10% "tax" of the total biomass
86!! production.\n
87!  This 'tax' could be made PFT-dependent as higher values are
88!  likely for crops. Further it should depend on the limitations that the plant
89!  experiences. If the PFT fares well, it should grow more fruits. However, this
90!  could only work if PFTs are also rewarded for having grown fruits i.e. by making
91!  the reproduction rate depend on the fruit growth of the past years. Otherwise,
92!  the fruit allocation would be a punishment for plants that are doing well.
93!! Following carbohydrate use to support budburst and initial growth, the
94!! carbohydrate reserve is refilled. The daily amount of carbon allocated to the
95!! reserve pool is proportional to leaf+root allocation (::LtoLSR and ::RtoLSR).\n
96!! Sapwood and root allocation (respectively ::StoLSR and ::RtoLSR) are proportional
97!! to the estimated light and soil (water and nitrogen) stress (::Limit_L and
98!! ::Limit_NtoW). Further, Sapwood allocation is separated in belowground sapwood
99!! and aboveground sapwood making use of the parameter (:: alloc_sap_above_tree
100!! or ::alloc_sap_above_grass). For trees partitioning between above and
101!! belowground compartments is a function of PFT age.\n
102!! Leaf allocation (::LtoLSR) is calculated as the residual of root and sapwood
103!! allocation (LtoLSR(:) = 1. - RtoLSR(:) - StoLSR(:).\n
104!!
105!! MAIN OUTPUT VARIABLE(S): :: f_alloc; fraction of NPP that is allocated to the
106!! six different biomass compartments (leaves, roots, above and belowground wood,
107!! carbohydrate reserves and fruits). DIMENSION(npts,nvm,nparts).\n
108!!
109!! REFERENCES   :
110!! - Friedlingstein, P., G. Joel, C.B. Field, and Y. Fung (1999), Towards an allocation
111!! scheme for global terrestrial carbon models, Global Change Biology, 5, 755-770.
112!! - Sharpe, P.J.H., and Rykiel, E.J. (1991), Modelling integrated response of plants
113!! to multiple stresses. In: Response of Plants to Multiple Stresses (eds Mooney, H.A.,
114!! Winner, W.E., Pell, E.J.), pp. 205-224, Academic Press, San Diego, CA.
115!!
116!! FLOWCHART    :
117!! \latexonly
118!! \includegraphics(scale=0.5){allocflow.jpg}
119!! \endlatexonly
120!! ==============================================================================
121  SUBROUTINE alloc (npts, dt, &
122       lai, veget_max, senescence, when_growthinit, &
123       moiavail_week, tsoil_month, soilhum_month, &
124       biomass, age, leaf_age, leaf_frac, rprof, f_alloc)
125
126!! Input variables
127
128!! Domain size - number of pixels (dimensionless)
129    INTEGER(i_std), INTENT(in)                                 :: npts
130!! Time step of the simulations for stomate (days)
131    REAL(r_std), INTENT(in)                                    :: dt
132!! PFT leaf area index (m**2/m**2)
133    REAL(r_std), DIMENSION(npts,nvm), INTENT(in)               :: lai
134!! PFT "Maximal" coverage fraction of a PFT (= ind*cn_ind) (m**2/m**2)
135    REAL(r_std), DIMENSION(npts,nvm), INTENT(in)               :: veget_max
136!! Is the PFT senescent? (true/false) - only for deciduous trees
137    LOGICAL, DIMENSION(npts,nvm), INTENT(in)                   :: senescence
138!! Days since beginning of growing season (days)
139    REAL(r_std), DIMENSION(npts,nvm), INTENT(in)               :: when_growthinit
140!! PFT moisture availability (??units??) - integrated over a week
141    REAL(r_std), DIMENSION(npts,nvm), INTENT(in)               :: moiavail_week
142!! PFT soil temperature (K) - integrated over a month
143    REAL(r_std), DIMENSION(npts,nbdl), INTENT(in)              :: tsoil_month
144!! PFT soil humidity (??units??) - integrated over a month
145    REAL(r_std), DIMENSION(npts,nbdl), INTENT(in)              :: soilhum_month
146!! PFT age (days)
147    REAL(r_std), DIMENSION(npts,nvm), INTENT(in)               :: age
148!! PFT total biomass (gC/m**2)
149    REAL(r_std), DIMENSION(npts,nvm,nparts), INTENT(inout)     :: biomass
150!! PFT age of different leaf classes (days)
151    REAL(r_std), DIMENSION(npts,nvm,nleafages), INTENT(inout)  :: leaf_age
152!! PFT fraction of leaves in leaf age class (dimensionless; 1)
153    REAL(r_std), DIMENSION(npts,nvm,nleafages), INTENT(inout)  :: leaf_frac
154
155!! Output variables
156
157!! [DISPENSABLE] PFT rooting depth - not calculated in the current version of the model (m)
158    REAL(r_std), DIMENSION(npts,nvm), INTENT(inout)            :: rprof
159!! PFT fraction of NPP that is allocated to the different components (dimensionless; 1)
160    REAL(r_std), DIMENSION(npts,nvm,nparts), INTENT(out)       :: f_alloc
161
162!! Local variables
163
164!! Lai threshold below which carbohydrate reserve may be used (m**2/m**2)
165    REAL(r_std), DIMENSION(nvm)                                :: lai_happy
166!! Lights stress (dimensionless; 1)
167    REAL(r_std), DIMENSION(npts)                               :: limit_L
168!! Total nitrogen stress (dimensionless; 1)
169    REAL(r_std), DIMENSION(npts)                               :: limit_N
170!! Stress from soil temperature on nitrogen mineralisation (dimensionless; 1)
171    REAL(r_std), DIMENSION(npts)                               :: limit_N_temp
172!! Stress from soil humidity on nitrogen mineralisation (dimensionless; 1)
173    REAL(r_std), DIMENSION(npts)                               :: limit_N_hum
174!! Soil water stress (dimensionless; 1)
175    REAL(r_std), DIMENSION(npts)                               :: limit_W
176!! Most limiting factor in the soil: nitrogen or water (dimensionless; 1)
177    REAL(r_std), DIMENSION(npts)                               :: limit_WorN
178!! Most limiting factor: amongst limit_N, limit_W and limit_L
179    REAL(r_std), DIMENSION(npts)                               :: limit
180!! Preliminairy soil temperature stress used as a proxy for niterogen stress (K)
181    REAL(r_std), DIMENSION(npts)                               :: t_nitrogen
182!! Preliminairy soil humidity stress used as a proxy for nitrogen stress (??units??)
183    REAL(r_std), DIMENSION(npts)                               :: h_nitrogen
184!! Scaling factor for integrating vertical soil profiles (-)
185    REAL(r_std), DIMENSION(npts)                               :: rpc
186!! Ratio between leaf-allocation and (leaf+sapwood+root)-allocation (dimensionless; 1)
187    REAL(r_std), DIMENSION(npts)                               :: LtoLSR
188!! Ratio between sapwood-allocation and (leaf+sapwood+root)-allocation (dimensionless; 1)
189    REAL(r_std), DIMENSION(npts)                               :: StoLSR
190!! Ratio between root-allocation and (leaf+sapwood+root)-allocation (dimensionless; 1)
191    REAL(r_std), DIMENSION(npts)                               :: RtoLSR
192!! Rescaling factor for allocation factors is carbon is allocated to carbohydrate reserve (dimensionless; 1)
193    REAL(r_std), DIMENSION(npts)                               :: carb_rescale
194!! Mass of carbohydrate reserve used to support growth (gC/m**2)
195    REAL(r_std), DIMENSION(npts)                               :: use_reserve
196!! Fraction of carbohydrate reserve used (::use_reserve) to support leaf growth (gC/m**2)
197    REAL(r_std), DIMENSION(npts)                               :: transloc_leaf
198!! Leaf biomass in youngest leaf age class (gC/m**2)
199    REAL(r_std), DIMENSION(npts)                               :: leaf_mass_young
200!! Variable to store leaf biomass from previous time step (gC/m**2)
201    REAL(r_std), DIMENSION(npts,nvm)                           :: lm_old
202!! Maximum number of days during which carbohydrate reserve may be used (days)
203    REAL(r_std)                                                :: reserve_time
204!! lai on natural part of the grid cell, or of agricultural PFTs
205    REAL(r_std), DIMENSION(npts,nvm)                           :: lai_around
206!! Vegetation cover of natural PFTs on the grid cell (agriculture masked)
207    REAL(r_std), DIMENSION(npts,nvm)                           :: veget_max_nat 
208!! Total natural vegetation cover on natural part of the grid cell
209    REAL(r_std), DIMENSION(npts)                               :: natveg_tot
210!! Average LAI on natural part of the grid cell
211    REAL(r_std), DIMENSION(npts)                               :: lai_nat
212!! [DISPENSABLE] intermediate array for looking for minimum
213    REAL(r_std), DIMENSION(npts)                               :: zdiff_min
214!! Prescribed fraction of sapwood allocation to above ground sapwood (dimensionless; 1)
215    REAL(r_std), DIMENSION(npts)                               :: alloc_sap_above
216!! Variable to store depth of the different soil layers (m)
217    REAL(r_std), SAVE, DIMENSION(0:nbdl)                       :: z_soil
218!! Indeces
219    INTEGER(i_std)                                             :: i,j,l,m
220
221!! Define local parameters
222
223!! [DISPENSABLE] Do we try to reach a minimum reservoir even if we are severely stressed? (true/false)
224    ! Duplicates bavard - no other function in this routine
225    LOGICAL, PARAMETER                                         :: ok_minres = .TRUE.
226!! Time required to develop a minimal LAI using the carbohydrate reserve (days)
227    REAL(r_std), PARAMETER                                     :: tau_leafinit = 10.
228!! Maximum number of days during which carbohydrate reserve may be used for trees (days)
229    REAL(r_std), PARAMETER                                     :: reserve_time_tree = 30.
230!! Maximum number of days during which carbohydrate reserve may be used for grasses (days)
231    REAL(r_std), PARAMETER                                     :: reserve_time_grass = 20.
232!! Default root allocation (dimensionless; 1)
233    REAL(r_std), PARAMETER                                     :: R0 = 0.3
234!! Default sapwood allocation (dimensionless; 1)
235    REAL(r_std), PARAMETER                                     :: S0 = 0.3
236!! Default leaf allocation (dimensionless; 1)
237    REAL(r_std), PARAMETER                                     :: L0 = 1. - R0 - S0
238!! Default fruit allocation (dimensionless; 1)
239    REAL(r_std), PARAMETER                                     :: f_fruit = 0.1
240!! Fraction of sapwood allocation above ground for trees (dimensionless; 1)
241    REAL(r_std), PARAMETER                                     :: alloc_sap_above_tree = 0.5
242!! Fraction of sapwood allocation above ground for trees (dimensionless; 1)
243    REAL(r_std), PARAMETER                                     :: alloc_sap_above_grass = 1.0
244!! Prescribed lower and upper bounds for leaf allocation (dimensionless; 1)
245    REAL(r_std), PARAMETER                                     :: min_LtoLSR = 0.2
246    REAL(r_std), PARAMETER                                     :: max_LtoLSR = 0.5
247!! Curvature of the root profile (m)
248    REAL(r_std), PARAMETER                                     :: z_nitrogen = 0.2
249!! =========================================================================
250
251
252!> @addtogroup Allocation Resource limitation and C-allocation
253!> @{
254!> @}
255    IF (bavard.GE.3) WRITE(numout,*) 'Entering alloc'
256!> @addtogroup Allocation
257!> 1. Initialize\n
258    !> 1.1 First call only\n
259    IF ( firstcall ) THEN
260       
261       !> @addtogroup Allocation
262       !> 1.1.1 Copy the depth of the different soil layers from ??constantes??\n
263       z_soil(0) = 0.
264       z_soil(1:nbdl) = diaglev(1:nbdl)
265
266!> @addtogroup Allocation
267       ! 1.1.2 Print flags and parameter settings\n
268       !!?? Should we mention the variable neame between brackets? Could help to link message
269       !!?? to variables??
270       WRITE(numout,*) 'alloc:'
271       WRITE(numout,'(a,$)') '    > We'
272       IF ( .NOT. ok_minres ) WRITE(numout,'(a,$)') ' do NOT'
273       WRITE(numout,*) 'try to reach a minumum reservoir when severely stressed.'
274       WRITE(numout,*) '   > Time delay (days) to build leaf mass (::tau_leafinit): ', &
275            tau_leafinit
276       WRITE(numout,*) '   > Curvature of root mass with increasing soil depth (::z_nitrogen): ', &
277            z_nitrogen
278       WRITE(numout,*) '   > Sap allocation above the ground / total sap allocation (dimensionless; 1): '
279       WRITE(numout,*) '       trees (::alloc_sap_above_tree):', alloc_sap_above_tree
280       WRITE(numout,*) '       grasses (::alloc_sap_above_grass) :', alloc_sap_above_grass
281       WRITE(numout,*) '   > Default root alloc fraction (1; ::R0): ', R0
282       WRITE(numout,*) '   > Default sapwood alloc fraction (1; ::S0): ', S0
283       WRITE(numout,*) '   > Default fruit allocation (1, ::f_fruit): ', f_fruit
284       WRITE(numout,*) '   > Minimum (min_LtoLSR)/maximum (::max_LtoLSR)leaf alloc fraction (dimensionless; 1): ',&
285            min_LtoLSR,max_LtoLSR
286       WRITE(numout,*) '   > Maximum time (days) the carbon reserve can be used:'
287       WRITE(numout,*) '       trees (reserve_time_tree):',reserve_time_tree
288       WRITE(numout,*) '       grasses (reserve_time_grass):',reserve_time_grass
289       firstcall = .FALSE.
290
291    ENDIF
292
293!> @addtogroup Allocation
294    !> 1.2 Every call\n
295    !> 1.2.1 Reset output variable (::f_alloc)\n
296    f_alloc(:,:,:) = 0.0
297    f_alloc(:,:,icarbres) = 1.0
298
299 !> @addtogroup Allocation
300    !> 1.2.2 Proxy for soil nitrogen stress\n
301    !> Nitrogen availability and thus N-stress can not be calculated by the model. Water and
302    !> temperature stress are used as proxy under the assumption that microbial activity is
303    !> determined by soil temperature and water availability. In turn, microbial activity is
304    !> assumed to be an indicator for nitrogen mineralisation and thus its availability.\n
305
306    !> 1.2.2.1 Convolution of nitrogen stress with root profile\n
307    !> The capacity of roots to take up nitrogen is assumed to decrease exponentially with
308    !> increasing soil depth. The curvature of the exponential function describing the
309    !> nitrogen-uptake capacity of roots (= root mass * uptake capacity) is given by
310    !> ::z_nitrogen. Strictly speaking its unit is meters (m). Despite its units this parameter
311    !> has no physical meaning.\n
312    !> Because the roots are described by an exponential function but the soil depth is limited to
313    !> ::z_soil(nbdl), the root profile is truncated at ::z_soil(nbdl). For numerical reasons,
314    !> the total capacity of the soil profile for nitrogen uptake should be 1. To this aim a scaling
315    !> factor (::rpc) is calculated as follows:\n
316    !! \latexonly
317    !! \input{alloc2.tex}
318    !! \endlatexonly
319    !> \n
320    !> Following, mean weighted (weighted by nitrogen uptake capacity) soil temperature
321    !> (::tsoil_month) or soil moisture (::t_soil_month) are calculated. In this calculation
322    !> ::rpc is used for scaling.\n
323       
324    ! Scaling factor for integration
325    rpc(:) = 1. / ( 1. - EXP( -z_soil(nbdl) / z_nitrogen ) )
326
327    ! Integrate over # soil layers
328    t_nitrogen(:) = 0.
329
330    DO l = 1, nbdl ! Loop over # soil layers
331
332       t_nitrogen(:) = &
333            t_nitrogen(:) + tsoil_month(:,l) * rpc(:) * &
334            ( EXP( -z_soil(l-1)/z_nitrogen ) - EXP( -z_soil(l)/z_nitrogen ) )
335
336    ENDDO ! Loop over # soil layers
337
338 !> @addtogroup Allocation
339   !> 1.2.2.2 Convolution for soil moisture
340
341    ! Scaling factor for integration
342    rpc(:) = 1. / ( 1. - EXP( -z_soil(nbdl) / z_nitrogen ) )
343
344    ! Integrate over # soil layers
345
346    h_nitrogen(:) = 0.0
347
348    DO l = 1, nbdl ! Loop over # soil layers
349
350       h_nitrogen(:) = &
351            h_nitrogen(:) + soilhum_month(:,l) * rpc(:) * &
352            ( EXP( -z_soil(l-1)/z_nitrogen ) - EXP( -z_soil(l)/z_nitrogen ) )
353
354    ENDDO ! Loop over # soil layers
355
356!> @addtogroup Allocation
357    !> 1.2.3 Separate between natural and agrigultural LAI\n
358    !> The models distinguishes different natural PFTs but does not contain information
359    !> on whether these PFTs are spatially separated or mixed. In line with the DGVM the
360    !> models treats the natural PFT's as mixed. Therefore, the average LAI over the
361    !> natural PFTs is calculated to estimate light stress. Agricultural PFTs are spatially
362    !> separated.
363    natveg_tot(:) = 0.0
364    lai_nat(:) = 0.0
365
366    DO j = 2, nvm ! Loop over # PFTs
367
368       IF ( natural(j) ) THEN
369          ! Mask agricultural vegetation
370          veget_max_nat(:,j) = veget_max(:,j)
371       ELSE
372          ! Mask natural vegetation
373          veget_max_nat(:,j) = 0.0
374       ENDIF
375
376       ! Sum up fraction of natural space covered by vegetation
377       natveg_tot(:) = natveg_tot(:) + veget_max_nat(:,j)
378
379       ! Sum up lai
380       lai_nat(:) = lai_nat(:) + veget_max_nat(:,j) * lai(:,j)
381
382    ENDDO ! Loop over # PFTs
383
384    DO j = 2, nvm ! Loop over # PFTs
385
386       IF ( natural(j) ) THEN
387          ! Use the mean LAI over all natural PFTs when estimating light stress
388          ! on a specific natural PFT
389          lai_around(:,j) = lai_nat(:)
390       ELSE
391          ! Use the actual LAI (specific for that PFT) when estimating light
392          ! stress on a specific agricultural PFT
393          lai_around(:,j) = lai(:,j)
394       ENDIF
395
396    ENDDO ! Loop over # PFTs
397
398!> @addtogroup Allocation
399    !> 1.2.4 Calculate LAI threshold (m**2/m**2) below which carbohydrate reserve is used.
400   
401    ! Lai_max is a PFT-dependent parameter specified in ??WHERE??
402    lai_happy(:) = lai_max(:) * 0.5
403
404!> @addtogroup Allocation
405!> 2. Use carbohydrate reserve to support growth and update leaf age\n
406
407    ! Save old leaf mass, biomass got last updated in ??WHERE??
408    lm_old(:,:) = biomass(:,:,ileaf)
409
410    DO j = 2, nvm ! Loop over # PFTs
411
412!> @addtogroup Allocation
413       !> 2.1 Calculate demand for carbonhydrate reserve to support leaf and root growth.\n
414
415       ! Maximum time (days) since start of the growing season during which carbohydrate
416       ! may be used
417       IF ( tree(j) ) THEN
418          reserve_time = reserve_time_tree
419       ELSE
420          reserve_time = reserve_time_grass
421       ENDIF
422
423!> @addtogroup Allocation
424       !> Growth is only supported by the use of carbohydrate reserves if the following
425       !> conditions are  statisfied:\n
426       !> - PFT is not senescent;\n
427       !> - LAI must be low (i.e. below ::lai_happy) and\n
428       !> - Day of year of the simulation is in the beginning of the growing season.\n
429       WHERE ( ( biomass(:,j,ileaf) .GT. 0.0 ) .AND. & 
430            ( .NOT. senescence(:,j) ) .AND. &
431            ( lai(:,j) .LT. lai_happy(j) ) .AND. &
432            ( when_growthinit(:,j) .LT. reserve_time ) ) 
433
434          ! Determine the mass from the carbohydrate reserve that can be used (gC/m**2).
435          ! Satisfy the demand or use everything that is available
436          ! (i.e. ::biomass(:,j,icarbres)). Distribute the demand evenly over the time
437          ! required (::tau_leafinit) to develop a minimal canopy from reserves (::lai_happy).
438          use_reserve(:) = &
439               MIN( biomass(:,j,icarbres), &
440               2._r_std * dt/tau_leafinit * lai_happy(j)/ sla(j) )
441
442          ! Distribute the reserve over leaves and fine roots
443          transloc_leaf(:) = L0/(L0+R0) * use_reserve(:)
444          biomass(:,j,ileaf) = biomass(:,j,ileaf) + transloc_leaf(:)
445          biomass(:,j,iroot) = biomass(:,j,iroot) + ( use_reserve(:) - transloc_leaf(:) )
446
447          ! Adjust the carbohydrate reserve mass by accounting for the reserves used during
448          ! this time step
449          biomass(:,j,icarbres) = biomass(:,j,icarbres) - use_reserve(:)
450
451       ELSEWHERE
452
453          transloc_leaf(:) = 0.0
454
455       ENDWHERE
456   
457 !> @addtogroup Allocation
458       !> 2.2 Update leaf age\n
459       !> 2.2.1 Decrease leaf age in youngest class\n
460       ! Adjust the mass of the youngest leaves by the newly grown leaves (gC/m**2)
461       leaf_mass_young(:) = leaf_frac(:,j,1) * lm_old(:,j) + transloc_leaf(:)
462
463       WHERE ( ( transloc_leaf(:) .GT. min_stomate ) .AND. ( leaf_mass_young(:) .GT. min_stomate ) )
464         
465          ! Adjust leaf age by the ratio of leaf_mass_young (t-1)/leaf_mass_young (t)
466          leaf_age(:,j,1) = MAX( zero, leaf_age(:,j,1) * ( leaf_mass_young(:) - transloc_leaf(:) ) / &
467               leaf_mass_young(:) )
468
469       ENDWHERE
470
471 !> @addtogroup Allocation
472       !> 2.2.2 Update leaf mass fraction for the different age classes\n
473       !> Mass fraction in the youngest age class is calculated as the ratio between
474       !> the new mass in the youngest class and the total leaf biomass
475       !> (inc. the new leaves)\n
476       
477       WHERE ( biomass(:,j,ileaf) .GT. min_stomate )
478         
479          leaf_frac(:,j,1) = leaf_mass_young(:) / biomass(:,j,ileaf)
480
481       ENDWHERE
482
483!> @addtogroup Allocation
484       !> Mass fraction in the other classes is calculated as the ratio bewteen
485       !> the current mass in that age and the total leaf biomass
486       !> (inc. the new leaves)\n
487
488       DO m = 2, nleafages ! Loop over # leaf age classes
489
490          WHERE ( biomass(:,j,ileaf) .GT. min_stomate )
491
492             leaf_frac(:,j,m) = leaf_frac(:,j,m) * lm_old(:,j) / biomass(:,j,ileaf)
493
494          ENDWHERE
495
496       ENDDO ! Loop over # leaf age classes
497
498    ENDDO ! loop over # PFTs
499
500!> @addtogroup Allocation
501!> 3. Calculate fractions of biomass production (NPP) to be allocated to the different
502!> biomass components.\n
503!> The fractions of NPP allocated (dimensionless; 1) to the different compartments depend on the
504!> availability of light, water, and nitrogen.
505
506    DO j = 2, nvm ! Loop over # PFTs
507
508       ! Reset values
509       RtoLSR(:)=0
510       LtoLSR(:)=0
511       StoLSR(:)=0
512
513       ! For trees, partitioning between above and belowground sapwood biomass is a function
514       ! of age. For the other PFTs it is prescribed.
515       ! ?? WHERE DO ::alloc_min, ::alloc_mx and ::demi_alloc come from??
516       IF ( tree(j) ) THEN
517
518          alloc_sap_above (:) = alloc_min(j)+(alloc_max(j)-alloc_min(j))*(1.-EXP(-age(:,j)/demi_alloc(j)))
519          !IF (j .EQ. 3) WRITE(*,*) '%allocated above = 'alloc_sap_above(dimensionless; 1),'age = ',age(1,j)
520       
521       ELSE
522         
523          alloc_sap_above(:) = alloc_sap_above_grass
524       
525       ENDIF
526
527!> @addtogroup Allocation
528       !> 3.1 Calculate light, water and proxy for nitrogen stress.\n
529       ! For the limiting factors a low value indicates a strong limitation
530
531       WHERE ( biomass(:,j,ileaf) .GT. min_stomate )
532
533          !> 3.1.1 Light stress\n
534          !> Light stress is a function of the mean lai on the natural part of the grid box
535          !> and of the PFT-specific LAI for agricultural crops. In line with the DGVM, natural
536          !> PFTs in the same gridbox are treated as if they were spatially mixed whereas
537          !> agricultural PFTs are considered to be spatially separated.
538          WHERE( lai_around(:,j) < 10 )
539             
540             limit_L(:) = MAX( 0.1_r_std, EXP( -0.5_r_std * lai_around(:,j) ) )
541         
542          ELSEWHERE
543             
544             limit_L(:) = 0.1_r_std
545         
546          ENDWHERE
547
548!> @addtogroup Allocation
549          !> 3.1.2 Water stress\n
550          !> Water stress is calculated as the weekly moisture availability.\n
551          limit_W(:) = MAX( 0.1_r_std, MIN( 1._r_std, moiavail_week(:,j) ) )
552
553!> @addtogroup Allocation
554          !> 3.1.3 Proxy for nitrogen stress
555          !> The proxy for nitrogen stress depends on monthly soil water availability
556          !> (::soilhum_month) and monthly soil temperature (::tsoil_month). See section
557          !> 1.2.2 for details on how ::t_nitrogen and ::h_nitrogen were calculated.\n
558          !  Currently nitrogen-stress is calculated for both natural and agricultural PFTs.
559          !  Due to intense fertilization of agricultural PFTs this is a strong
560          !  assumption for several agricultural regions in the world (US, Europe, India, ...)
561           
562          ! Water stress on nitrogen mineralisation
563          limit_N_hum(:) = MAX( 0.5_r_std, MIN( 1._r_std, h_nitrogen(:) ) )
564
565          ! Temperature stress on nitrogen mineralisation using a Q10 decomposition model
566          ! where Q10 was set to 2
567          limit_N_temp(:) = 2.**((t_nitrogen(:)-ZeroCelsius-25.)/10.)
568          limit_N_temp(:) = MAX( 0.1_r_std, MIN( 1._r_std, limit_N_temp(:) ) )
569
570          ! Combine water and temperature factors to get total nitrogen stress
571          limit_N(:) = MAX( 0.1_r_std, MIN( 1._r_std, limit_N_hum(:) * limit_N_temp(:) ) )
572
573          ! Take the most limiting factor among soil water and nitrogen
574          limit_WorN(:) = MIN( limit_W(:), limit_N(:) )
575
576          ! Take the most limiting factor among aboveground (i.e. light) and belowground
577          ! (i.e. water & nitrogen) limitations
578          limit(:) = MIN( limit_WorN(:), limit_L(:) )
579
580
581!> @addtogroup Allocation
582          !> 3.2 Calculate ratio between allocation to leaves, sapwood and roots\n
583          !> Partitioning between belowground and aboveground biomass components is assumed
584          !> to be proportional to the ratio of belowground and aboveground stresses.\n
585          !! \latexonly
586          !! \input{alloc1.tex}
587          !! \endlatexonly
588          !! \n
589
590          ! Root allocation
591          RtoLSR(:) = &
592               MAX( .15_r_std, &
593               R0 * 3._r_std * limit_L(:) / ( limit_L(:) + 2._r_std * limit_WorN(:) ) )
594
595          ! Sapwood allocation
596          StoLSR(:) = S0 * 3. * limit_WorN(:) / ( 2._r_std * limit_L(:) + limit_WorN(:) )
597
598          ! Leaf allocation
599          LtoLSR(:) = 1. - RtoLSR(:) - StoLSR(:)
600          LtoLSR(:) = MAX( min_LtoLSR, MIN( max_LtoLSR, LtoLSR(:) ) )
601
602          ! Recalculate roots allocation as the residual carbon
603          RtoLSR(:) = 1. - LtoLSR(:) - StoLSR(:)
604
605       ENDWHERE
606
607       ! Check whether allocation needs to be adjusted. If LAI exceeds maximum LAI
608       ! (::lai_max), no addition carbon should be allocated to leaf biomass. Allocation is
609       ! then partioned between root and sapwood biomass.
610
611       WHERE ( (biomass(:,j,ileaf) .GT. min_stomate) .AND. (lai(:,j) .GT. lai_max(j)) )
612
613          StoLSR(:) = StoLSR(:) + LtoLSR(:)
614          LtoLSR(:) = 0.0
615
616       ENDWHERE
617
618 !> @addtogroup Allocation
619       !> 3.3 Calculate the allocation factors
620       !> The allocation factors (::f_alloc) are an output variable (dimensionless; 1). f_alloc
621       !> has three dimensions (npts,nvm,nparts). Where ::npts is the # of pixels, ::nvm is the
622       !> number of PFTs and ::nparts the number of biomass comoponents. Currently six biomass
623       !> are distinhuished: (dimensionless; 1) Carbon reserves, (2) Aboveground sapwood, (3) Belowground
624       !> sapwood, (4) Roots, (5) fruits/seeds and (6) Leaves.\n
625
626       DO i = 1, npts ! Loop over # pixels - domain size
627
628          IF ( biomass(i,j,ileaf) .GT. min_stomate ) THEN
629     
630             IF ( senescence(i,j) ) THEN
631               
632!> @addtogroup Allocation
633                !> 3.3.1 If the PFT is senescent allocate all C to carbohydrate reserve.\n
634                f_alloc(i,j,icarbres) = 1.0
635
636             ELSE
637
638!> @addtogroup Allocation
639                !> 3.3.2 Allocation during the growing season 
640                f_alloc(i,j,ifruit) = f_fruit
641
642!> @addtogroup Allocation
643                !> Allocation to the carbohydrate reserve is proportional to leaf and root
644                !> allocation. If carbon is allocated to the carbohydtae reserve, rescaling
645                !> of allocation factors is required to ensure carbon mass preservation.\n
646                !> Carbon is allocated to the carbohydrate reserve when the pool size of the
647                !> reserve is less than the carbon needed to grow a canopy twice the size of
648                !> the maximum LAI (::lai_max). Twice the size was used as a threshold because
649                !> the reserves needs to be sufficiently to grow a canopy and roots. In case
650                !> the carbohydrate pool is full, there is no need to rescale the other
651                !> allocation factors.
652                !!?? Where does ecureuil come from??
653                !!?? current equation of (f_alloc(i,j,icarbes) is equivalent to:
654                !!?? reserve alloc = ecureuil*(LtoLSR+StoLSR)*(1-fruit_alloc)*carb_rescale
655                IF ( ( biomass(i,j,icarbres)*sla(j) ) .LT. 2*lai_max(j) ) THEN
656                   !!??NEED to know what is in ecureuil to understand this line??
657                   carb_rescale(i) = 1. / ( 1. + ecureuil(j) * ( LtoLSR(i) + RtoLSR(i) ) )
658                ELSE
659                   carb_rescale(i) = 1.
660                ENDIF
661
662                f_alloc(i,j,ileaf) = LtoLSR(i) * ( 1.-f_alloc(i,j,ifruit) ) * carb_rescale(i)
663                f_alloc(i,j,isapabove) = StoLSR(i) * alloc_sap_above(i) * &
664                     ( 1. - f_alloc(i,j,ifruit) ) * carb_rescale(i)
665                f_alloc(i,j,isapbelow) = StoLSR(i) * ( 1. - alloc_sap_above(i) ) * &
666                     ( 1. - f_alloc(i,j,ifruit) ) * carb_rescale(i)
667                f_alloc(i,j,iroot) = RtoLSR(i) * ( 1.-f_alloc(i,j,ifruit) ) * carb_rescale(i)
668                f_alloc(i,j,icarbres) = ( 1. - carb_rescale(i) ) * ( 1.-f_alloc(i,j,ifruit) )
669
670             ENDIF  ! Is senescent?
671
672          ENDIF  ! There are leaves
673
674       ENDDO  ! Loop over # pixels - domain size
675
676    ENDDO  ! loop over # PFTs
677
678    IF (bavard.GE.3) WRITE(numout,*) 'Leaving alloc'
679
680  END SUBROUTINE alloc
681
682END MODULE stomate_alloc