Context Navigation

source: branches/ORCHIDEE_EXT/ORCHIDEE/src_stomate/stomate_lpj.f90 @ 257

Last change on this file since 257 was 257, checked in by didier.solyga, 13 years ago
Externalized version merged with the trunk
File size: 48.6 KB

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

Note: See TracBrowser for help on using the repository browser.

Download in other formats: