source: tags/ORCHIDEE/src_stomate/stomate_npp.f90 @ 6

Last change on this file since 6 was 6, checked in by orchidee, 14 years ago

import first tag equivalent to CVS orchidee_1_9_5 + OOL_1_9_5
File size: 15.9 KB
Line 
1! Npp: Maintenance and growth respiration
2! We calculte first the maintenance rspiration. This is substracted from the
3! allocatable biomass (and from the present biomass if the GPP is too low).
4! Of the rest, a part is lost as growth respiration, while the other part is
5! effectively allocated.
6!
7! $Header: /home/ssipsl/CVSREP/ORCHIDEE/src_stomate/stomate_npp.f90,v 1.14 2010/04/20 14:12:04 ssipsl Exp $
8! IPSL (2006)
9!  This software is governed by the CeCILL licence see ORCHIDEE/ORCHIDEE_CeCILL.LIC
10!
11MODULE stomate_npp
12
13  ! modules used:
14
15  USE ioipsl
16  USE stomate_constants
17  USE constantes_veg
18
19  IMPLICIT NONE
20
21  ! private & public routines
22
23  PRIVATE
24  PUBLIC npp_calc,npp_calc_clear
25
26  ! first call
27  LOGICAL, SAVE                                              :: firstcall = .TRUE.
28
29CONTAINS
30
31  SUBROUTINE npp_calc_clear
32    firstcall=.TRUE.
33  END SUBROUTINE npp_calc_clear
34
35  SUBROUTINE npp_calc (npts, dt, &
36       PFTpresent, &
37       tlong_ref, t2m, tsoil, lai, rprof, &
38       gpp, f_alloc, bm_alloc, resp_maint_part,&
39       biomass, leaf_age, leaf_frac, age, &
40       resp_maint, resp_growth, npp)
41
42    !
43    ! 0 declarations
44    !
45
46    ! 0.1 input
47
48    ! Domain size
49    INTEGER(i_std), INTENT(in)                                        :: npts
50    ! time step (days)
51    REAL(r_std), INTENT(in)                                     :: dt
52    ! PFT exists
53    LOGICAL, DIMENSION(npts,nvm), INTENT(in)                  :: PFTpresent
54    ! long term annual mean 2 meter reference temperature
55    REAL(r_std), DIMENSION(npts), INTENT(in)                    :: tlong_ref
56    ! 2 meter temperature
57    REAL(r_std), DIMENSION(npts), INTENT(in)                    :: t2m
58    ! soil temperature (K)
59    REAL(r_std), DIMENSION(npts,nbdl), INTENT(in)               :: tsoil
60    ! leaf area index
61    REAL(r_std), DIMENSION(npts,nvm), INTENT(in)               :: lai
62    ! root depth (m)
63    REAL(r_std), DIMENSION(npts,nvm), INTENT(in)               :: rprof
64    ! gross primary productivity (gC/days/(m**2 of ground))
65    REAL(r_std), DIMENSION(npts,nvm), INTENT(in)               :: gpp
66    ! fraction that goes into plant part
67    REAL(r_std), DIMENSION(npts,nvm,nparts), INTENT(in)        :: f_alloc
68    ! maintenance respiration of different plant parts (gC/day/m**2 of ground)
69    REAL(r_std), DIMENSION(npts,nvm,nparts), INTENT(in)             :: resp_maint_part
70
71    ! 0.2 modified fields
72
73    ! biomass (gC/(m**2 of ground))
74    REAL(r_std), DIMENSION(npts,nvm,nparts), INTENT(inout)     :: biomass
75    ! leaf age (days)
76    REAL(r_std), DIMENSION(npts,nvm,nleafages), INTENT(inout)  :: leaf_age
77    ! fraction of leaves in leaf age class
78    REAL(r_std), DIMENSION(npts,nvm,nleafages), INTENT(inout)  :: leaf_frac
79    ! age (years)
80    REAL(r_std), DIMENSION(npts,nvm), INTENT(inout)            :: age
81
82    ! 0.3 output
83
84    ! maintenance respiration (gC/day/m**2 of total ground)
85    REAL(r_std), DIMENSION(npts,nvm), INTENT(out)              :: resp_maint
86    ! autotrophic respiration (gC/day/m**2 of total ground)
87    REAL(r_std), DIMENSION(npts,nvm), INTENT(out)              :: resp_growth
88    ! net primary productivity (gC/day/m**2 of ground)
89    REAL(r_std), DIMENSION(npts,nvm), INTENT(out)              :: npp
90    ! biomass increase, i.e. NPP per plant part
91    REAL(r_std), DIMENSION(npts,nvm,nparts), INTENT(out)       :: bm_alloc
92
93    ! 0.4 local
94
95    ! maximum fraction of allocatable biomass used for maintenance respiration
96    REAL(r_std), PARAMETER                                      :: tax_max = 0.8
97    ! soil levels (m)
98    REAL(r_std), SAVE, DIMENSION(0:nbdl)                        :: z_soil
99    ! root temperature (convolution of root and soil temperature profiles)
100    REAL(r_std), DIMENSION(npts,nvm)                           :: t_root
101    ! maintenance respiration coefficients at 0 deg C (g/g d**-1)
102    REAL(r_std), DIMENSION(npts,nvm,nparts)                    :: coeff_maint
103    ! temperature which is pertinent for maintenance respiration (K)
104    REAL(r_std), DIMENSION(npts,nparts)                         :: t_maint
105    ! integration constant for root profile
106    REAL(r_std), DIMENSION(npts)                                :: rpc
107    ! long term annual mean temperature, C
108    REAL(r_std), DIMENSION(npts)                                :: tl
109    ! slope of maintenance respiration coefficient (1/K)
110    REAL(r_std), DIMENSION(npts)                                :: slope
111    ! growth respiration of different plant parts (gC/day/m**2 of ground)
112    REAL(r_std), DIMENSION(npts,nparts)                         :: resp_growth_part
113    ! allocatable biomass (gC/m**2 of ground) for the whole plant
114    REAL(r_std), DIMENSION(npts,nvm)                           :: bm_alloc_tot
115    ! biomass increase
116    REAL(r_std), DIMENSION(npts)                                :: bm_add
117    ! new biomass
118    REAL(r_std), DIMENSION(npts)                                :: bm_new
119    ! leaf mass in youngest age class (gC/m**2 of ground)
120    REAL(r_std), DIMENSION(npts,nvm)                           :: leaf_mass_young
121    ! leaf mass after maintenance respiration
122    REAL(r_std), DIMENSION(npts,nvm)                           :: lm_old
123    ! biomass created when biomass<0 because of dark respiration (gC/m**2 of ground)
124    REAL(r_std), DIMENSION(npts,nvm)                           :: bm_create
125    ! maximum part of allocatable biomass used for respiration
126    REAL(r_std), DIMENSION(npts)                                :: bm_tax_max
127    ! biomass that remains to be taken away
128    REAL(r_std), DIMENSION(npts)                                :: bm_pump
129    ! Index
130    INTEGER(i_std)                                             :: i,j,k,l,m
131
132    ! =========================================================================
133
134    IF (bavard.GE.3) WRITE(numout,*) 'Entering npp'
135    !
136    ! 1 Initializations
137    !
138
139    !
140    ! 1.1 first call
141    !
142
143    IF ( firstcall ) THEN
144
145       ! 1.1.1 soil levels
146
147       z_soil(0) = 0.
148       z_soil(1:nbdl) = diaglev(1:nbdl)
149
150       ! 1.1.2 messages
151
152       WRITE(numout,*) 'npp:'
153
154       WRITE(numout,*) '   > max. fraction of allocatable biomass used for'// &
155            ' maint. resp.:', tax_max
156
157       firstcall = .FALSE.
158
159    ENDIF
160
161    !
162    ! 1.2 set output to zero
163    !
164
165    bm_alloc(:,:,:) = zero
166    resp_maint(:,:) = zero
167    resp_growth(:,:) = zero
168    npp(:,:) = zero
169
170    !
171    ! 1.3 root temperature: convolution of root and temperature profiles
172    !     suppose exponential root profile.
173    !
174
175    DO j = 2,nvm
176
177       ! 1.3.1 rpc is an integration constant such that the integral of the root profile is 1.
178       rpc(:) = 1. / ( 1. - EXP( -z_soil(nbdl) / rprof(:,j) ) )
179
180       ! 1.3.2 integrate over the nbdl levels
181
182       t_root(:,j) = zero
183
184       DO l = 1, nbdl
185
186          t_root(:,j) = &
187               t_root(:,j) + tsoil(:,l) * rpc(:) * &
188               ( EXP( -z_soil(l-1)/rprof(:,j) ) - EXP( -z_soil(l)/rprof(:,j) ) )
189
190       ENDDO
191
192    ENDDO
193
194    !
195    ! 1.4 total allocatable biomass
196    !
197
198    bm_alloc_tot(:,:) = gpp(:,:) * dt
199
200    !
201    ! 2 define maintenance respiration coefficients
202    !
203
204    DO j = 2,nvm
205
206       !
207       ! 2.1 temperature which is taken for the plant part we are talking about
208       !
209
210       ! 2.1.1 parts above the ground
211
212       t_maint(:,ileaf) = t2m(:)
213       t_maint(:,isapabove) = t2m(:)
214       t_maint(:,ifruit) = t2m(:)
215
216       ! 2.1.2 parts below the ground
217
218       t_maint(:,isapbelow) = t_root(:,j)
219       t_maint(:,iroot) = t_root(:,j)
220
221       ! 2.1.3 heartwood: does not respire. Any temperature
222
223       t_maint(:,iheartbelow) = t_root(:,j)
224       t_maint(:,iheartabove) = t2m(:)
225
226       ! 2.1.4 reserve: above the ground for trees, below for grasses
227
228       IF ( tree(j) ) THEN
229          t_maint(:,icarbres) = t2m(:)
230       ELSE
231          t_maint(:,icarbres) = t_root(:,j)
232       ENDIF
233
234       !
235       ! 2.2 calculate coefficient
236       !
237
238       tl(:) = tlong_ref(:) - ZeroCelsius
239       slope(:) = maint_resp_slope(j,1) + tl(:) * maint_resp_slope(j,2) + &
240            tl(:)*tl(:) * maint_resp_slope(j,3)
241
242       DO k = 1, nparts
243
244          coeff_maint(:,j,k) = &
245               MAX( coeff_maint_zero(j,k) * &
246               ( 1. + slope(:) * (t_maint(:,k)-ZeroCelsius) ), zero )
247
248       ENDDO
249
250    ENDDO
251
252    !
253    ! 3 calculate maintenance and growth respiration.
254    !   NPP = GPP - maintenance resp - growth resp.
255    !
256    !
257    DO j = 2,nvm
258
259       !
260       ! 3.1 maintenance respiration of the different plant parts
261       !
262       !
263       ! 3.2 Total maintenance respiration of the plant
264       !     VPP killer:
265       !     resp_maint(:,j) = SUM( resp_maint_part(:,:), DIM=2 )
266       !
267
268       resp_maint(:,j) = zero
269
270       ! with the new calculation of hourly respiration, we must verify that
271       ! PFT has not been killed after calcul of resp_maint_part in stomate
272       DO k= 1, nparts
273          WHERE (PFTpresent(:,j))
274             resp_maint(:,j) = resp_maint(:,j) + resp_maint_part(:,j,k)
275          ENDWHERE
276       ENDDO
277       !
278       ! 3.3 This maintenance respiration is taken away from the newly produced
279       !     allocatable biomass. However, we avoid that no allocatable biomass remains.
280       !     If the respiration is larger than a given fraction of the allocatable biomass,
281       !     the rest is taken from the tissues themselves.
282       !     We suppose that respiration is not dependent on leaf age ->
283       !     do not change age structure.
284       !
285
286       ! maximum part of allocatable biomass used for respiration
287       bm_tax_max(:) = tax_max * bm_alloc_tot(:,j)
288
289       DO i = 1, npts
290
291          IF ( ( bm_alloc_tot(i,j) .GT. zero ) .AND. &
292               ( ( resp_maint(i,j) * dt ) .LT. bm_tax_max(i) ) ) THEN
293
294             bm_alloc_tot(i,j) = bm_alloc_tot(i,j) - resp_maint(i,j) * dt
295
296             !Shilong        ELSEIF ( resp_maint(i,j) .GT. zero ) THEN
297          ELSEIF ( resp_maint(i,j) .GT. min_stomate ) THEN
298
299             ! remaining allocatable biomass
300             bm_alloc_tot(i,j) = bm_alloc_tot(i,j) - bm_tax_max(i)
301
302             ! biomass that remains to be taken away from tissues
303             bm_pump(i) = resp_maint(i,j) * dt - bm_tax_max(i)
304
305             ! take biomass from tissues
306
307             biomass(i,j,ileaf) = biomass(i,j,ileaf) - &
308                  bm_pump(i) * resp_maint_part(i,j,ileaf) / resp_maint(i,j)
309             biomass(i,j,isapabove) = biomass(i,j,isapabove) - &
310                  bm_pump(i) * resp_maint_part(i,j,isapabove) / resp_maint(i,j)
311             biomass(i,j,isapbelow) = biomass(i,j,isapbelow) - &
312                  bm_pump(i) * resp_maint_part(i,j,isapbelow) / resp_maint(i,j)
313             biomass(i,j,iroot) = biomass(i,j,iroot) - &
314                  bm_pump(i) * resp_maint_part(i,j,iroot) / resp_maint(i,j)
315             biomass(i,j,ifruit) = biomass(i,j,ifruit) - &
316                  bm_pump(i) * resp_maint_part(i,j,ifruit) / resp_maint(i,j)
317             biomass(i,j,icarbres) = biomass(i,j,icarbres) - &
318                  bm_pump(i) * resp_maint_part(i,j,icarbres) / resp_maint(i,j)
319
320          ENDIF
321
322       ENDDO   ! Fortran95: WHERE - ELSEWHERE construct
323
324       !
325       ! 3.4 dispatch allocatable biomass
326       !
327
328       DO k = 1, nparts
329          bm_alloc(:,j,k) = f_alloc(:,j,k) * bm_alloc_tot(:,j)
330       ENDDO
331
332       !
333       ! 3.5 growth respiration of a plant part is a given fraction of the
334       !     remaining allocatable biomass.
335       !
336
337       resp_growth_part(:,:) = frac_growthresp * bm_alloc(:,j,:) / dt
338
339       bm_alloc(:,j,:) = ( 1. - frac_growthresp ) * bm_alloc(:,j,:)
340
341       !
342       ! 3.6 Total growth respiration of the plant
343       !     VPP killer:
344       !     resp_growth(:,j) = SUM( resp_growth_part(:,:), DIM=2 )
345       !
346
347       resp_growth(:,j) = zero
348
349       DO k = 1, nparts
350          resp_growth(:,j) = resp_growth(:,j) + resp_growth_part(:,k)
351       ENDDO
352
353    ENDDO
354
355    !
356    ! 4 update the biomass, but save the old leaf mass for later
357    !   "old" leaf mass is leaf mass after maintenance respiration
358    !
359
360    lm_old(:,:) = biomass(:,:,ileaf)
361
362    biomass(:,:,:) = biomass(:,:,:) + bm_alloc(:,:,:)
363
364    !
365    ! 5 biomass can become negative in some rare cases, as the GPP can be negative
366    !   (dark respiration).
367    !   In this case, set biomass to some small value. This creation of matter is taken into
368    !   account by decreasing the autotrophic respiration. In this case, maintenance respiration
369    !   can become negative !!!
370    !
371
372    DO k = 1, nparts
373
374       DO j = 2,nvm
375
376          WHERE ( biomass(:,j,k) .LT. zero )
377
378             bm_create(:,j) = min_stomate - biomass(:,j,k)
379             
380             biomass(:,j,k) = biomass(:,j,k) + bm_create(:,j)
381             
382             resp_maint(:,j) = resp_maint(:,j) - bm_create(:,j) / dt
383
384          ENDWHERE
385
386       ENDDO
387
388    ENDDO
389
390    !
391    ! 6 Calculate the NPP (gC/(m**2 of ground/day)
392    !
393
394    DO j = 2,nvm
395       npp(:,j) = gpp(:,j) - resp_growth(:,j) - resp_maint(:,j)
396    ENDDO
397
398    !
399    ! 7 leaf age
400    !
401
402    !
403    ! 7.1 Decrease leaf age in youngest class if new leaf biomass is higher than old one.
404    !
405
406    DO j = 2,nvm
407       leaf_mass_young(:,j) = leaf_frac(:,j,1) * lm_old(:,j) + bm_alloc(:,j,ileaf)
408    ENDDO
409
410    DO j = 2,nvm
411       WHERE ( ( bm_alloc(:,j,ileaf) .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) ) / &
417               & leaf_mass_young(:,j) )
418         
419       ENDWHERE
420    ENDDO
421
422    !
423    ! 7.2 new age class fractions (fraction in youngest class increases)
424    !
425
426    ! 7.2.1 youngest class: new mass in youngest class divided by total new mass
427
428    DO j = 2,nvm
429       WHERE ( biomass(:,j,ileaf) .GT. min_stomate )
430
431          leaf_frac(:,j,1) = leaf_mass_young(:,j) / biomass(:,j,ileaf)
432
433       ENDWHERE
434    ENDDO
435
436    ! 7.2.2 other classes: old mass in leaf age class divided by new mass
437
438    DO m = 2, nleafages
439
440       DO j = 2,nvm
441          WHERE ( biomass(:,j,ileaf) .GT. min_stomate )
442
443             leaf_frac(:,j,m) = leaf_frac(:,j,m) * lm_old(:,j) / biomass(:,j,ileaf)
444
445          ENDWHERE
446       ENDDO
447
448    ENDDO
449
450    !
451    ! 8 Plant age (years)
452    !
453
454    !
455    ! 8.1 Increase age at every time step
456    !
457
458    WHERE ( PFTpresent(:,:) )
459
460       age(:,:) = age(:,:) + dt/one_year
461
462    ELSEWHERE
463
464       age(:,:) = zero
465
466    ENDWHERE
467
468    !
469    ! 8.2 For grasses, decrease age
470    !     if new biomass is higher than old one.
471    !     For trees, age is treated in "establish" if vegetation is dynamic,
472    !     and in turnover routines if it is static (in this case, only take
473    !     into account the age of the heartwood).
474    !
475
476    DO j = 2,nvm
477
478       IF ( .NOT. tree(j) ) THEN
479
480          ! Only four compartments for grasses
481          ! VPP killer:
482          ! bm_new(:) = SUM( biomass(:,j,:), DIM=2 )
483          ! bm_add(:) = SUM( bm_alloc(:,j,:), DIM=2 )
484
485          bm_new(:) = biomass(:,j,ileaf) + biomass(:,j,isapabove) + &
486               biomass(:,j,iroot) + biomass(:,j,ifruit)
487          bm_add(:) = bm_alloc(:,j,ileaf) + bm_alloc(:,j,isapabove) + &
488               bm_alloc(:,j,iroot) + bm_alloc(:,j,ifruit)
489
490          WHERE ( ( bm_new(:) .GT. zero ) .AND. ( bm_add(:) .GT. zero ) )
491             age(:,j) = age(:,j) * ( bm_new(:) - bm_add(:) ) / bm_new(:)
492          ENDWHERE
493
494       ENDIF
495
496    ENDDO
497
498    !
499    ! 9 history
500    !
501
502    CALL histwrite (hist_id_stomate, 'BM_ALLOC_LEAF', itime, &
503         bm_alloc(:,:,ileaf), npts*nvm, horipft_index)
504    CALL histwrite (hist_id_stomate, 'BM_ALLOC_SAP_AB', itime, &
505         bm_alloc(:,:,isapabove), npts*nvm, horipft_index)
506    CALL histwrite (hist_id_stomate, 'BM_ALLOC_SAP_BE', itime, &
507         bm_alloc(:,:,isapbelow), npts*nvm, horipft_index)
508    CALL histwrite (hist_id_stomate, 'BM_ALLOC_ROOT', itime, &
509         bm_alloc(:,:,iroot), npts*nvm, horipft_index)
510    CALL histwrite (hist_id_stomate, 'BM_ALLOC_FRUIT', itime, &
511         bm_alloc(:,:,ifruit), npts*nvm, horipft_index)
512    CALL histwrite (hist_id_stomate, 'BM_ALLOC_RES', itime, &
513         bm_alloc(:,:,icarbres), npts*nvm, horipft_index)
514
515    IF (bavard.GE.4) WRITE(numout,*) 'Leaving npp'
516
517  END SUBROUTINE npp_calc
518
519END MODULE stomate_npp
Note: See TracBrowser for help on using the repository browser.