source: tags/ORCHIDEE_1_9_5_1/ORCHIDEE/src_stomate/stomate_lpj.f90 @ 4114

Last change on this file since 4114 was 119, checked in by martial.mancip, 13 years ago

PC, CE, MM : corrections for IPSLCM5A_C couple carbon model.
File size: 44.2 KB
RevLine 
[8]1! Stomate: phenology, allocation, etc.
2!
3! authors: A. Botta, P. Friedlingstein, C. Morphopoulos, N. Viovy, et al.
4!
5! bits and pieces put together by G. Krinner
6!
7! version 0.0: August 1998
8!
9! $Header: /home/ssipsl/CVSREP/ORCHIDEE/src_stomate/stomate_lpj.f90,v 1.26 2010/08/05 16:02:37 ssipsl Exp $
10! IPSL (2006)
11!  This software is governed by the CeCILL licence see ORCHIDEE/ORCHIDEE_CeCILL.LIC
12!
13MODULE stomate_lpj
14
15  ! modules used:
16
17  USE ioipsl
18  USE grid
19  USE stomate_constants
20  USE lpj_constraints
21  USE lpj_pftinout
22  USE lpj_kill
23  USE lpj_crown
24  USE lpj_fire
25  USE lpj_gap
26  USE lpj_light
27  USE lpj_establish
28  USE lpj_cover
29  USE stomate_prescribe
30  USE stomate_phenology
31  USE stomate_alloc
32  USE stomate_npp
33  USE stomate_turnover
34  USE stomate_litter
35  USE stomate_soilcarbon
36  USE stomate_vmax
37  USE stomate_assimtemp
38  USE constantes_veg
39  USE stomate_lcchange
40  !  USE Write_Field_p
41
42  IMPLICIT NONE
43
44  ! private & public routines
45
46  PRIVATE
47  PUBLIC StomateLpj,StomateLpj_clear
48
49  ! first call
50  LOGICAL, SAVE                                           :: firstcall = .TRUE.
51
52CONTAINS
53
54  SUBROUTINE StomateLpj_clear
55
56    CALL prescribe_clear
57    CALL phenology_clear
58    CALL npp_calc_clear
59    CALL turn_clear
60    CALL soilcarbon_clear
61    CALL constraints_clear
62    CALL establish_clear
63    CALL fire_clear
64    CALL gap_clear
65    CALL light_clear
66    CALL pftinout_clear
67    CALL alloc_clear
68  END SUBROUTINE StomateLpj_clear
69
[46]70  SUBROUTINE StomateLpj (npts, dt_days, EndOfYear, EndOfMonth, &
[8]71       neighbours, resolution, &
72       clay, herbivores, &
73       tsurf_daily, tsoil_daily, t2m_daily, t2m_min_daily, &
74       litterhum_daily, soilhum_daily, &
75       maxmoiavail_lastyear, minmoiavail_lastyear, &
76       gdd0_lastyear, precip_lastyear, &
77       moiavail_month, moiavail_week, tlong_ref, t2m_month, t2m_week, &
78       tsoil_month, soilhum_month, &
79       gdd_m5_dormance, gdd_midwinter, ncd_dormance, ngd_minus5, &
80       turnover_longterm, gpp_daily, time_lowgpp, &
81       time_hum_min, maxfpc_lastyear, resp_maint_part, &
82       PFTpresent, age, fireindex, firelitter, &
83       leaf_age, leaf_frac, biomass, ind, adapted, regenerate, &
84       senescence, when_growthinit, &
85       litterpart, litter, dead_leaves, carbon, black_carbon, lignin_struc, &
86       veget_max, veget, npp_longterm, lm_lastyearmax, veget_lastlight, &
87       everywhere, need_adjacent, RIP_time, &
88       lai, rprof,npp_daily, turnover_daily, turnover_time,&
89       control_moist, control_temp, soilcarbon_input, &
90       co2_to_bm, co2_fire, resp_hetero, resp_maint, resp_growth, &
91       height, deadleaf_cover, vcmax, vjmax, &
92       t_photo_min, t_photo_opt, t_photo_max,bm_to_litter, &
93       prod10,prod100,flux10, flux100, veget_max_new, &
94       convflux,cflux_prod10,cflux_prod100, harvest_above, lcchange)
95
96    !
97    ! 0 declarations
98    !
99
100    ! 0.1 input
101
102    ! Domain size
103    INTEGER(i_std), INTENT(in)                                           :: npts
104    ! time step of Stomate in days
105    REAL(r_std), INTENT(in)                                        :: dt_days
106    ! indices of the 8 neighbours of each grid point (1=N, 2=NE, 3=E, 4=SE, 5=S, 6=SW, 7=W, 8=NW)
107    INTEGER(i_std), DIMENSION(npts,8), INTENT(in)                 :: neighbours
108    ! resolution at each grid point in m (1=E-W, 2=N-S)
109    REAL(r_std), DIMENSION(npts,2), INTENT(in)                     :: resolution
110    ! clay fraction
111    REAL(r_std), DIMENSION(npts), INTENT(in)                       :: clay
112    ! time constant of probability of a leaf to be eaten by a herbivore (days)
113    REAL(r_std), DIMENSION(npts,nvm), INTENT(in)                  :: herbivores
114    ! daily surface temperatures (K)
115    REAL(r_std), DIMENSION(npts), INTENT(in)                       :: tsurf_daily
116    ! daily soil temperatures (K)
117    REAL(r_std), DIMENSION(npts,nbdl), INTENT(in)                  :: tsoil_daily
118    ! daily 2 meter temperatures (K)
119    REAL(r_std), DIMENSION(npts), INTENT(in)                       :: t2m_daily
120    ! daily minimum 2 meter temperatures (K)
121    REAL(r_std), DIMENSION(npts), INTENT(in)                       :: t2m_min_daily
122    ! daily litter humidity
123    REAL(r_std), DIMENSION(npts), INTENT(in)                       :: litterhum_daily
124    ! daily soil humidity
125    REAL(r_std), DIMENSION(npts,nbdl), INTENT(in)                  :: soilhum_daily
126    ! last year's maximum moisture availability
127    REAL(r_std), DIMENSION(npts,nvm), INTENT(in)                  :: maxmoiavail_lastyear
128    ! last year's minimum moisture availability
129    REAL(r_std), DIMENSION(npts,nvm), INTENT(in)                  :: minmoiavail_lastyear
130    ! last year's GDD0
131    REAL(r_std), DIMENSION(npts), INTENT(in)                       :: gdd0_lastyear
132    ! lastyear's precipitation (mm/year)
133    REAL(r_std), DIMENSION(npts), INTENT(in)                       :: precip_lastyear
134    ! "monthly" moisture availability
135    REAL(r_std), DIMENSION(npts,nvm), INTENT(in)                  :: moiavail_month
136    ! "weekly" moisture availability
137    REAL(r_std), DIMENSION(npts,nvm), INTENT(in)                  :: moiavail_week
138    ! "long term" 2 meter reference temperatures (K)
139    REAL(r_std), DIMENSION(npts), INTENT(in)                       :: tlong_ref
140    ! "monthly" 2-meter temperatures (K)
141    REAL(r_std), DIMENSION(npts), INTENT(in)                       :: t2m_month
142    ! "weekly" 2-meter temperatures (K)
143    REAL(r_std), DIMENSION(npts), INTENT(in)                       :: t2m_week
144    ! "monthly" soil temperatures (K)
145    REAL(r_std), DIMENSION(npts,nbdl), INTENT(in)                  :: tsoil_month
146    ! "monthly" soil humidity
147    REAL(r_std), DIMENSION(npts,nbdl), INTENT(in)                  :: soilhum_month
148    ! growing degree days, threshold -5 deg C (for phenology)
149    REAL(r_std), DIMENSION(npts,nvm), INTENT(in)                   :: gdd_m5_dormance
150    ! growing degree days, since midwinter (for phenology)
151    REAL(r_std), DIMENSION(npts,nvm), INTENT(inout)                   :: gdd_midwinter
152    ! number of chilling days, since leaves were lost (for phenology)
153    REAL(r_std), DIMENSION(npts,nvm), INTENT(in)                  :: ncd_dormance
154    ! number of growing days, threshold -5 deg C (for phenology)
155    REAL(r_std), DIMENSION(npts,nvm), INTENT(in)                  :: ngd_minus5
156    ! "long term" turnover rate (gC/(m**2 of ground)/year)
157    REAL(r_std), DIMENSION(npts,nvm,nparts), INTENT(in)           :: turnover_longterm
158    ! daily gross primary productivity (gC/(m**2 of ground)/day)
159    REAL(r_std), DIMENSION(npts,nvm), INTENT(in)                  :: gpp_daily
160    ! duration of dormance (d)
161    REAL(r_std), DIMENSION(npts,nvm), INTENT(in)                  :: time_lowgpp
162    ! time elapsed since strongest moisture availability (d)
163    REAL(r_std), DIMENSION(npts,nvm), INTENT(in)                  :: time_hum_min
164    ! last year's maximum fpc for each natural PFT, on ground
165    REAL(r_std), DIMENSION(npts,nvm), INTENT(in)                  :: maxfpc_lastyear
166    ! maintenance respiration of different plant parts (gC/day/m**2 of ground)
167    REAL(r_std), DIMENSION(npts,nvm,nparts), INTENT(in)             :: resp_maint_part
168
169    ! 0.2 modified fields
170
171    ! PFT exists
172    LOGICAL, DIMENSION(npts,nvm), INTENT(inout)                  :: PFTpresent
173    ! age (years)
174    REAL(r_std), DIMENSION(npts,nvm), INTENT(inout)               :: age
175    ! Probability of fire
176    REAL(r_std), DIMENSION(npts,nvm), INTENT(inout)          :: fireindex
177    ! Longer term litter above the ground, gC/m**2 of ground
178    REAL(r_std), DIMENSION(npts,nvm), INTENT(inout)          :: firelitter
179    ! leaf age (days)
180    REAL(r_std), DIMENSION(npts,nvm,nleafages), INTENT(inout)     :: leaf_age
181    ! fraction of leaves in leaf age class
182    REAL(r_std), DIMENSION(npts,nvm,nleafages), INTENT(inout)     :: leaf_frac
183    ! biomass (gC/(m**2 of ground))
184    REAL(r_std), DIMENSION(npts,nvm,nparts), INTENT(inout)        :: biomass
185    ! density of individuals (1/(m**2 of ground))
186    REAL(r_std), DIMENSION(npts,nvm), INTENT(inout)               :: ind
187    ! Winter too cold? between 0 and 1
188    REAL(r_std), DIMENSION(npts,nvm), INTENT(inout)               :: adapted
189    ! Winter sufficiently cold? between 0 and 1
190    REAL(r_std), DIMENSION(npts,nvm), INTENT(inout)               :: regenerate
191    ! is the plant senescent? (only for deciduous trees - carbohydrate reserve)
192    LOGICAL, DIMENSION(npts,nvm), INTENT(inout)                  :: senescence
193    ! how many days ago was the beginning of the growing season
194    REAL(r_std), DIMENSION(npts,nvm), INTENT(inout)               :: when_growthinit
195    ! fraction of litter above the ground belonging to different PFTs
196    REAL(r_std), DIMENSION(npts,nvm,nlitt), INTENT(inout)         :: litterpart
197    ! metabolic and structural litter, above and below ground (gC/(m**2 of ground))
198    REAL(r_std), DIMENSION(npts,nlitt,nvm,nlevs), INTENT(inout)  :: litter
199    ! dead leaves on ground, per PFT, metabolic and structural,
200    !   in gC/(m**2 of ground)
201    REAL(r_std), DIMENSION(npts,nvm,nlitt), INTENT(inout)         :: dead_leaves
202    ! carbon pool: active, slow, or passive,(gC/(m**2 of ground))
203    REAL(r_std), DIMENSION(npts,ncarb,nvm), INTENT(inout)    :: carbon
204    ! black carbon on the ground (gC/(m**2 of total ground))
205    REAL(r_std), DIMENSION(npts), INTENT(inout)                    :: black_carbon
206    ! ratio Lignine/Carbon in structural litter, above and below ground,
207    ! (gC/(m**2 of ground))
208    REAL(r_std), DIMENSION(npts,nvm,nlevs), INTENT(inout)    :: lignin_struc
209    ! "maximal" coverage fraction of a PFT (LAI -> infinity) on ground
210    REAL(r_std), DIMENSION(npts,nvm), INTENT(inout)               :: veget_max
211    ! fractional coverage on ground, taking into account LAI (=grid-scale fpc)
212    REAL(r_std), DIMENSION(npts,nvm), INTENT(inout)               :: veget
213    ! "long term" net primary productivity (gC/(m**2 of ground)/year)
214    REAL(r_std), DIMENSION(npts,nvm), INTENT(inout)               :: npp_longterm
215    ! last year's maximum leaf mass, for each PFT (gC/(m**2 of ground))
216    REAL(r_std), DIMENSION(npts,nvm), INTENT(inout)               :: lm_lastyearmax
217    ! vegetation fractions (on ground) after last light competition
218    REAL(r_std), DIMENSION(npts,nvm), INTENT(inout)               :: veget_lastlight
219    ! is the PFT everywhere in the grid box or very localized (after its introduction)
220    REAL(r_std), DIMENSION(npts,nvm), INTENT(inout)               :: everywhere
221    ! in order for this PFT to be introduced, does it have to be present in an
222    !   adjacent grid box?
223    LOGICAL, DIMENSION(npts,nvm), INTENT(inout)                  :: need_adjacent
224    ! How much time ago was the PFT eliminated for the last time (y)
225    REAL(r_std), DIMENSION(npts,nvm), INTENT(inout)               :: RIP_time
226    ! Turnover_time of leaves for grasses (d)
227    REAL(r_std), DIMENSION(npts,nvm), INTENT(inout)          :: turnover_time
228
229    ! 0.3 output
230
231    ! leaf area index OF AN INDIVIDUAL PLANT
232    REAL(r_std), DIMENSION(npts,nvm), INTENT(inout)                 :: lai
233    ! root depth. This will, one day, be a prognostic variable. It will be calculated by
234    ! STOMATE (save in restart file & give to hydrology module!), probably somewhere
235    ! in the allocation routine. For the moment, it is prescribed.
236    REAL(r_std), DIMENSION(npts,nvm), INTENT(inout)                 :: rprof
237    ! net primary productivity (gC/day/(m**2 of ground))
238    REAL(r_std), DIMENSION(npts,nvm), INTENT(out)                 :: npp_daily
239    ! Turnover rates (gC/(m**2 of ground)/day)
240    REAL(r_std), DIMENSION(npts,nvm,nparts), INTENT(out)          :: turnover_daily
241    ! moisture control of heterotrophic respiration
242    REAL(r_std), DIMENSION(npts,nlevs), INTENT(inout)                :: control_moist
243    ! temperature control of heterotrophic respiration, above and below
244    REAL(r_std), DIMENSION(npts,nlevs), INTENT(inout)                :: control_temp
245    ! quantity of carbon going into carbon pools from litter decomposition
246    !   (gC/(m**2 of ground)/day)
247    REAL(r_std), DIMENSION(npts,ncarb,nvm), INTENT(inout)      :: soilcarbon_input
248    ! co2 taken up (gC/(m**2 of total ground)/day)
249    !NV devient 2D
250    REAL(r_std), DIMENSION(npts,nvm), INTENT(out)                      :: co2_to_bm
251    ! carbon emitted into the atmosphere by fire (living and dead biomass)
252    ! (in gC/m**2 of total ground/day)
253    !NV devient 2D
254    REAL(r_std), DIMENSION(npts,nvm), INTENT(out)                      :: co2_fire
255    ! heterotrophic respiration (gC/day/m**2 of total ground)
256    REAL(r_std), DIMENSION(npts,nvm), INTENT(inout)            :: resp_hetero
257    ! maintenance respiration (gC/day/(m**2 of total ground))
258    REAL(r_std), DIMENSION(npts,nvm), INTENT(out)                 :: resp_maint
259    ! growth respiration (gC/day/(m**2 of total ground))
260    REAL(r_std), DIMENSION(npts,nvm), INTENT(out)                 :: resp_growth
261    ! height of vegetation (m)
262    REAL(r_std), DIMENSION(npts,nvm), INTENT(inout)                 :: height
263    ! fraction of soil covered by dead leaves
[46]264    REAL(r_std), DIMENSION(npts), INTENT(inout)                      :: deadleaf_cover
[8]265    ! Maximum rate of carboxylation
266    REAL(r_std), DIMENSION(npts,nvm), INTENT(out)                 :: vcmax
267    ! Maximum rate of RUbp regeneration
268    REAL(r_std), DIMENSION(npts,nvm), INTENT(out)                 :: vjmax
269    ! Minimum temperature for photosynthesis (deg C)
270    REAL(r_std), DIMENSION(npts,nvm), INTENT(out)                 :: t_photo_min
271    ! Optimum temperature for photosynthesis (deg C)
272    REAL(r_std), DIMENSION(npts,nvm), INTENT(out)                 :: t_photo_opt
273    ! Maximum temperature for photosynthesis (deg C)
274    REAL(r_std), DIMENSION(npts,nvm), INTENT(out)                 :: t_photo_max
275    ! conversion of biomass to litter (gC/(m**2 of ground)) / day
276    REAL(r_std), DIMENSION(npts,nvm,nparts), INTENT(out)          :: bm_to_litter
277    !
278    ! new "maximal" coverage fraction of a PFT (LAI -> infinity)
279    REAL(r_std), DIMENSION(npts,nvm),INTENT(inout)                 :: veget_max_new 
280
281    ! products remaining in the 10/100 year-turnover pool after the annual release for each compartment
282    ! (10 or 100 + 1 : input from year of land cover change)
283    REAL(r_std),DIMENSION(npts,0:10), INTENT(inout)                        :: prod10
284    REAL(r_std),DIMENSION(npts,0:100), INTENT(inout)                       :: prod100
285    ! annual release from the 10/100 year-turnover pool compartments
286    REAL(r_std),DIMENSION(npts,10), INTENT(inout)                       :: flux10
287    REAL(r_std),DIMENSION(npts,100), INTENT(inout)                      :: flux100
288    ! release during first year following land cover change
289    REAL(r_std),DIMENSION(npts), INTENT(inout)                          :: convflux
290    ! total annual release from the 10/100 year-turnover pool
291    REAL(r_std),DIMENSION(npts), INTENT(inout)                          :: cflux_prod10, cflux_prod100
292    ! harvest above ground biomass for agriculture
293    REAL(r_std), DIMENSION(npts), INTENT(inout)                       :: harvest_above
294
295    ! land cover change flag
296    LOGICAL, INTENT(in)                                                :: lcchange
297
298    ! Land cover change variables + EndOfYear
299    ! Do update of yearly variables? This variable must be .TRUE. once a year
300    LOGICAL, INTENT(in)                                                :: EndOfYear
[46]301    ! Do update of monthly variables ? This variable must be .TRUE. once a month
302    LOGICAL, INTENT(in)                                                :: EndOfMonth
[8]303
304
305    ! 0.4 local
306
307    ! total conversion of biomass to litter (gC/(m**2)) / day
308    REAL(r_std), DIMENSION(npts,nvm)                                   :: tot_bm_to_litter
309    ! total living biomass (gC/(m**2))
310    REAL(r_std), DIMENSION(npts,nvm)                                   :: tot_live_biomass
311    ! biomass increase, i.e. NPP per plant part
312    REAL(r_std), DIMENSION(npts,nvm,nparts)                            :: bm_alloc
313    ! total turnover rate (gC/(m**2)) / day
314    REAL(r_std), DIMENSION(npts,nvm)                                   :: tot_turnover
315    ! total soil and litter carbon (gC/(m**2))
316    REAL(r_std), DIMENSION(npts,nvm)                                   :: tot_litter_soil_carb
317    ! total litter carbon (gC/(m**2))
318    REAL(r_std), DIMENSION(npts,nvm)                                   :: tot_litter_carb
319    ! total soil carbon (gC/(m**2))
320    REAL(r_std), DIMENSION(npts,nvm)                                   :: tot_soil_carb
321    ! crown area of individuals (m**2)
322    REAL(r_std), DIMENSION(npts,nvm)                               :: cn_ind
323    ! fraction that goes into plant part
324    REAL(r_std), DIMENSION(npts,nvm,nparts)                        :: f_alloc
325    ! space availability for trees
326    REAL(r_std), DIMENSION(npts)                                   :: avail_tree
327    ! space availability for grasses
328    REAL(r_std), DIMENSION(npts)                                   :: avail_grass
329
330    INTEGER                                                       :: j
331
332    ! total products remaining in the pool after the annual release
333    REAL(r_std),DIMENSION(npts)                                   :: prod10_total, prod100_total 
334    ! total flux from conflux and the 10/100 year-turnover pool
335    REAL(r_std),DIMENSION(npts)                                       :: cflux_prod_total
336
337    ! "maximal" coverage fraction of a PFT (LAI -> infinity) on ground
338    REAL(r_std),DIMENSION(npts,nvm)                                :: veget_max_old
339
340    REAL(r_std), DIMENSION(npts)                                   :: vartmp
341
342    REAL(r_std), DIMENSION(npts,nvm)                          :: histvar
343
344    ! =========================================================================
345
346    IF (bavard.GE.3) WRITE(numout,*) 'Entering stomate_lpj'
347
348    !
349    ! 1 Initializations
350    !
351
352    !
353    ! 1.1 set outputs to zero
354    !
355    co2_to_bm(:,:) = zero
356    co2_fire(:,:) = zero
357    npp_daily(:,:) = zero
358    turnover_daily(:,:,:) = zero
359    resp_maint(:,:) = zero
360    resp_growth(:,:) = zero
361    harvest_above(:) = zero
362
363    !
364    ! 1.2  initialize some variables
365    !
366
367    bm_to_litter(:,:,:) = zero
368    cn_ind(:,:) = zero
369    veget_max_old(:,:) = veget_max(:,:)
370
371    !
372    ! 1.3 Prescribe some vegetation characteristics if the vegetation is not dynamic
373    !     IF the DGVM is not activated, the density of individuals and their crown
374    !     areas don't matter, but they should be defined for the case we switch on
375    !     the DGVM afterwards.
376    !     At first call, if the DGVM is not activated, impose a minimum biomass for
377    !       prescribed PFTs and declare them present.
378    !
379
380    CALL prescribe (npts, &
381         veget_max, PFTpresent, everywhere, when_growthinit, &
382         biomass, leaf_frac, ind, cn_ind)
383
384    !
385    ! 2 climatic constraints for PFT presence and regenerativeness
386    !   call this even when DGVM is not activated so that "adapted" and "regenerate"
387    !   are kept up to date for the moment when the DGVM is activated.
388    !
389
390    CALL constraints (npts, dt_days, &
391         t2m_month, t2m_min_daily, when_growthinit, &
392         adapted, regenerate)
393
394    !
395    ! 3 PFTs in and out, based on climate criteria
396    !
397
398    IF ( control%ok_dgvm ) THEN
399
400       !
401       ! 3.1 do introduction/elimination
402       !
403
404       CALL pftinout (npts, dt_days, adapted, regenerate, &
405            neighbours, veget, veget_max, &
406            biomass, ind, age, leaf_frac, npp_longterm, lm_lastyearmax, senescence, &
407            PFTpresent, everywhere, when_growthinit, need_adjacent, RIP_time, &
408            co2_to_bm, &
409            avail_tree, avail_grass)
410
411       !
412       ! 3.2 reset attributes for eliminated PFTs.
413       !     This also kills PFTs that had 0 leafmass during the last year. The message
414       !     "... after pftinout" is misleading in this case.
415       !
416
417       CALL kill (npts, 'pftinout  ', lm_lastyearmax, &
418            ind, PFTpresent, cn_ind, biomass, senescence, RIP_time, &
419            lai, age, leaf_age, leaf_frac, &
420            when_growthinit, everywhere, veget, veget_max, bm_to_litter)
421
422       !
423       ! 3.3 calculate new crown area and maximum vegetation cover
424       !
425
426       CALL crown (npts, PFTpresent, &
427            ind, biomass, &
428            veget_max, cn_ind, height)
429
430    ENDIF
431
432    !
433    ! 4 phenology
434    !
435    CALL histwrite (hist_id_stomate, 'WHEN_GROWTHINIT', itime, when_growthinit, npts*nvm, horipft_index)
436    CALL histwrite (hist_id_stomate, 'TIME_LOWGPP', itime, time_lowgpp, npts*nvm, horipft_index)
437
438    WHERE(PFTpresent)
439       histvar=un
440    ELSEWHERE
441       histvar=zero
442    ENDWHERE
443    CALL histwrite (hist_id_stomate, 'PFTPRESENT', itime, histvar, npts*nvm, horipft_index)
444
445    WHERE(gdd_midwinter.EQ.undef)
446       histvar=val_exp
447    ELSEWHERE
448       histvar=gdd_midwinter
449    ENDWHERE
450    CALL histwrite (hist_id_stomate, 'GDD_MIDWINTER', itime, histvar, npts*nvm, horipft_index)
451
452    WHERE(ncd_dormance.EQ.undef)
453       histvar=val_exp
454    ELSEWHERE
455       histvar=ncd_dormance
456    ENDWHERE
457    CALL histwrite (hist_id_stomate, 'NCD_DORMANCE', itime, histvar, npts*nvm, horipft_index)
458
459    CALL phenology (npts, dt_days, PFTpresent, &
460         veget_max, &
461         tlong_ref, t2m_month, t2m_week, gpp_daily, &
462         maxmoiavail_lastyear, minmoiavail_lastyear, &
463         moiavail_month, moiavail_week, &
464         gdd_m5_dormance, gdd_midwinter, ncd_dormance, ngd_minus5, &
465         senescence, time_lowgpp, time_hum_min, &
466         biomass, leaf_frac, leaf_age, &
467         when_growthinit, co2_to_bm, lai)
468
469    !
470    ! 5 allocation
471    !
472
473    CALL alloc (npts, dt_days, &
474         lai, veget_max, senescence, when_growthinit, &
475         moiavail_week, tsoil_month, soilhum_month, &
476         biomass, age, leaf_age, leaf_frac, rprof, f_alloc)
477
478    !
479    ! 6 maintenance and growth respiration. NPP
480    !
481
482    CALL npp_calc (npts, dt_days, &
483         PFTpresent, &
484         tlong_ref, t2m_daily, tsoil_daily, lai, rprof, &
485         gpp_daily, f_alloc, bm_alloc, resp_maint_part,&
486         biomass, leaf_age, leaf_frac, age, &
487         resp_maint, resp_growth, npp_daily)
488
489    IF ( control%ok_dgvm ) THEN
490
491       ! new provisional crown area and maximum vegetation cover after growth
492
493       CALL crown (npts, PFTpresent, &
494            ind, biomass, &
495            veget_max, cn_ind, height)
496
497    ENDIF
498
499    !
500    ! 7 fire.
501    !
502
503    CALL fire (npts, dt_days, litterpart, &
504         litterhum_daily, t2m_daily, lignin_struc, &
505         fireindex, firelitter, biomass, ind, &
506         litter, dead_leaves, bm_to_litter, black_carbon, &
507         co2_fire)
508
509    IF ( control%ok_dgvm ) THEN
510
511       ! reset attributes for eliminated PFTs
512
513       CALL kill (npts, 'fire      ', lm_lastyearmax, &
514            ind, PFTpresent, cn_ind, biomass, senescence, RIP_time, &
515            lai, age, leaf_age, leaf_frac, &
516            when_growthinit, everywhere, veget, veget_max, bm_to_litter)
517
518    ENDIF
519
520    !
521    ! 8 tree mortality. Does not depend on age, therefore does not change crown area.
522    !
523
524    CALL gap (npts, dt_days, &
525         npp_longterm, turnover_longterm, lm_lastyearmax, &
526         PFTpresent, biomass, ind, bm_to_litter)
527
528    IF ( control%ok_dgvm ) THEN
529
530       ! reset attributes for eliminated PFTs
531
532       CALL kill (npts, 'gap       ', lm_lastyearmax, &
533            ind, PFTpresent, cn_ind, biomass, senescence, RIP_time, &
534            lai, age, leaf_age, leaf_frac, &
535            when_growthinit, everywhere, veget, veget_max, bm_to_litter)
536
537    ENDIF
538
539    !
540    ! 9 calculate vcmax, vjmax and photosynthesis temperatures
541    !
542
543    CALL vmax (npts, dt_days, &
544         leaf_age, leaf_frac, &
545         vcmax, vjmax)
546
547    CALL assim_temp (npts, tlong_ref, t2m_month, &
548         t_photo_min, t_photo_opt, t_photo_max)
549
550    !
551    ! 10 leaf senescence and other turnover processes. New lai
552    !
553
554    CALL turn (npts, dt_days, PFTpresent, &
555         herbivores, &
556         maxmoiavail_lastyear, minmoiavail_lastyear, &
557         moiavail_week,  moiavail_month,tlong_ref, t2m_month, t2m_week, veget_max, &
558         leaf_age, leaf_frac, age, lai, biomass, &
559         turnover_daily, senescence,turnover_time)
560
561    !
562    ! 11 light competition
563    !
564
565    IF ( control%ok_dgvm ) THEN
566
567       !
568       ! 11.1 do light competition
569       !
570
571       CALL light (npts, dt_days, &
572            PFTpresent, cn_ind, lai, maxfpc_lastyear, &
573            ind, biomass, veget_lastlight, bm_to_litter)
574
575       !
576       ! 11.2 reset attributes for eliminated PFTs
577       !
578
579       CALL kill (npts, 'light     ', lm_lastyearmax, &
580            ind, PFTpresent, cn_ind, biomass, senescence, RIP_time, &
581            lai, age, leaf_age, leaf_frac, &
582            when_growthinit, everywhere, veget, veget_max, bm_to_litter)
583
584    ENDIF
585
586    !
587    ! 12 establishment of saplings
588    !
589
590    IF ( control%ok_dgvm ) THEN
591
592       !
593       ! 12.1 do establishment
594       !
595
596       CALL establish (npts, dt_days, PFTpresent, regenerate, &
597            neighbours, resolution, need_adjacent, herbivores, &
598            precip_lastyear, gdd0_lastyear, lm_lastyearmax, &
599            cn_ind, lai, avail_tree, avail_grass, &
600            leaf_age, leaf_frac, &
601            ind, biomass, age, everywhere, co2_to_bm, veget_max)
602
603       !
604       ! 12.2 calculate new crown area (and maximum vegetation cover)
605       !
606
607       CALL crown (npts, PFTpresent, &
608            ind, biomass, &
609            veget_max, cn_ind, height)
610
611    ENDIF
612
613    !
614    ! 13 calculate final LAI and vegetation cover.
615    !
616
617    CALL cover (npts, cn_ind, ind, biomass, &
618         veget_max, veget_max_old, veget, &
619         lai, litter, carbon)
620
621    !
622    ! 14 the whole litter stuff:
623    !    litter update, lignin content, PFT parts, litter decay,
624    !    litter heterotrophic respiration, dead leaf soil cover.
625    !    No vertical discretisation in the soil for litter decay.
626    !
627    ! added by shilong for harvest
628    IF(harvest_agri) THEN
629       CALL harvest(npts, dt_days, veget_max, veget, &
630            bm_to_litter, turnover_daily, &
631            harvest_above)
632    ENDIF
633
634    ! 15.1  Land cover change
635
636    !shilong adde turnover_daily
637    IF(EndOfYear) THEN
638       IF (lcchange) THEN
639          CALL lcchange_main (npts, dt_days, veget_max, veget_max_new, &
640               biomass, ind, age, PFTpresent, senescence, when_growthinit, &
641               everywhere, veget, & 
642               co2_to_bm, bm_to_litter, turnover_daily, bm_sapl, tree, cn_ind,flux10,flux100, &
643!!$               prod10,prod100,prod10_total,prod100_total,&
644!!$               convflux,cflux_prod_total,cflux_prod10,cflux_prod100,leaf_frac,&
645               prod10,prod100,convflux,cflux_prod10,cflux_prod100,leaf_frac,&
646               npp_longterm, lm_lastyearmax, litter, carbon)
647       ENDIF
648    ENDIF
649!MM déplacement pour initialisation correcte des grandeurs cumulées :
650    cflux_prod_total(:) = convflux(:) + cflux_prod10(:) + cflux_prod100(:)
651    prod10_total(:)=SUM(prod10,dim=2)
652    prod100_total(:)=SUM(prod100,dim=2)
653    !
654    ! 16 total heterotrophic respiration
655
656    tot_soil_carb=zero
657    tot_litter_carb=zero
658    DO j=2,nvm
659
660       tot_litter_carb(:,j) = tot_litter_carb(:,j) + (litter(:,istructural,j,iabove) + &
661            &          litter(:,imetabolic,j,iabove) + &
662            &          litter(:,istructural,j,ibelow) + litter(:,imetabolic,j,ibelow))
663
664       tot_soil_carb(:,j) = tot_soil_carb(:,j) + (carbon(:,iactive,j) + &
665            &          carbon(:,islow,j)+  carbon(:,ipassive,j))
666
667    ENDDO
668    tot_litter_soil_carb = tot_litter_carb + tot_soil_carb
669
670    tot_live_biomass = biomass(:,:,ileaf) + biomass(:,:,isapabove) + biomass(:,:,isapbelow) +&
671         &             biomass(:,:,iheartabove) + biomass(:,:,iheartbelow) + &
672         &             biomass(:,:,iroot)+ biomass(:,:,ifruit)+ biomass(:,:,icarbres)
673
674    tot_turnover = turnover_daily(:,:,ileaf) + turnover_daily(:,:,isapabove) + &
675         &         turnover_daily(:,:,isapbelow) + turnover_daily(:,:,iheartabove) + &
676         &         turnover_daily(:,:,iheartbelow) + turnover_daily(:,:,iroot) + &
677         &         turnover_daily(:,:,ifruit) + turnover_daily(:,:,icarbres)
678
679    tot_bm_to_litter = bm_to_litter(:,:,ileaf) + bm_to_litter(:,:,isapabove) +&
680         &             bm_to_litter(:,:,isapbelow) + bm_to_litter(:,:,iheartbelow) +&
681         &             bm_to_litter(:,:,iheartabove) + bm_to_litter(:,:,iroot) + &
682         &             bm_to_litter(:,:,ifruit) + bm_to_litter(:,:,icarbres)
683
684    !
685    ! 17 history
686    !
687
688    ! 2d
689
690    CALL histwrite (hist_id_stomate, 'RESOLUTION_X', itime, &
691         resolution(:,1), npts, hori_index)
692    CALL histwrite (hist_id_stomate, 'RESOLUTION_Y', itime, &
693         resolution(:,2), npts, hori_index)
694    CALL histwrite (hist_id_stomate, 'CONTFRAC', itime, &
695         contfrac(:), npts, hori_index)
696
697    CALL histwrite (hist_id_stomate, 'LITTER_STR_AB', itime, &
698         litter(:,istructural,:,iabove), npts*nvm, horipft_index)
699    CALL histwrite (hist_id_stomate, 'LITTER_MET_AB', itime, &
700         litter(:,imetabolic,:,iabove), npts*nvm, horipft_index)
701    CALL histwrite (hist_id_stomate, 'LITTER_STR_BE', itime, &
702         litter(:,istructural,:,ibelow), npts*nvm, horipft_index)
703    CALL histwrite (hist_id_stomate, 'LITTER_MET_BE', itime, &
704         litter(:,imetabolic,:,ibelow), npts*nvm, horipft_index)
705
706    CALL histwrite (hist_id_stomate, 'DEADLEAF_COVER', itime, &
707         deadleaf_cover, npts, hori_index)
708
709    CALL histwrite (hist_id_stomate, 'TOTAL_SOIL_CARB', itime, &
710         tot_litter_soil_carb, npts*nvm, horipft_index)
711    CALL histwrite (hist_id_stomate, 'CARBON_ACTIVE', itime, &
712         carbon(:,iactive,:), npts*nvm, horipft_index)
713    CALL histwrite (hist_id_stomate, 'CARBON_SLOW', itime, &
714         carbon(:,islow,:), npts*nvm, horipft_index)
715    CALL histwrite (hist_id_stomate, 'CARBON_PASSIVE', itime, &
716         carbon(:,ipassive,:), npts*nvm, horipft_index)
717
718    CALL histwrite (hist_id_stomate, 'T2M_MONTH', itime, &
719         t2m_month, npts, hori_index)
720    CALL histwrite (hist_id_stomate, 'T2M_WEEK', itime, &
721         t2m_week, npts, hori_index)
722
723    CALL histwrite (hist_id_stomate, 'HET_RESP', itime, &
724         resp_hetero(:,:), npts*nvm, horipft_index)
725
726    CALL histwrite (hist_id_stomate, 'BLACK_CARBON', itime, &
727         black_carbon, npts, hori_index)
728
729    CALL histwrite (hist_id_stomate, 'FIREINDEX', itime, &
730         fireindex(:,:), npts*nvm, horipft_index)
731    CALL histwrite (hist_id_stomate, 'LITTERHUM', itime, &
732         litterhum_daily, npts, hori_index)
733    CALL histwrite (hist_id_stomate, 'CO2_FIRE', itime, &
734         co2_fire, npts*nvm, horipft_index)
735    CALL histwrite (hist_id_stomate, 'CO2_TAKEN', itime, &
736         co2_to_bm, npts*nvm, horipft_index)
737    ! land cover change
738    CALL histwrite (hist_id_stomate, 'CONVFLUX', itime, &
739         convflux, npts, hori_index)
740    CALL histwrite (hist_id_stomate, 'CFLUX_PROD10', itime, &
741         cflux_prod10, npts, hori_index)
742    CALL histwrite (hist_id_stomate, 'CFLUX_PROD100', itime, &
743         cflux_prod100, npts, hori_index)
744    CALL histwrite (hist_id_stomate, 'HARVEST_ABOVE', itime, &
745         harvest_above, npts, hori_index)
746
747    ! 3d
748
749    CALL histwrite (hist_id_stomate, 'LAI', itime, &
750         lai, npts*nvm, horipft_index)
751    CALL histwrite (hist_id_stomate, 'VEGET', itime, &
752         veget, npts*nvm, horipft_index)
753    CALL histwrite (hist_id_stomate, 'VEGET_MAX', itime, &
754         veget_max, npts*nvm, horipft_index)
755    CALL histwrite (hist_id_stomate, 'NPP', itime, &
756         npp_daily, npts*nvm, horipft_index)
757    CALL histwrite (hist_id_stomate, 'GPP', itime, &
758         gpp_daily, npts*nvm, horipft_index)
759    CALL histwrite (hist_id_stomate, 'IND', itime, &
760         ind, npts*nvm, horipft_index)
761    CALL histwrite (hist_id_stomate, 'TOTAL_M', itime, &
762         tot_live_biomass, npts*nvm, horipft_index)
763    CALL histwrite (hist_id_stomate, 'LEAF_M', itime, &
764         biomass(:,:,ileaf), npts*nvm, horipft_index)
765    CALL histwrite (hist_id_stomate, 'SAP_M_AB', itime, &
766         biomass(:,:,isapabove), npts*nvm, horipft_index)
767    CALL histwrite (hist_id_stomate, 'SAP_M_BE', itime, &
768         biomass(:,:,isapbelow), npts*nvm, horipft_index)
769    CALL histwrite (hist_id_stomate, 'HEART_M_AB', itime, &
770         biomass(:,:,iheartabove), npts*nvm, horipft_index)
771    CALL histwrite (hist_id_stomate, 'HEART_M_BE', itime, &
772         biomass(:,:,iheartbelow), npts*nvm, horipft_index)
773    CALL histwrite (hist_id_stomate, 'ROOT_M', itime, &
774         biomass(:,:,iroot), npts*nvm, horipft_index)
775    CALL histwrite (hist_id_stomate, 'FRUIT_M', itime, &
776         biomass(:,:,ifruit), npts*nvm, horipft_index)
777    CALL histwrite (hist_id_stomate, 'RESERVE_M', itime, &
778         biomass(:,:,icarbres), npts*nvm, horipft_index)
779    CALL histwrite (hist_id_stomate, 'TOTAL_TURN', itime, &
780         tot_turnover, npts*nvm, horipft_index)
781    CALL histwrite (hist_id_stomate, 'LEAF_TURN', itime, &
782         turnover_daily(:,:,ileaf), npts*nvm, horipft_index)
783    CALL histwrite (hist_id_stomate, 'SAP_AB_TURN', itime, &
784         turnover_daily(:,:,isapabove), npts*nvm, horipft_index)
785    CALL histwrite (hist_id_stomate, 'ROOT_TURN', itime, &
786         turnover_daily(:,:,iroot), npts*nvm, horipft_index)
787    CALL histwrite (hist_id_stomate, 'FRUIT_TURN', itime, &
788         turnover_daily(:,:,ifruit), npts*nvm, horipft_index)
789    CALL histwrite (hist_id_stomate, 'TOTAL_BM_LITTER', itime, &
790         tot_bm_to_litter, npts*nvm, horipft_index)
791    CALL histwrite (hist_id_stomate, 'LEAF_BM_LITTER', itime, &
792         bm_to_litter(:,:,ileaf), npts*nvm, horipft_index)
793    CALL histwrite (hist_id_stomate, 'SAP_AB_BM_LITTER', itime, &
794         bm_to_litter(:,:,isapabove), npts*nvm, horipft_index)
795    CALL histwrite (hist_id_stomate, 'SAP_BE_BM_LITTER', itime, &
796         bm_to_litter(:,:,isapbelow), npts*nvm, horipft_index)
797    CALL histwrite (hist_id_stomate, 'HEART_AB_BM_LITTER', itime, &
798         bm_to_litter(:,:,iheartabove), npts*nvm, horipft_index)
799    CALL histwrite (hist_id_stomate, 'HEART_BE_BM_LITTER', itime, &
800         bm_to_litter(:,:,iheartbelow), npts*nvm, horipft_index)
801    CALL histwrite (hist_id_stomate, 'ROOT_BM_LITTER', itime, &
802         bm_to_litter(:,:,iroot), npts*nvm, horipft_index)
803    CALL histwrite (hist_id_stomate, 'FRUIT_BM_LITTER', itime, &
804         bm_to_litter(:,:,ifruit), npts*nvm, horipft_index)
805    CALL histwrite (hist_id_stomate, 'RESERVE_BM_LITTER', itime, &
806         bm_to_litter(:,:,icarbres), npts*nvm, horipft_index)
807    CALL histwrite (hist_id_stomate, 'MAINT_RESP', itime, &
808         resp_maint, npts*nvm, horipft_index)
809    CALL histwrite (hist_id_stomate, 'GROWTH_RESP', itime, &
810         resp_growth, npts*nvm, horipft_index)
811    CALL histwrite (hist_id_stomate, 'AGE', itime, &
812         age, npts*nvm, horipft_index)
813    CALL histwrite (hist_id_stomate, 'HEIGHT', itime, &
814         height, npts*nvm, horipft_index)
815    CALL histwrite (hist_id_stomate, 'MOISTRESS', itime, &
816         moiavail_week, npts*nvm, horipft_index)
817    CALL histwrite (hist_id_stomate, 'VCMAX', itime, &
818         vcmax, npts*nvm, horipft_index)
819    CALL histwrite (hist_id_stomate, 'TURNOVER_TIME', itime, &
820         turnover_time, npts*nvm, horipft_index)
821    ! land cover change
822    CALL histwrite (hist_id_stomate, 'PROD10', itime, &
823         prod10, npts*11, horip11_index)
824    CALL histwrite (hist_id_stomate, 'PROD100', itime, &
825         prod100, npts*101, horip101_index)
826    CALL histwrite (hist_id_stomate, 'FLUX10', itime, &
827         flux10, npts*10, horip10_index)
828    CALL histwrite (hist_id_stomate, 'FLUX100', itime, &
829         flux100, npts*100, horip100_index)
830
831    IF ( hist_id_stomate_IPCC > 0 ) THEN
832       vartmp(:)=SUM(tot_live_biomass*veget_max,dim=2)/1e3*contfrac
833       CALL histwrite (hist_id_stomate_IPCC, "cVeg", itime, &
834         vartmp, npts, hori_index)
835       vartmp(:)=SUM(tot_litter_carb*veget_max,dim=2)/1e3*contfrac
836       CALL histwrite (hist_id_stomate_IPCC, "cLitter", itime, &
837         vartmp, npts, hori_index)
838       vartmp(:)=SUM(tot_soil_carb*veget_max,dim=2)/1e3*contfrac
839       CALL histwrite (hist_id_stomate_IPCC, "cSoil", itime, &
840         vartmp, npts, hori_index)
841       vartmp(:)=(prod10_total + prod100_total)/1e3
842       CALL histwrite (hist_id_stomate_IPCC, "cProduct", itime, &
843         vartmp, npts, hori_index)
844       vartmp(:)=SUM(lai*veget_max,dim=2)*contfrac
845       CALL histwrite (hist_id_stomate_IPCC, "lai", itime, &
846         vartmp, npts, hori_index)
847       vartmp(:)=SUM(gpp_daily*veget_max,dim=2)/1e3/one_day*contfrac
848       CALL histwrite (hist_id_stomate_IPCC, "gpp", itime, &
849         vartmp, npts, hori_index)
850       vartmp(:)=SUM((resp_maint+resp_growth)*veget_max,dim=2)/1e3/one_day*contfrac
851       CALL histwrite (hist_id_stomate_IPCC, "ra", itime, &
852         vartmp, npts, hori_index)
853       vartmp(:)=SUM(npp_daily*veget_max,dim=2)/1e3/one_day*contfrac
854       CALL histwrite (hist_id_stomate_IPCC, "npp", itime, &
855         vartmp, npts, hori_index)
856       vartmp(:)=SUM(resp_hetero*veget_max,dim=2)/1e3/one_day*contfrac
857       CALL histwrite (hist_id_stomate_IPCC, "rh", itime, &
858         vartmp, npts, hori_index)
859       vartmp(:)=SUM(co2_fire*veget_max,dim=2)/1e3/one_day*contfrac
860       CALL histwrite (hist_id_stomate_IPCC, "fFire", itime, &
861         vartmp, npts, hori_index)
862       vartmp(:)=harvest_above/1e3/one_day*contfrac
863       CALL histwrite (hist_id_stomate_IPCC, "fHarvest", itime, &
864         vartmp, npts, hori_index)
865       vartmp(:)=cflux_prod_total/1e3/one_day*contfrac
866       CALL histwrite (hist_id_stomate_IPCC, "fLuc", itime, &
867         vartmp, npts, hori_index)
868       vartmp(:)=(SUM((gpp_daily-(resp_maint+resp_growth+resp_hetero)-co2_fire) &
869            &        *veget_max,dim=2)-cflux_prod_total-harvest_above)/1e3/one_day*contfrac
870       CALL histwrite (hist_id_stomate_IPCC, "nbp", itime, &
871         vartmp, npts, hori_index)
872       vartmp(:)=SUM(tot_bm_to_litter*veget_max,dim=2)/1e3/one_day*contfrac
873       CALL histwrite (hist_id_stomate_IPCC, "fVegLitter", itime, &
874         vartmp, npts, hori_index)
875       vartmp(:)=SUM(SUM(soilcarbon_input,dim=2)*veget_max,dim=2)/1e3/one_day*contfrac
876       CALL histwrite (hist_id_stomate_IPCC, "fLitterSoil", itime, &
877         vartmp, npts, hori_index)
878       vartmp(:)=SUM(biomass(:,:,ileaf)*veget_max,dim=2)/1e3*contfrac
879       CALL histwrite (hist_id_stomate_IPCC, "cLeaf", itime, &
880         vartmp, npts, hori_index)
881       vartmp(:)=SUM((biomass(:,:,isapabove)+biomass(:,:,iheartabove))*veget_max,dim=2)/1e3*contfrac
882       CALL histwrite (hist_id_stomate_IPCC, "cWood", itime, &
883         vartmp, npts, hori_index)
884       vartmp(:)=SUM(( biomass(:,:,iroot) + biomass(:,:,isapbelow) + biomass(:,:,iheartbelow) ) &
885            &        *veget_max,dim=2)/1e3*contfrac
886       CALL histwrite (hist_id_stomate_IPCC, "cRoot", itime, &
887         vartmp, npts, hori_index)
888       vartmp(:)=SUM(( biomass(:,:,icarbres) + biomass(:,:,ifruit))*veget_max,dim=2)/1e3*contfrac
889       CALL histwrite (hist_id_stomate_IPCC, "cMisc", itime, &
890         vartmp, npts, hori_index)
891       vartmp(:)=SUM((litter(:,istructural,:,iabove)+litter(:,imetabolic,:,iabove))*veget_max,dim=2)/1e3*contfrac
892       CALL histwrite (hist_id_stomate_IPCC, "cLitterAbove", itime, &
893         vartmp, npts, hori_index)
894       vartmp(:)=SUM((litter(:,istructural,:,ibelow)+litter(:,imetabolic,:,ibelow))*veget_max,dim=2)/1e3*contfrac
895       CALL histwrite (hist_id_stomate_IPCC, "cLitterBelow", itime, &
896         vartmp, npts, hori_index)
897       vartmp(:)=SUM(carbon(:,iactive,:)*veget_max,dim=2)/1e3*contfrac
898       CALL histwrite (hist_id_stomate_IPCC, "cSoilFast", itime, &
899         vartmp, npts, hori_index)
900       vartmp(:)=SUM(carbon(:,islow,:)*veget_max,dim=2)/1e3*contfrac
901       CALL histwrite (hist_id_stomate_IPCC, "cSoilMedium", itime, &
902         vartmp, npts, hori_index)
903       vartmp(:)=SUM(carbon(:,ipassive,:)*veget_max,dim=2)/1e3*contfrac
904       CALL histwrite (hist_id_stomate_IPCC, "cSoilSlow", itime, &
905         vartmp, npts, hori_index)
906       DO j=1,nvm
907          histvar(:,j)=veget_max(:,j)*contfrac(:)*100
908       ENDDO
909       CALL histwrite (hist_id_stomate_IPCC, "landCoverFrac", itime, &
910         histvar, npts*nvm, horipft_index)
911       vartmp(:)=(veget_max(:,3)+veget_max(:,6)+veget_max(:,8)+veget_max(:,9))*contfrac*100
912       CALL histwrite (hist_id_stomate_IPCC, "treeFracPrimDec", itime, &
913          vartmp, npts, hori_index)
914       vartmp(:)=(veget_max(:,2)+veget_max(:,4)+veget_max(:,5)+veget_max(:,7))*contfrac*100
915       CALL histwrite (hist_id_stomate_IPCC, "treeFracPrimEver", itime, &
916         vartmp, npts, hori_index)
917       vartmp(:)=(veget_max(:,10)+veget_max(:,12))*contfrac*100
918       CALL histwrite (hist_id_stomate_IPCC, "c3PftFrac", itime, &
919         vartmp, npts, hori_index)
920       vartmp(:)=(veget_max(:,11)+veget_max(:,13))*contfrac*100
921       CALL histwrite (hist_id_stomate_IPCC, "c4PftFrac", itime, &
922         vartmp, npts, hori_index)
923       vartmp(:)=SUM(resp_growth*veget_max,dim=2)/1e3/one_day*contfrac
924       CALL histwrite (hist_id_stomate_IPCC, "rGrowth", itime, &
925         vartmp, npts, hori_index)
926       vartmp(:)=SUM(resp_maint*veget_max,dim=2)/1e3/one_day*contfrac
927       CALL histwrite (hist_id_stomate_IPCC, "rMaint", itime, &
928         vartmp, npts, hori_index)
929       vartmp(:)=SUM(bm_alloc(:,:,ileaf)*veget_max,dim=2)/1e3/one_day*contfrac
930       CALL histwrite (hist_id_stomate_IPCC, "nppLeaf", itime, &
931         vartmp, npts, hori_index)
932       vartmp(:)=SUM(bm_alloc(:,:,isapabove)*veget_max,dim=2)/1e3/one_day*contfrac
933       CALL histwrite (hist_id_stomate_IPCC, "nppWood", itime, &
934         vartmp, npts, hori_index)
935       vartmp(:)=SUM(( bm_alloc(:,:,isapbelow) + bm_alloc(:,:,iroot) )*veget_max,dim=2)/1e3/one_day*contfrac
936       CALL histwrite (hist_id_stomate_IPCC, "nppRoot", itime, &
937         vartmp, npts, hori_index)
938
939       CALL histwrite (hist_id_stomate_IPCC, 'RESOLUTION_X', itime, &
940            resolution(:,1), npts, hori_index)
941       CALL histwrite (hist_id_stomate_IPCC, 'RESOLUTION_Y', itime, &
942            resolution(:,2), npts, hori_index)
943       CALL histwrite (hist_id_stomate_IPCC, 'CONTFRAC', itime, &
944            contfrac(:), npts, hori_index)
945
946    ENDIF
947
948    IF (bavard.GE.4) WRITE(numout,*) 'Leaving stomate_lpj'
949
950  END SUBROUTINE StomateLpj
951
952  SUBROUTINE harvest(npts, dt_days, veget_max, veget, &
953       bm_to_litter, turnover_daily, &
954       harvest_above)
955    ! 0.1 input
956
957    ! Domain size
958    INTEGER, INTENT(in)                                            :: npts
959
960    ! Time step (days)
961    REAL(r_std), INTENT(in)                                         :: dt_days
962
963    ! new "maximal" coverage fraction of a PFT (LAI -> infinity) on ground
964    REAL(r_std), DIMENSION(npts,nvm), INTENT(in)                 :: veget_max
965    ! 0.2 modified fields
966
967    ! fractional coverage on natural/agricultural ground, taking into
968    !   account LAI (=grid-scale fpc)
969    REAL(r_std), DIMENSION(npts,nvm), INTENT(inout)                    :: veget
970
971    ! conversion of biomass to litter (gC/(m**2 of nat/agri ground)) / day
972    REAL(r_std), DIMENSION(npts,nvm,nparts), INTENT(inout)           :: bm_to_litter
973
974    ! Turnover rates (gC/(m**2 of ground)/day)
975    REAL(r_std), DIMENSION(npts,nvm,nparts), INTENT(inout)          :: turnover_daily
976    ! harvest above ground biomass for agriculture
977    REAL(r_std), DIMENSION(npts), INTENT(inout)                       :: harvest_above
978    ! 0.4 local
979
980    ! indices
981    INTEGER(i_std)                                                     :: i, j, k, l, m
982
983    ! biomass increase (gC/(m**2 of ground))
984    REAL(r_std)                                                        :: above_old
985    ! yearly initialisation
986    above_old             = zero
987    harvest_above         = zero
988    DO i = 1, npts
989       DO j=1, nvm
990          IF (j > 11) THEN
991             above_old = turnover_daily(i,j,ileaf) + turnover_daily(i,j,isapabove) + &
992                  &       turnover_daily(i,j,iheartabove) + turnover_daily(i,j,ifruit) + &
993                  &       turnover_daily(i,j,icarbres) + turnover_daily(i,j,isapbelow) + &
994                  &       turnover_daily(i,j,iheartbelow) + turnover_daily(i,j,iroot)
995
996             turnover_daily(i,j,ileaf) = turnover_daily(i,j,ileaf)*0.55
997             turnover_daily(i,j,isapabove) = turnover_daily(i,j,isapabove)*0.55
998             turnover_daily(i,j,isapbelow) = turnover_daily(i,j,isapbelow)*0.55
999             turnover_daily(i,j,iheartabove) = turnover_daily(i,j,iheartabove)*0.55
1000             turnover_daily(i,j,iheartbelow) = turnover_daily(i,j,iheartbelow)*0.55
1001             turnover_daily(i,j,iroot) = turnover_daily(i,j,iroot)*0.55
1002             turnover_daily(i,j,ifruit) = turnover_daily(i,j,ifruit)*0.55
1003             turnover_daily(i,j,icarbres) = turnover_daily(i,j,icarbres)*0.55
1004             harvest_above(i)  = harvest_above(i) + veget_max(i,j) * above_old *0.45
1005
1006          ENDIF
1007       ENDDO
1008    ENDDO
1009
1010!!$    harvest_above = harvest_above
1011  END SUBROUTINE harvest
1012END MODULE stomate_lpj
Note: See TracBrowser for help on using the repository browser.