source: branches/publications/ORCHIDEE-LEAK-r5919/src_stomate/stomate_npp.f90 @ 5925

Last change on this file since 5925 was 2917, checked in by josefine.ghattas, 9 years ago

Vertical soil discretization change: ticket #190
Done mainly by Fuxing Wang and F Cheruy, J Polcher, JL Dufresnes

The parameter HYDROL_SOIL_DEPTH changed into DEPTH_WMAX but it may change again in a later cleaning commit. No change in the discretization for Choisnel.

  • Property svn:keywords set to HeadURL Date Author Revision
File size: 28.3 KB
Line 
1! =================================================================================================================================
2! MODULE          : stomate_npp
3!
4! CONTACT         : orchidee-help _at_ ipsl.jussieu.fr
5!
6! LICENCE         : IPSL (2006)
7!                 This software is governed by the CeCILL licence see ORCHIDEE/ORCHIDEE_CeCILL.LIC
8!
9!>\BRIEF          This modules calculates NPP: Maintenance and growth respiration
10!!
11!!\n DESCRIPTION: We calculate first the maintenance respiration. This is substracted from the
12!!                allocatable biomass (and from the present biomass if the GPP is too low).\n
13!!                Of the rest, a part is lost as growth respiration, while the other part is
14!!                effectively allocated.
15!!
16!! RECENT CHANGE(S): None
17!!
18!! REFERENCE(S) :
19!!
20!! SVN          :
21!! $HeadURL$
22!! $Date$
23!! $Revision$
24!! \n
25!_ ================================================================================================================================
26
27MODULE stomate_npp
28
29  ! modules used:
30  USE xios_orchidee
31  USE ioipsl_para
32  USE stomate_data
33  USE constantes
34  USE constantes_soil
35  USE pft_parameters
36
37  IMPLICIT NONE
38
39  ! private & public routines
40
41  PRIVATE
42  PUBLIC npp_calc,npp_calc_clear
43
44  LOGICAL, SAVE                                              :: firstcall_npp = .TRUE.         !! first call
45!$OMP THREADPRIVATE(firstcall_npp)
46
47CONTAINS
48
49!! ================================================================================================================================
50!! SUBROUTINE   : npp_calc_clear
51!!
52!>\BRIEF        : Set the flag ::firstcall_npp to .TRUE. and as such activate section
53!! 1.1 of the subroutine npp_calc (see below).\n
54!_ ================================================================================================================================
55
56  SUBROUTINE npp_calc_clear
57    firstcall_npp=.TRUE.
58  END SUBROUTINE npp_calc_clear
59
60
61
62
63
64!! ================================================================================================================================
65!! SUBROUTINE   : npp_calc
66!!
67!>\BRIEF        Calculate NPP as the difference between GPP and respiration (= growth + maintenance respiration).
68!!              Update biomass of all compartments after calculating respiration and allocation.
69!!
70!!
71!! DESCRIPTION  : NPP is calculated from three components: Gross Primary Productivity (GPP), maintenance respiration
72!! and growth respiration (all in @tex $ gC.m^{-2}dt^{-1} $ @endtex), following the convention that positive fluxes denote
73!! fluxes plants to the atmosphere. GPP is the input variable from which, in the end, NPP or total allocatable biomass
74!! @tex $(gC.m^{-2}dt^{-1}))$ @endtex is calculated. Net primary production is then calculated as:\n   
75!! NPP = GPP - growth_resp - maint-resp   [eq. 1]\n   
76!!     
77!! The calculation of maintenance respiration is done in routine stomate_resp.f90. Maintenance respiration is calculated for
78!! the whole plant and is therefore removed from the total allocatable biomass. In order to prevent all allocatable biomass
79!! from being used for maintenance respiration, a limit fraction of total allocatable biomass, tax_max, is defined (in
80!! variables declaration). If maintenance respiration exceeds tax_max (::bm_tax_max), the maximum allowed allocatable biomass
81!! will be respired and the remaining respiration, required in excess of tax_max, is taken out from tissues already present in
82!! the plant (biomass).\n 
83!!
84!! After total allocatable biomass has been updated by removing maintenance respiration, total allocatable biomass is distributed
85!! to all plant compartments according to the f_alloc fractions calculated in stomate_alloc.f90.\n
86!!
87!! Growth respiration is calculated as a fraction of allocatable biomass for each part of the plant. The fraction coefficient
88!! ::frac_growth_resp is defined in stomate_constants.f90 and is currently set to be the same for all plant compartments.
89!! Allocatable biomass of all plant compartments are updated by removing what is lost through growth respiration. Net allocatable
90!! biomass (total allocatable biomass after maintenance and growth respiration) is added to the current biomass for  each plant
91!! compartment.
92!!
93!! Finally, leaf age and plant age are updated. Leaf age is described with the concept of "leaf age classes". A number of leaf
94!! age classes (nleafages) is defined in stomate_constants.f90. Each leaf age class contains a fraction (::leaf_frac) of the
95!! total leaf biomass. When new biomass is added to leaves, the age of the biomass in the youngest leaf age class is decreased.
96!! The fractions of leaves in the other leaf ages classes are also updated as the total biomass has increased. Plant age is
97!! updated first by increasing the age of the previous biomass by one time step, and then by adjusting this age as the average
98!! of the ages of the previous and the new biomass.
99!!
100!! RECENT CHANGE(S): None
101!!
102!! MAIN OUTPUT VARIABLE(S): ::npp
103!!
104!! REFERENCE(S) :
105!! - F.W.T.Penning De Vries, A.H.M. Brunsting, H.H. Van Laar. 1974. Products, requirements and efficiency of biosynthesis a
106!! quantitative approach. Journal of Theoretical Biology, Volume 45, Issue 2, June 1974, Pages 339-377.
107!!
108!! FLOWCHART :
109!! \latexonly
110!! \includegraphics[scale=0.14]{stomate_npp_flow.jpg}
111!! \endlatexonly
112!! \n
113!_ ================================================================================================================================
114
115  SUBROUTINE npp_calc (npts, dt, &
116       PFTpresent, &
117       t2m, tsoil, lai, rprof, &
118       gpp, f_alloc, bm_alloc, resp_maint_part,&
119       biomass, leaf_age, leaf_frac, age, &
120       resp_maint, resp_growth, npp)
121   
122!! 0 Variable and parameter declaration
123
124    !! 0.1 Input variables
125
126    INTEGER(i_std), INTENT(in)                                :: npts             !! Domain size - number of pixels (unitless)
127    REAL(r_std), INTENT(in)                                   :: dt               !! Time step (days)
128    LOGICAL, DIMENSION(npts,nvm), INTENT(in)                  :: PFTpresent       !! PFT exists (true/false)
129    REAL(r_std), DIMENSION(npts), INTENT(in)                  :: t2m              !! Temperature at 2 meter (K)
130    REAL(r_std), DIMENSION(npts,nbdl), INTENT(in)             :: tsoil            !! Soil temperature of each soil layer (K)
131    REAL(r_std), DIMENSION(npts,nvm), INTENT(in)              :: lai              !! PFT leaf area index (unitless)
132    REAL(r_std), DIMENSION(npts,nvm), INTENT(in)              :: rprof            !! PFT root depth as calculated in stomate.f90
133                                                                                  !! from root profile parameter humcste (m)
134    REAL(r_std), DIMENSION(npts,nvm), INTENT(in)              :: gpp              !! PFT gross primary productivity
135                                                                                  !! @tex $(gC.m^{-2}dt^{-1})$ @endtex
136    REAL(r_std), DIMENSION(npts,nvm,nparts), INTENT(in)       :: f_alloc          !! Fraction of total allocatable biomass that
137                                                                                  !! goes into each plant part (unitless)
138    REAL(r_std), DIMENSION(npts,nvm,nparts), INTENT(in)       :: resp_maint_part  !! Maintenance respiration of different plant
139                                                                                  !! parts @tex $(gC.m^{-2}dt^{-1})$ @endtex
140    !! 0.2 Output variables
141
142    REAL(r_std), DIMENSION(npts,nvm), INTENT(out)             :: resp_maint       !! PFT maintenance respiration
143                                                                                  !! @tex $(gC.m^{-2}dt^{-1})$ @endtex             
144    REAL(r_std), DIMENSION(npts,nvm), INTENT(out)             :: resp_growth      !! PFT growth respiration
145                                                                                  !! @tex $(gC.m^{-2}dt^{-1})$ @endtex                         
146    REAL(r_std), DIMENSION(npts,nvm), INTENT(out)             :: npp              !! PFT net primary productivity
147                                                                                  !! @tex $(gC.m^{-2}dt^{-1})$ @endtex         
148    REAL(r_std), DIMENSION(npts,nvm,nparts,nelements), INTENT(out) :: bm_alloc    !! PFT biomass increase, i.e. NPP per plant part
149                                                                                  !! @tex $(gC.m^{-2}dt^{-1})$ @endtex         
150
151    !! 0.3 Modified variables
152
153    REAL(r_std), DIMENSION(npts,nvm,nparts,nelements), INTENT(inout) :: biomass   !! PFT total biomass of each plant part
154                                                                                  !! @tex $(gC.m^{-2})$ @endtex
155    REAL(r_std), DIMENSION(npts,nvm,nleafages), INTENT(inout) :: leaf_age         !! PFT age of different leaf age classes (days)
156    REAL(r_std), DIMENSION(npts,nvm,nleafages), INTENT(inout) :: leaf_frac        !! PFT fraction of total leaves in leaf age
157                                                                                  !! class (unitless)
158    REAL(r_std), DIMENSION(npts,nvm), INTENT(inout)           :: age              !! PFT age (years)
159
160    !! 0.4 Local variables
161
162    REAL(r_std), SAVE, ALLOCATABLE, DIMENSION(:)              :: z_soil           !! Soil levels  representing soil depth (m)
163!$OMP THREADPRIVATE(z_soil)
164    REAL(r_std), DIMENSION(npts,nvm)                          :: t_root           !! Root temperature (convolution of root and
165                                                                                  !! soil temperature profiles)(K)
166    REAL(r_std), DIMENSION(npts,nvm,nparts)                   :: coeff_maint      !! PFT maintenance respiration coefficients of
167                                                                                  !! different plant compartments at 0 deg C
168                                                                                  !! @tex $(g.g^{-1}dt^{-1})$ @endtex
169    REAL(r_std), DIMENSION(npts,nparts)                       :: t_maint          !! Temperature which is pertinent for maintenance
170                                                                                  !! respiration, which is air/root temperature for
171                                                                                  !! above/below-ground compartments (K)
172    REAL(r_std), DIMENSION(npts)                              :: rpc              !! Scaling factor for integrating vertical soil
173                                                                                  !! profiles (unitless)
174    REAL(r_std), DIMENSION(npts)                              :: tl               !! Long term annual mean temperature (C)
175    REAL(r_std), DIMENSION(npts)                              :: slope            !! Slope of maintenance respiration coefficient
176                                                                                  !! (1/K)
177    REAL(r_std), DIMENSION(npts,nparts)                       :: resp_growth_part !! Growth respiration of different plant parts
178                                                                                  !! @tex $(gC.m^{-2}dt^{-1})$ @endtex         
179    REAL(r_std), DIMENSION(npts,nvm)                          :: bm_alloc_tot     !! Allocatable biomass for the whole plant
180                                                                                  !! @tex $(gC.m^{-2})$ @endtex
181    REAL(r_std), DIMENSION(npts)                              :: bm_add           !! Biomass increase @tex $(gC.m^{-2})$ @endtex               
182    REAL(r_std), DIMENSION(npts)                              :: bm_new           !! New biomass @tex $(gC.m^{-2})$ @endtex     
183    REAL(r_std), DIMENSION(npts,nvm)                          :: leaf_mass_young  !! Leaf mass in youngest age class
184                                                                                  !! @tex $(gC.m^{-2})$ @endtex         
185    REAL(r_std), DIMENSION(npts,nvm)                          :: lm_old           !! Leaf mass after maintenance respiration
186                                                                                  !! @tex $(gC.m^{-2})$ @endtex                 
187    REAL(r_std), DIMENSION(npts,nvm)                          :: bm_create        !! Biomass created when biomass<0 because of dark
188                                                                                  !! respiration @tex $(gC.m^{-2})$ @endtex
189    REAL(r_std), DIMENSION(npts)                              :: bm_tax_max       !! Maximum part of allocatable biomass used for
190                                                                                  !! respiration @tex $(gC.m^{-2})$ @endtex     
191    REAL(r_std), DIMENSION(npts)                              :: bm_pump          !! Biomass that remains to be taken away
192                                                                                  !! @tex $(gC.m^{-2})$ @endtex
193    INTEGER(i_std)                                            :: i,j,k,l,m        !! Indeces(unitless)
194    INTEGER(i_std)                                            :: ier              !! Error handling
195
196!_ ================================================================================================================================
197
198    IF (printlev>=3) WRITE(numout,*) 'Entering npp'
199   
200 !! 1. Initializations
201   
202    !! 1.1 First call
203    IF ( firstcall_npp ) THEN
204
205       !! 1.1.1 Soil levels
206       !  Get the depth of the different soil layers (number of layers=nbdl)
207       !  previously calculated as variable diaglev in routines sechiba.f90 and slowproc.f90 
208       ALLOCATE(z_soil(0:nbdl), stat=ier)
209       IF ( ier /= 0 ) CALL ipslerr_p(3,'npp_calc','Pb in allocate of z_soil','','')
210
211       z_soil(0) = zero
212       z_soil(1:nbdl) = diaglev(1:nbdl)
213
214       !! 1.1.2 Output message
215       !  Write message including value used for tax_max       
216       WRITE(numout,*) 'npp:'
217
218       WRITE(numout,*) '   > max. fraction of allocatable biomass used for'// &
219            ' maint. resp.:', tax_max
220
221       firstcall_npp = .FALSE.
222
223    ENDIF ! End if first call
224
225    !! 1.2 Set output variables to zero
226    bm_alloc(:,:,:,:) = zero
227    resp_maint(:,:) = zero
228    resp_growth(:,:) = zero
229    npp(:,:) = zero
230
231    !! 1.3 Total allocatable biomass
232    ! total allocatable biomass during this time step determined from GPP.
233    ! GPP was calculated as CO2 assimilation in enerbil.f90
234    bm_alloc_tot(:,:) = gpp(:,:) * dt
235
236   
237 
238    !! 3. Calculate maintenance and growth respiration
239    ! First, total maintenance respiration for the whole plant is calculated by summing maintenance
240    ! respiration of the different plant compartments. Then, maintenance respiration is subtracted
241    ! from whole-plant allocatable biomass (up to a maximum fraction of the total allocatable biomass).
242    ! Growth respiration is then calculated for each plant compartment as a fraction of remaining
243    ! allocatable biomass for this compartment. NPP is calculated by substracting total autotrophic
244    ! respiration from GPP i.e. NPP = GPP - maintenance resp - growth resp.
245    DO j = 2,nvm        ! Loop over # of PFTs
246
247       !! 3.1 Maintenance respiration of the different plant parts
248       !      Maintenance respiration of the different plant parts is calculated in
249       !      stomate_resp.f90 as a function of the plant's temperature,
250       !      the long term temperature and plant coefficients
251       !      VPP killer:
252       resp_maint(:,j) = zero
253
254       !  Following the calculation of hourly maintenance respiration, verify that
255       !  the PFT has not been killed after calcul of resp_maint_part in stomate.
256       DO k= 1, nparts
257          WHERE (PFTpresent(:,j))
258             resp_maint(:,j) = resp_maint(:,j) + resp_maint_part(:,j,k)
259          ENDWHERE
260       ENDDO
261       
262       !! 3.2 Substract maintenance respiration from allocatable biomass
263       !      The total maintenance respiration calculated in 3.2 is substracted  from the newly
264       !      produced allocatable biomass (bm_alloc_tot). However, ensure that not all allocatable
265       !      biomass is removed by setting a maximum to the fraction of allocatable biomass used
266       !      for maintenance respiration: tax_max. If the maintenance respiration is larger than
267       !      tax_max,the amount tax_max is taken from allocatable biomass, and the remaining of
268       !      maintenance respiration is taken from the tissues themselves (biomass). We suppose
269       !      that respiration is not dependent on leaf age -> therefore the leaf age structure is
270       !      not changed.
271       !      The maximum fraction of allocatable biomass used for respiration is defined as tax_max.
272       !      The value of tax_max is set in the declarations section (0.4 Local variables) of this
273       !      routine
274       bm_tax_max(:) = tax_max * bm_alloc_tot(:,j)
275
276       DO i = 1, npts   ! Loop over # of pixels
277
278          ! If there is enough allocatable biomass to cover maintenance respiration,
279          ! then biomass associated with maintenance respiration is removed from allocatable biomass
280          IF ( bm_alloc_tot(i,j) .GT. zero ) THEN
281               IF ( ( resp_maint(i,j) * dt ) .LT. bm_tax_max(i) )  THEN
282       
283                  bm_alloc_tot(i,j) = bm_alloc_tot(i,j) - resp_maint(i,j) * dt
284
285                  ! If there is not enough allocatable biomass to cover maintenance respiration, the 
286                  ! - maximum allowed allocatable biomass (bm_tax_max) is removed from allocatable biomass.
287               ELSE
288             
289                  bm_alloc_tot(i,j) = bm_alloc_tot(i,j) - bm_tax_max(i)
290
291                  ! ::bm_pump is the amount of maintenance respiration that exceeds the maximum allocatable biomass
292                  ! This amount of biomass still needs to be respired and will be removed from tissues biomass of each
293                  ! plant compartment
294                  bm_pump(i) = resp_maint(i,j) * dt - bm_tax_max(i)
295
296                  ! The biomass is removed from each plant compartment tissues as the ratio of the maintenance         
297                  ! respiration of the plant compartment to the total maintenance respiration (resp_maint_part/resp_maint)
298                  biomass(i,j,ileaf,icarbon) = biomass(i,j,ileaf,icarbon) - &
299                       bm_pump(i) * resp_maint_part(i,j,ileaf) / resp_maint(i,j)
300                  biomass(i,j,isapabove,icarbon) = biomass(i,j,isapabove,icarbon) - &
301                       bm_pump(i) * resp_maint_part(i,j,isapabove) / resp_maint(i,j)
302                  biomass(i,j,isapbelow,icarbon) = biomass(i,j,isapbelow,icarbon) - &
303                       bm_pump(i) * resp_maint_part(i,j,isapbelow) / resp_maint(i,j)
304                  biomass(i,j,iroot,icarbon) = biomass(i,j,iroot,icarbon) - &
305                       bm_pump(i) * resp_maint_part(i,j,iroot) / resp_maint(i,j)
306                  biomass(i,j,ifruit,icarbon) = biomass(i,j,ifruit,icarbon) - &
307                       bm_pump(i) * resp_maint_part(i,j,ifruit) / resp_maint(i,j)
308                  biomass(i,j,icarbres,icarbon) = biomass(i,j,icarbres,icarbon) - &
309                       bm_pump(i) * resp_maint_part(i,j,icarbres) / resp_maint(i,j)
310               ENDIF
311          ELSE
312             biomass(i,j,icarbres,icarbon) = biomass(i,j,icarbres,icarbon) - & 
313                  bm_alloc_tot(i,j) - resp_maint(i,j) * dt 
314             bm_alloc_tot(i,j) = 0. 
315          ENDIF ! End if there is enough allocatable biomass to cover maintenance respiration
316
317       ENDDO   ! Fortran95: WHERE - ELSEWHERE construct
318
319       
320       !! 3.3 Allocate allocatable biomass to different plant compartments.
321       !      The amount of allocatable biomass of each compartment is a fraction according f_alloc of total
322       !      allocatable biomass (the f_alloc of the different plant parts are calculated in stomate_alloc.f90)
323       DO k = 1, nparts
324          bm_alloc(:,j,k,icarbon) = f_alloc(:,j,k) * bm_alloc_tot(:,j)
325       ENDDO
326
327       
328       !! 3.4 Calculate growth respiration of each plant compartment.
329       !      Growth respiration of a plant compartment is a fraction of the allocatable biomass remaining after
330       !      maintenance respiration losses have been taken into account. The fraction of allocatable biomass
331       !      removed for growth respiration is the same for all plant compartments and is defined by the parameter
332       !      frac_growth_resp in stomate_constants.f90. Allocatable biomass ::bm_alloc is updated as a result of
333       !      the removal of growth resp.
334       resp_growth_part(:,:) = frac_growthresp(j) * bm_alloc(:,j,:,icarbon) / dt
335       bm_alloc(:,j,:,icarbon) = ( un - frac_growthresp(j) ) * bm_alloc(:,j,:,icarbon)
336
337       
338       !! 3.5 Total growth respiration
339       !      Calculate total growth respiration of the plant as the sum of growth respiration of all plant parts       
340       resp_growth(:,j) = zero
341
342       DO k = 1, nparts
343          resp_growth(:,j) = resp_growth(:,j) + resp_growth_part(:,k)
344       ENDDO
345
346    ENDDO ! # End Loop over # of PFTs
347
348   
349 !! 4. Update the biomass with newly allocated biomass after respiration
350 
351    !  Save the old leaf biomass for later. "old" leaf mass is leaf mass after maintenance respiration in the case
352    !  where maintenance respiration has required taking biomass from tissues in section 3.3
353    lm_old(:,:) = biomass(:,:,ileaf,icarbon)
354    biomass(:,:,:,:) = biomass(:,:,:,:) + bm_alloc(:,:,:,:)
355
356   
357 !! 5. Deal with negative biomasses
358   
359    !  Biomass can become negative in some rare cases, as the GPP can be negative. This corresponds to very
360    !  situations that can be seen as the 'creation' of a seed ('virtual photosynthesis'). In this case, we set
361    !  biomass to some small value min_stomate. For carbon budget to remain balanced, this creation of matter (carbon)
362    !  is taken into account by decreasing the autotrophic respiration by the same amount that has been added to biomass
363    !  for it to become positive. In this case, maintenance respiration can become negative in extreme cases (deserts)!!
364
365    DO k = 1, nparts    ! Loop over # of plant parts
366
367       DO j = 2,nvm     ! Loop over # of PFTs
368
369          WHERE ( biomass(:,j,k,icarbon) .LT. zero )
370
371             bm_create(:,j) = min_stomate - biomass(:,j,k,icarbon)
372             
373             biomass(:,j,k,icarbon) = biomass(:,j,k,icarbon) + bm_create(:,j)
374             
375             resp_maint(:,j) = resp_maint(:,j) - bm_create(:,j) / dt
376
377          ENDWHERE
378
379       ENDDO    ! Loop over # of PFTs
380
381    ENDDO       ! Loop over # plant parts
382
383   
384 !! 6. Calculate NPP (See Eq 1 in header)
385   
386    !  Calculate the NPP @tex $(gC.m^{-2}dt^{-1})$ @endtex as the difference between GPP
387    !  and autotrophic respiration (maintenance and growth respirations)
388    DO j = 2,nvm        ! Loop over # PFTs
389       npp(:,j) = gpp(:,j) - resp_growth(:,j) - resp_maint(:,j)
390    ENDDO       ! Loop over # PFTs
391
392   
393 !! 7. Update leaf age
394
395    !  Leaf age is needed for calculation of turnover and vmax in stomate_turnover.f90 and stomate_vmax.f90 routines.
396    !  Leaf biomass is distributed according to its age into several "age classes" with age class=1 representing the
397    !  youngest class, and consisting of the most newly allocated leaf biomass
398   
399    !! 7.1 Update quantity and age of the leaf biomass in the youngest class
400    !      The new amount of leaf biomass in the youngest age class (leaf_mass_young) is the sum of :
401    !      - the leaf biomass that was already in the youngest age class (leaf_frac(:,j,1) * lm_old(:,j)) with the
402    !        leaf age given in leaf_age(:,j,1)
403    !      - and the new biomass allocated to leaves (bm_alloc(:,j,ileaf)) with a leaf age of zero.
404    DO j = 2,nvm
405       leaf_mass_young(:,j) = leaf_frac(:,j,1) * lm_old(:,j) + bm_alloc(:,j,ileaf,icarbon)
406    ENDDO
407
408    ! The age of the updated youngest age class is the average of the ages of its 2 components: bm_alloc(leaf) of age
409    ! '0', and leaf_frac*lm_old(=leaf_mass_young-bm_alloc) of age 'leaf_age(:,j,1)'
410    DO j = 2,nvm
411       WHERE ( ( bm_alloc(:,j,ileaf,icarbon) .GT. zero ) .AND. &
412         ( leaf_mass_young(:,j) .GT. zero ) )
413
414          leaf_age(:,j,1) = MAX ( zero, &
415               & leaf_age(:,j,1) * &
416               & ( leaf_mass_young(:,j) - bm_alloc(:,j,ileaf,icarbon) ) / &
417               & leaf_mass_young(:,j) )
418         
419       ENDWHERE
420    ENDDO
421
422    !! 7.2 Update leaf age
423    !      Update fractions of leaf biomass in each age class (fraction in youngest class increases)
424
425    !! 7.2.1 Update age of youngest leaves
426    !        For age class 1 (youngest class), because we have added biomass to the youngest class, we need to update
427    !        the fraction of total leaf biomass that belongs to the youngest age class : updated mass in class divided
428    !        by new total leaf mass
429    DO j = 2,nvm
430       WHERE ( biomass(:,j,ileaf,icarbon) .GT. min_stomate )
431
432          leaf_frac(:,j,1) = leaf_mass_young(:,j) / biomass(:,j,ileaf,icarbon)
433
434       ENDWHERE
435    ENDDO
436
437    !! 7.2.2 Update age of other age classes
438    !        Because the total leaf biomass has changed, we need to update the fraction of leaves in each age class:
439    !        mass in leaf age class (from previous fraction of leaves in this class and previous total leaf biomass)
440    !        divided by new total mass
441    DO m = 2, nleafages ! Loop over # leaf age classes
442
443       DO j = 2,nvm     ! Loop over # PFTs
444          WHERE ( biomass(:,j,ileaf,icarbon) .GT. min_stomate )
445
446             leaf_frac(:,j,m) = leaf_frac(:,j,m) * lm_old(:,j) / biomass(:,j,ileaf,icarbon)
447
448          ENDWHERE
449       ENDDO
450
451    ENDDO       ! Loop over # leaf age classes
452
453 !! 8. Update whole-plant age
454   
455    !! 8.1 PFT age
456    !      At every time step, increase age of the biomass that was already present at previous time step.
457    !      Age is expressed in years, and the time step 'dt' in days so age increase is: dt divided by number
458    !      of days in a year.
459    WHERE ( PFTpresent(:,:) )
460
461       age(:,:) = age(:,:) + dt/one_year
462
463    ELSEWHERE
464
465       age(:,:) = zero
466
467    ENDWHERE
468
469    !! 8.2 Age of grasses and crops
470    !  For grasses and crops, biomass with age 0 has been added to the whole plant with age 'age'. New biomass is the sum of
471    !  the current total biomass in all plant parts (bm_new), bm_new(:) = SUM( biomass(:,j,:), DIM=2 ). The biomass that has
472    !  just been added is the sum of the allocatable biomass of all plant parts (bm_add), its age is zero. bm_add(:) =
473    !  SUM( bm_alloc(:,j,:), DIM=2 ). Before allocation, the plant biomass is bm_new-bm_add, its age is "age(:,j)". The age of
474    !  the new biomass is the average of the ages of previous and added biomass.
475    !  For trees, age is treated in "establish" if vegetation is dynamic, and in turnover routines if it is static (in this
476    !  case, only the age of the heartwood is accounted for).
477    DO j = 2,nvm
478
479       IF ( .NOT. is_tree(j) ) THEN
480
481          bm_new(:) = biomass(:,j,ileaf,icarbon) + biomass(:,j,isapabove,icarbon) + &
482               biomass(:,j,iroot,icarbon) + biomass(:,j,ifruit,icarbon)
483          bm_add(:) = bm_alloc(:,j,ileaf,icarbon) + bm_alloc(:,j,isapabove,icarbon) + &
484               bm_alloc(:,j,iroot,icarbon) + bm_alloc(:,j,ifruit,icarbon)
485
486          WHERE ( ( bm_new(:) .GT. zero ) .AND. ( bm_add(:) .GT. zero ) )
487             age(:,j) = age(:,j) * ( bm_new(:) - bm_add(:) ) / bm_new(:)
488          ENDWHERE
489
490       ENDIF
491
492    ENDDO
493
494 !! 9. Write history files
495
496    CALL xios_orchidee_send_field("BM_ALLOC_LEAF",bm_alloc(:,:,ileaf,icarbon))
497    CALL xios_orchidee_send_field("BM_ALLOC_SAP_AB",bm_alloc(:,:,isapabove,icarbon))
498    CALL xios_orchidee_send_field("BM_ALLOC_SAP_BE",bm_alloc(:,:,isapbelow,icarbon))
499    CALL xios_orchidee_send_field("BM_ALLOC_ROOT",bm_alloc(:,:,iroot,icarbon))
500    CALL xios_orchidee_send_field("BM_ALLOC_FRUIT",bm_alloc(:,:,ifruit,icarbon))
501    CALL xios_orchidee_send_field("BM_ALLOC_RES",bm_alloc(:,:,icarbres,icarbon))
502
503
504    ! Save in history file the variables describing the biomass allocated to the plant parts
505    CALL histwrite_p (hist_id_stomate, 'BM_ALLOC_LEAF', itime, &
506         bm_alloc(:,:,ileaf,icarbon), npts*nvm, horipft_index)
507    CALL histwrite_p (hist_id_stomate, 'BM_ALLOC_SAP_AB', itime, &
508         bm_alloc(:,:,isapabove,icarbon), npts*nvm, horipft_index)
509    CALL histwrite_p (hist_id_stomate, 'BM_ALLOC_SAP_BE', itime, &
510         bm_alloc(:,:,isapbelow,icarbon), npts*nvm, horipft_index)
511    CALL histwrite_p (hist_id_stomate, 'BM_ALLOC_ROOT', itime, &
512         bm_alloc(:,:,iroot,icarbon), npts*nvm, horipft_index)
513    CALL histwrite_p (hist_id_stomate, 'BM_ALLOC_FRUIT', itime, &
514         bm_alloc(:,:,ifruit,icarbon), npts*nvm, horipft_index)
515    CALL histwrite_p (hist_id_stomate, 'BM_ALLOC_RES', itime, &
516         bm_alloc(:,:,icarbres,icarbon), npts*nvm, horipft_index)
517
518
519    IF (printlev>=4) WRITE(numout,*) 'Leaving npp'
520
521  END SUBROUTINE npp_calc
522
523END MODULE stomate_npp
Note: See TracBrowser for help on using the repository browser.