Context Navigation

source: tags/ORCHIDEE_4_1/ORCHIDEE/src_stomate/lpj_establish.f90 @ 7761

Last change on this file since 7761 was 7273, checked in by sebastiaan.luyssaert, 3 years ago
Contributes to tickets #633 and #653. Adding a period for unmanaged forest to decay after reaching the min density. Made crown dimensions PFT-dependent
Property svn:keywords set to `HeadURL Date Author Revision`
File size: 46.3 KB

Line
1	! =================================================================================================================================
2	! MODULE : lpj_establish
3	!
4	! CONTACT : orchidee-help _at_ listes.ipsl.fr
5	!
6	! LICENCE : IPSL (2006)
7	! This software is governed by the CeCILL licence see ORCHIDEE/ORCHIDEE_CeCILL.LIC
8	!
9	!>\BRIEF Establish pft's
10	!!
11	!!\n DESCRIPTION: None
12	!!
13	!! RECENT CHANGE(S): None
14	!!
15	!! REFERENCE(S) :
16	!! - Sitch, S., B. Smith, et al. (2003), Evaluation of ecosystem dynamics,
17	!! plant geography and terrestrial carbon cycling in the LPJ dynamic
18	!! global vegetation model, Global Change Biology, 9, 161-185.\n
19	!! - Haxeltine, A. and I. C. Prentice (1996), BIOME3: An equilibrium
20	!! terrestrial biosphere model based on ecophysiological constraints,
21	!! resource availability, and competition among plant functional types,
22	!! Global Biogeochemical Cycles, 10(4), 693-709.\n
23	!! - Smith, B., I. C. Prentice, et al. (2001), Representation of vegetation
24	!! dynamics in the modelling of terrestrial ecosystems: comparing two
25	!! contrasting approaches within European climate space,
26	!! Global Ecology and Biogeography, 10, 621-637.\n
27	!!
28	!! SVN :
29	!! $HeadURL$
30	!! $Date$
31	!! $Revision$
32	!! \n
33	!_ ================================================================================================================================
34
35	MODULE lpj_establish
36
37	! modules used:
38	USE xios_orchidee
39	USE ioipsl_para
40	USE stomate_data
41	USE constantes
42	USE grid
43	USE function_library, ONLY: biomass_to_lai
44
45	IMPLICIT NONE
46
47	! private & public routines
48	PRIVATE
49	PUBLIC establish,establish_clear
50
51	LOGICAL, SAVE :: firstcall_establish = .TRUE. !! first call
52	!$OMP THREADPRIVATE(firstcall_establish)
53	CONTAINS
54
55
56	!! ================================================================================================================================
57	!! SUBROUTINE : fire_clear
58	!!
59	!>\BRIEF Set the firstcall_establish flag to .TRUE. and activate initialization
60	!_ ================================================================================================================================
61
62	SUBROUTINE establish_clear
63	firstcall_establish = .TRUE.
64	END SUBROUTINE establish_clear
65
66
67	! =================================================================================================================================
68	! SUBROUTINE : establish
69	!
70	!>\BRIEF Calculate sstablishment of new woody PFT and herbaceous PFTs
71	!!
72	!! DESCRIPTION : Establishments of new woody and herbaceous PFT are simulated.
73	!! Maximum establishment rate (0.12) declines due to competition for light (space).
74	!! There are two establishment estimates: one for the for DGVM and one for the
75	!! static cases.\n
76	!! In the case of DGVM, competitive process of establishment for the area of
77	!! available space is represented using more detailed description compared with static
78	!! one. Biomass and distribution of plant age are updated on the basis of changes
79	!! in number of individuals. Finally excess sapwood of is converted to heartwood.
80	!!
81	!! \latexonly
82	!! \input{equation_lpj_establish.tex}
83	!! \endlatexonly
84	!! \n
85	!!
86	!! RECENT CHANGE(S): None
87	!!
88	!! REFERENCE(S) :
89	!! Smith, B., I. C. Prentice, et al. (2001), Representation of vegetation
90	!! dynamics in the modelling of terrestrial ecosystems: comparing two
91	!! contrasting approaches within European climate space,
92	!! Global Ecology and Biogeography, 10, 621-637.
93	!!
94	!! FLOWCHART :
95	!! \latexonly
96	!! \includegraphics[scale = 0.7]{establish.png}
97	!! \endlatexonly
98	!! \n
99	!_ ================================================================================================================================
100
101	SUBROUTINE establish (npts, dt, PFTpresent, regenerate, &
102	neighbours, resolution, need_adjacent, herbivores, &
103	precip_annual, gdd0, lm_lastyearmax, &
104	cn_ind, lai, avail_tree, avail_grass, npp_longterm, &
105	leaf_age, leaf_frac, &
106	ind, biomass, age, everywhere, co2_to_bm, &
107	soil_n_min, n_uptake, nstress_season, veget_cov_max, woodmass_ind, &
108	mortality, bm_sapl_2D, bm_to_litter)
109
110	!! 0. Variable and parameter declaration
111
112	!! 0.1 Input variables
113
114	INTEGER(i_std), INTENT(in) :: npts !! Domain size - number of pixels (dimensionless)
115	REAL(r_std), INTENT(in) :: dt !! Time step of vegetation dynamics for stomate
116	!! (days)
117	LOGICAL, DIMENSION(npts,nvm), INTENT(in) :: PFTpresent !! Is pft there (unitless)
118	REAL(r_std), DIMENSION(npts,nvm), INTENT(in) :: regenerate !! Winter sufficiently cold (unitless)
119	INTEGER(i_std), DIMENSION(npts,NbNeighb), INTENT(in) :: neighbours !! indices of the neighbours of each grid point
120	!! (unitless);
121	!! [1=North and then clockwise]
122	REAL(r_std), DIMENSION(npts,2), INTENT(in) :: resolution !! resolution at each grid point (m); 1=E-W, 2=N-S
123	LOGICAL, DIMENSION(npts,nvm), INTENT(in) :: need_adjacent !! in order for this PFT to be introduced, does it
124	!! have to be present in an adjacent grid box?
125	REAL(r_std), DIMENSION(npts,nvm), INTENT(in) :: herbivores !! time constant of probability of a leaf to
126	!! be eaten by a herbivore (days)
127	REAL(r_std), DIMENSION(npts), INTENT(in) :: precip_annual !! annual precipitation (mm year^{-1})
128	REAL(r_std), DIMENSION(npts), INTENT(in) :: gdd0 !! growing degree days (degree C)
129	REAL(r_std), DIMENSION(npts,nvm), INTENT(in) :: lm_lastyearmax !! last year's maximum leaf mass for each PFT
130	!! (gC m^{-2 })
131	REAL(r_std), DIMENSION(npts,nvm), INTENT(in) :: cn_ind !! crown area of individuals (m^2)
132	REAL(r_std), DIMENSION(npts,nvm), INTENT(in) :: lai !! leaf area index OF an individual plant
133	!! (m^2 m^{-2})
134	REAL(r_std), DIMENSION(npts), INTENT(in) :: avail_tree !! space availability for trees (unitless)
135	REAL(r_std), DIMENSION(npts), INTENT(in) :: avail_grass !! space availability for grasses (unitless)
136	REAL(r_std), DIMENSION(npts,nvm), INTENT(in) :: npp_longterm !! longterm NPP, for each PFT (gC m^{-2})
137	REAL(r_std), DIMENSION(npts,nvm), INTENT(in) :: veget_cov_max !! "maximal" coverage fraction of a PFT
138	!! (LAI -> infinity) on ground (unitless)
139	REAL(r_std), DIMENSION(npts,nvm),INTENT(in) :: mortality !! Fraction of individual dying this time
140	!! step (0 to 1, unitless)
141	REAL(r_std), DIMENSION(npts,nvm,nparts,nelements),INTENT(in) :: bm_sapl_2D
142
143	!! 0.2 Output variables
144
145	!! 0.3 Modified variables
146
147	REAL(r_std), DIMENSION(npts,nvm,nleafages), INTENT(inout) :: leaf_age !! leaf age (days)
148	REAL(r_std), DIMENSION(npts,nvm,nleafages), INTENT(inout) :: leaf_frac !! fraction of leaves in leaf age class (unitless)
149	REAL(r_std), DIMENSION(npts,nvm), INTENT(inout) :: ind !! Number of individuals (individuals m^{-2})
150	REAL(r_std), DIMENSION(npts,nvm,nparts,nelements), INTENT(inout):: biomass !! biomass (gC m^{-2 })
151	REAL(r_std), DIMENSION(npts,nvm), INTENT(inout) :: age !! mean age (years)
152	REAL(r_std), DIMENSION(npts,nvm), INTENT(inout) :: everywhere !! is the PFT everywhere in the grid box or very
153	!! localized (unitless)
154	REAL(r_std), DIMENSION(npts,nvm), INTENT(inout) :: co2_to_bm !! biomass up take for establishment i.e.
155	!! pseudo-photosynthesis (gC m^{-2} day^{-1})
156	REAL(r_std), DIMENSION(npts,nvm,nnspec), INTENT(inout) :: soil_n_min !! Nitrogen in mineral soil,
157	!! (gN/(m**2 of PFT))
158	!! to correct for nitrogen creation associated with establishmnet
159	REAL(r_std), DIMENSION(npts,nvm,nionspec), INTENT(inout) :: n_uptake !! Plant uptake of ammonium and nitrate
160	!! (gN/m**2/day)
161	REAL(r_std), DIMENSION(npts,nvm), INTENT(in) :: nstress_season !! seasonal nitrogen stress
162	!! (unitless)
163	!! (get's initialised at first establishment)
164	REAL(r_std), DIMENSION(npts,nvm), INTENT(inout) :: woodmass_ind !! woodmass of the individual, needed to calculate
165	!! crownarea in lpj_crownarea (gC m^{-2 })
166	REAL(r_std), DIMENSION(npts,nvm,nparts,nelements), INTENT(inout) :: bm_to_litter !!Biomass transfer to litter
167
168	!! 0.4 Local variables
169
170	REAL(r_std) :: tau_eatup !! time during which a sapling can be entirely
171	!! eaten by herbivores (days)
172	REAL(r_std), DIMENSION(npts,nvm) :: fpc_nat !! new fpc, foliage projective cover: fractional
173	!! coverage (unitless)
174	REAL(r_std), DIMENSION(npts) :: estab_rate_max_climate_tree !! maximum tree establishment rate,
175	!! based on climate only (unitless)
176	REAL(r_std), DIMENSION(npts) :: estab_rate_max_climate_grass !! maximum grass establishment rate,
177	!! based on climate only (unitless)
178	REAL(r_std), DIMENSION(npts) :: estab_rate_max_tree !! maximum tree establishment rate,
179	!! based on climate and fpc
180	!! (unitless)
181	REAL(r_std), DIMENSION(npts) :: estab_rate_max_grass !! maximum grass establishment rate,
182	!! based on climate and fpc
183	!! (unitless)
184	REAL(r_std), DIMENSION(npts) :: sumfpc !! total natural fpc (unitless)
185	REAL(r_std), DIMENSION(npts) :: fracnat !! total fraction occupied by natural
186	!! vegetation (unitless)
187	REAL(r_std), DIMENSION(npts) :: sumfpc_wood !! total woody fpc (unitless)
188	REAL(r_std), DIMENSION(npts) :: spacefight_grass!! for grasses, measures the total concurrence
189	!! for available space (unitless)
190	REAL(r_std), DIMENSION(npts,nvm) :: d_ind !! change in number of individuals per time step
191	!! (individuals m^{-2} day{-1})
192	REAL(r_std), DIMENSION(npts,nelements) :: bm_new !! biomass increase (gC m^{-2 })
193	REAL(r_std), DIMENSION(npts,nvm,nparts,nelements) :: biomass_old !! Save the original biomass passed into the subroutine
194	REAL(r_std), DIMENSION(npts,nelements) :: bm_non !! Non-effective establishment: the "virtual" saplings that die instantly
195	REAL(r_std), DIMENSION(npts,nelements) :: bm_eff !! Effective (or real) establishment
196	REAL(r_std), DIMENSION(npts) :: dia !! stem diameter (m)
197	REAL(r_std), DIMENSION(npts) :: b1 !! temporary variable
198	REAL(r_std), DIMENSION(npts) :: woodmass !! woodmass of an individual (gC m^{-2})
199	REAL(r_std), DIMENSION(npts) :: leaf_mass_young !! carbon mass in youngest leaf age class
200	!! (gC m^{-2})
201	REAL(r_std), DIMENSION(npts) :: factor !! reduction factor for establishment if many
202	!! trees or grasses are present (unitless)
203	REAL(r_std), DIMENSION(npts,nelements) :: total_bm !! Total carbon (or nitrogen) mass for all pools (gC m^{-2})
204	REAL(r_std), DIMENSION(npts,nelements) :: total_bm_sapl !! Total sappling biomass for all pools
205	!! (gC(orN) m^{-2})
206	REAL(r_std), DIMENSION(npts,nelements) :: total_bm_sapl_non !! total non-effective sapling biomass
207
208	INTEGER(i_std) :: nfrontx !! from how many sides is the grid box invaded
209	!! (unitless?)
210	INTEGER(i_std) :: nfronty !! from how many sides is the grid box invaded
211	!! (unitless?)
212	REAL(r_std), DIMENSION(npts) :: vn !! flow due to new individuals veget_cov_max after
213	!! establishment, to get a proper estimate of
214	!! carbon and nitrogen
215	REAL(r_std), DIMENSION(npts) :: lai_ind !! lai on each PFT surface (m^2 m^{-2})
216	REAL(r_std), DIMENSION(npts) :: veget_cov_max_tree !! Sum of veget_cov_max for the pft's which are trees
217	INTEGER(i_std) :: nbtree !! Number of PFT's which are trees
218	INTEGER(i_std) :: i,j,k,m !! indices (unitless)
219	REAL(r_std), DIMENSION(npts) :: soil_n_constraint !! Soil inorganic N pool (gN m-2)
220	REAL(r_std), DIMENSION(npts) :: ntake !! total N associated to establishment (gN)
221	REAL(r_std), DIMENSION(npts) :: ntake_nitrate !! nitrate associated to establishment (gN)
222	!_ ================================================================================================================================
223
224	IF (printlev>=3) WRITE(numout,*) 'Entering establish'
225
226	!! 1. messages and initialization
227
228	! time during which young plants can be completely eaten by herbivores after germination
229	! (and then individual die) assume to be half year
230	! No reference
231	tau_eatup = one_year/2.
232
233	! Calculate the sum of the vegetation over the tree pft's and the number of pft's which are trees
234	veget_cov_max_tree(:) = 0.0
235	nbtree=0
236	DO j = 1, nvm
237	IF (is_tree(j)) THEN
238	veget_cov_max_tree(:) = veget_cov_max_tree(:) + veget_cov_max(:,j)
239	nbtree = nbtree + 1
240	END IF
241	END DO
242	! Set nbtree=1 to avoid zero division later if there are no tree PFT's.
243	! For that case veget_cov_max_tree=0 so there will not be any impact.
244	IF (nbtree == 0) nbtree=1
245
246	!! 1.1 First call only
247	IF ( firstcall_establish ) THEN
248
249	WRITE(numout,*) 'establish:'
250
251	WRITE(numout,*) ' > time during which a sapling can be entirely eaten by herbivores (d): ', &
252	tau_eatup
253
254	firstcall_establish = .FALSE.
255
256	ENDIF
257
258	!! 2. recalculate fpc
259
260	IF (ok_dgvm) THEN
261	fracnat(:) = un
262
263	!! 2.1 Only natural part of the grid cell
264	do j = 2,nvm ! Loop over # PFTs
265
266	IF ( .NOT. natural(j) ) THEN
267	fracnat(:) = fracnat(:) - veget_cov_max(:,j)
268	ENDIF
269	ENDDO ! Loop over # PFTs
270
271	sumfpc(:) = zero
272
273	!! 2.2 Total natural fpc on grid
274	! The overall fractional coverage of a PFT in a grid is calculated here.
275	! FPC is related to mean individual leaf area index by the Lambert-Beer law.
276	! See Eq. (1) in tex file.\n
277	DO j = 2,nvm ! Loop over # PFTs
278	IF ( natural(j) ) THEN
279	WHERE(fracnat(:).GT.min_stomate)
280	WHERE (lai(:,j) == val_exp)
281	fpc_nat(:,j) = cn_ind(:,j) * ind(:,j) / fracnat(:)
282	ELSEWHERE
283	fpc_nat(:,j) = cn_ind(:,j) * ind(:,j) / fracnat(:) &
284	* ( un - exp( - biomass_to_lai(lm_lastyearmax(:,j),npts,j)* ext_coeff(j) ) )
285	ENDWHERE
286	ENDWHERE
287
288	WHERE ( PFTpresent(:,j) )
289	sumfpc(:) = sumfpc(:) + fpc_nat(:,j)
290	ENDWHERE
291	ELSE
292
293	fpc_nat(:,j) = zero
294
295	ENDIF
296
297	ENDDO ! Loop over # PFTs
298
299	!! 2.3 Total woody fpc on grid and number of regenerative tree pfts
300	! Total woody FPC increases by adding new FPC.
301	! Under the condition that temperature in last winter is higher than a threshold,
302	! woody plants is exposed in higher competitive environment.
303	sumfpc_wood(:) = zero
304
305	DO j = 2,nvm ! Loop over # PFTs
306
307	IF ( is_tree(j) .AND. natural(j) ) THEN
308
309	! total woody fpc
310	WHERE ( PFTpresent(:,j) )
311	sumfpc_wood(:) = sumfpc_wood(:) + fpc_nat(:,j)
312	ENDWHERE
313
314	ENDIF
315
316	ENDDO ! Loop over # PFTs
317
318	!! 2.4 Total number of natural grasses on grid\n
319	! Grass increment equals 'everywhere'\n
320	spacefight_grass(:) = zero
321
322	DO j = 2,nvm ! Loop over # PFTs
323
324	IF ( .NOT. is_tree(j) .AND. natural(j) ) THEN
325
326	! Count a PFT fully only if it is present on a grid.
327	WHERE ( PFTpresent(:,j) )
328	spacefight_grass(:) = spacefight_grass(:) + everywhere(:,j)
329	ENDWHERE
330
331	ENDIF
332
333	ENDDO ! Loop over # PFTs
334
335	!! 2.5 Maximum establishment rate, based on climate only\n
336	WHERE ( ( precip_annual(:) .GE. precip_crit ) .AND. ( gdd0(:) .GE. gdd_crit_estab ) )
337
338	estab_rate_max_climate_tree(:) = estab_max_tree ! 'estab_max_*'; see 'stomate_constants.f90'
339	estab_rate_max_climate_grass(:) = estab_max_grass
340
341	ELSEWHERE
342
343	estab_rate_max_climate_tree(:) = zero
344	estab_rate_max_climate_grass(:) = zero
345
346	ENDWHERE
347
348	!! 2.6 Reduce maximum tree establishment rate if many trees are present.
349	! In the original DGVM, this is done using a step function which yields a
350	! reduction by factor 4 if sumfpc_wood(i) .GT. fpc_crit - 0.05.
351	! This can lead to small oscillations (without consequences however).
352	! Here, a steady linear transition is used between fpc_crit-0.075 and
353	! fpc_crit-0.025.
354	! factor(:) = 1. - 15. * ( sumfpc_wood(:) - (fpc_crit-.075))
355	! factor(:) = MAX( 0.25_r_std, MIN( 1._r_std, factor(:)))
356	! S. Zaehle modified according to Smith et al. 2001
357	! See Eq. (2) in header
358	factor(:)=(un - exp(- establish_scal_fact * (un - sumfpc_wood(:))))*(un - sumfpc_wood(:))
359	estab_rate_max_tree(:) = estab_rate_max_climate_tree(:) * factor(:)
360
361	!! 2.7 Modulate grass establishment rate.
362	! If canopy is not closed (fpc < fpc_crit-0.05), normal establishment.
363	! If canopy is closed, establishment is reduced by a factor 4.
364	! Factor is linear between these two bounds.
365	! This is different from the original DGVM where a step function is
366	! used at fpc_crit-0.05 (This can lead to small oscillations,
367	! without consequences however).
368	! factor(:) = 1. - 15. * ( sumfpc(:) - (fpc_crit-.05))
369	! factor(:) = MAX( 0.25_r_std, MIN( 1._r_std, factor(:)))
370	! estab_rate_max_grass(:) = estab_rate_max_climate_grass(:) * factor(:)
371	! S. Zaehle modified to true LPJ formulation, grasses are only allowed in the
372	! fpc fraction not occupied by trees..., 080806
373	! estab_rate_max_grass(:)=MAX(0.98-sumfpc(:),zero)
374	! See Eq. (3) in header
375	estab_rate_max_grass(:) = MAX(MIN(estab_rate_max_climate_grass(:), max_tree_coverage - sumfpc(:)),zero)
376
377	!! 2.8 Longterm grass NPP for competition between C4 and C3 grasses
378	! to avoid equal veget_cov_max, the idea is that more reestablishment
379	! is possible for the more productive PFT
380	factor(:) = min_stomate
381
382	DO j = 2,nvm ! Loop over # PFTs
383	IF ( natural(j) .AND. .NOT.is_tree(j)) &
384	factor(:) = factor(:) + npp_longterm(:,j) * &
385	biomass_to_lai(lm_lastyearmax(:,j),npts,j)
386	ENDDO ! Loop over # PFTs
387
388	!! 2.9 Establish natural PFTs
389	d_ind(:,:) = zero
390
391	IF ( NbNeighb /= 8 ) THEN
392	CALL ipslerr(3, "establish", "This routine needs to be adapted to non rectengular grids", "Talk to Jan Polcher", " ")
393	ENDIF
394
395	DO j = 2,nvm ! Loop over # PFTs
396
397	IF ( natural(j) ) THEN ! only for natural PFTs
398
399	!! 2.9.1 PFT expansion across the grid box. Not to be confused with areal coverage.
400	IF ( treat_expansion ) THEN
401
402	! only treat plants that are regenerative and present and still can expand
403	DO i = 1, npts ! Loop over # pixels - domain size
404
405	IF ( PFTpresent(i,j) .AND. &
406	( everywhere(i,j) .LT. un ) .AND. &
407	( regenerate(i,j) .GT. regenerate_crit ) ) THEN
408
409	! from how many sides is the grid box invaded (separate x and y directions
410	! because resolution may be strongly anisotropic)
411	! For the moment we only look into 4 direction but that can be expanded (JP)
412	nfrontx = 0
413	IF ( neighbours(i,3) .GT. 0 ) THEN
414	IF ( everywhere(neighbours(i,3),j) .GT. 1.-min_stomate ) nfrontx = nfrontx+1
415	ENDIF
416	IF ( neighbours(i,7) .GT. 0 ) THEN
417	IF ( everywhere(neighbours(i,7),j) .GT. 1.-min_stomate ) nfrontx = nfrontx+1
418	ENDIF
419
420	nfronty = 0
421	IF ( neighbours(i,1) .GT. 0 ) THEN
422	IF ( everywhere(neighbours(i,1),j) .GT. 1.-min_stomate ) nfronty = nfronty+1
423	ENDIF
424	IF ( neighbours(i,5) .GT. 0 ) THEN
425	IF ( everywhere(neighbours(i,5),j) .GT. 1.-min_stomate ) nfronty = nfronty+1
426	ENDIF
427
428	everywhere(i,j) = &
429	everywhere(i,j) + migrate(j) * dt/one_year * &
430	( nfrontx / resolution(i,1) + nfronty / resolution(i,2) )
431
432	IF ( .NOT. need_adjacent(i,j) ) THEN
433
434	! in that case, we also assume that the PFT expands from places within
435	! the grid box (e.g., oasis).
436	! What is this equation? No reference.
437	everywhere(i,j) = &
438	everywhere(i,j) + migrate(j) * dt/one_year * &
439	2. * SQRT( pi*everywhere(i,j)/(area(i)) )
440
441	ENDIF
442
443	everywhere(i,j) = MIN( everywhere(i,j), un )
444
445	ENDIF
446
447	ENDDO ! Loop over # pixels - domain size
448
449	ENDIF ! treat expansion?
450
451	!! 2.9.2 Establishment rate
452	! - Is lower if the PFT is only present in a small part of the grid box
453	! (after its introduction), therefore multiplied by "everywhere".
454	! - Is divided by the number of PFTs that compete ("spacefight").
455	! - Is modulated by space availability (avail_tree, avail_grass).
456
457	!! 2.9.2.1 present and regenerative trees
458	IF ( is_tree(j) ) THEN
459
460	WHERE ( PFTpresent(:,j) .AND. ( regenerate(:,j) .GT. regenerate_crit ) )
461	d_ind(:,j) = estab_rate_max_tree(:)everywhere(:,j) &
462	avail_tree(:) * dt/one_year
463	ENDWHERE
464
465	!! 2.9.2.2 present and regenerative grasses
466	ELSE
467	WHERE ( PFTpresent(:,j) .AND. ( regenerate(:,j) .GT. regenerate_crit ) &
468	.AND.factor(:).GT.min_stomate .AND. spacefight_grass(:).GT. min_stomate)
469
470	d_ind(:,j) = estab_rate_max_grass(:)everywhere(:,j)/spacefight_grass(:) &
471	MAX(min_stomate,npp_longterm(:,j)biomass_to_lai(lm_lastyearmax(:,j),npts,j)/factor(:)) &
472	fracnat(:) * dt/one_year
473	ENDWHERE
474
475	ENDIF ! tree/grass
476
477	ENDIF ! if natural
478	ENDDO ! Loop over # PFTs
479
480	!! 3. Lpj establishment in static case
481
482	! Lpj establishment in static case, S. Zaehle 080806, account for real LPJ dynamics in
483	! prescribed vegetation, i.e. population dynamics within a given area of the grid cell.
484	ELSE
485
486	d_ind(:,:) = zero
487
488	DO j = 2,nvm ! Loop over # PFTs
489	WHERE(ind(:,j)*cn_ind(:,j).GT.min_stomate)
490	lai_ind(:) = biomass_to_lai(lm_lastyearmax(:,j),npts,j)/(ind(:,j)*cn_ind(:,j))
491	ELSEWHERE
492	lai_ind(:) = zero
493	ENDWHERE
494
495	!! 3.1 For natural woody PFTs
496	IF ( natural(j) .AND. is_tree(j)) THEN
497
498	! See Eq. (4) in tex file.
499	fpc_nat(:,j) = MIN(un, cn_ind(:,j) * ind(:,j) * &
500	MAX( ( un - exp( - ext_coeff(j) * lai_ind(:) ) ), min_cover ) )
501
502
503	WHERE (veget_cov_max(:,j).GT.min_stomate.AND.ind(:,j).LE.2.)
504
505	!! 3.1.1 Only establish into growing stands
506	! Only establish into growing stands, ind can become very
507	! large in the static mode because LAI is very low in poor
508	! growing conditions, favouring continuous establishment.
509	! To avoid this a maximum IND is set. BLARPP: This should be
510	! replaced by a better stand density criteria.
511	factor(:)=(un - exp(-establish_scal_fact * (un - fpc_nat(:,j))))*(un - fpc_nat(:,j))
512
513	estab_rate_max_tree(:) = estab_max_tree * factor(:) * MAX(0.11 , nstress_season(:,j) )
514
515	!! 3.1.2 do establishment for natural PFTs\n
516	d_ind(:,j) = MAX( zero, estab_rate_max_tree(:) * dt/one_year)
517
518	ENDWHERE
519
520	!S. Zaehle: quickfix: to simulate even aged stand, uncomment the following lines...
521	!where (ind(:,j) .LE. min_stomate)
522	!d_ind(:,j) = 0.1 !MAX( 0.0, estab_rate_max_tree(:) * dt/one_year)
523	WHERE (veget_cov_max(:,j).GT.min_stomate .AND. ind(:,j).EQ.zero)
524	d_ind(:,j) = ind_0_estab
525	ENDWHERE
526
527	!! 3.2 For natural grass PFTs
528	ELSEIF ( natural(j) .AND. .NOT.is_tree(j)) THEN
529
530	WHERE (veget_cov_max(:,j).GT.min_stomate)
531
532	fpc_nat(:,j) = cn_ind(:,j) * ind(:,j) * &
533	MAX( ( un - exp( - ext_coeff(j) * lai_ind(:) ) ), min_cover )
534
535	d_ind(:,j) = MAX(zero , (un - fpc_nat(:,j)) * dt/one_year )
536
537	ENDWHERE
538
539	WHERE (veget_cov_max(:,j).GT.min_stomate .AND. ind(:,j).EQ. zero)
540	d_ind(:,j) = ind_0_estab
541	ENDWHERE
542
543	ENDIF
544
545	ENDDO ! Loop over # PFTs
546
547	ENDIF ! DGVM OR NOT
548
549	!! 4. Biomass calculation
550
551	DO j = 2,nvm ! Loop over # PFTs
552
553	IF ( natural(j) ) THEN ! only for natural PFTs
554
555	!! 4.1 Herbivores reduce establishment rate
556	! We suppose that saplings are vulnerable during a given time after establishment.
557	! This is taken into account by preventively reducing the establishment rate.
558	IF ( ok_herbivores ) THEN
559
560	d_ind(:,j) = d_ind(:,j) * EXP( - tau_eatup/herbivores(:,j) )
561
562	ENDIF
563
564	!! 4.2 Total biomass.
565	! Add biomass only if d_ind, over one year, is of the order of ind.
566	! save old leaf mass to calculate leaf age
567	leaf_mass_young(:) = leaf_frac(:,j,1) * biomass(:,j,ileaf,icarbon)
568
569	! total biomass of existing PFT to limit biomass added from establishment
570	total_bm(:,:) = zero
571
572	DO k = 1, nparts
573	total_bm(:,:) = total_bm(:,:) + biomass(:,j,k,:)
574	ENDDO
575	! soil nitrogen available to feed plant establishment.
576	! This needs to be taken from soil N. Reason: plants must not produce nitrogen by establishment
577	! (otherwise this creates N)
578	soil_n_constraint(:)=soil_n_min(:,j,iammonium)+soil_n_min(:,j,initrate)
579
580
581	IF(ok_dgvm) THEN
582	vn(:) = veget_cov_max(:,j)
583	ELSE
584	vn(:) = un
585	ENDIF
586
587	!! 4.3 Woodmass calculation
588
589	!! 4.3.1 with DGVM
590	IF(ok_dgvm) THEN
591
592	! S. Zaehle calculate new woodmass_ind and veget_cov_max after establishment (needed for correct scaling!)
593	! essential correction for MERGE!
594	IF(is_tree(j))THEN
595	DO i=1,npts ! Loop over # pixels - domain size
596	IF((d_ind(i,j)+ind(i,j)).GT.min_stomate) THEN
597
598	IF((total_bm(i,icarbon).LE.min_stomate) .OR. (veget_cov_max(i,j) .LE. min_stomate)) THEN
599
600	! new wood mass of PFT
601	woodmass_ind(i,j) = &
602	(((biomass(i,j,isapabove,icarbon) + biomass(i,j,isapbelow,icarbon) &
603	+ biomass(i,j,iheartabove,icarbon) + biomass(i,j,iheartbelow,icarbon))*veget_cov_max(i,j)) &
604	+ (bm_sapl_2D(i,j,isapabove,icarbon) + bm_sapl_2D(i,j,isapbelow,icarbon) &
605	+ bm_sapl_2D(i,j,iheartabove,icarbon) + &
606	bm_sapl_2D(i,j,iheartbelow,icarbon))*d_ind(i,j))/(ind(i,j) + d_ind(i,j))
607
608	ELSE
609
610	! new biomass is added to the labile pool, hence there is no change
611	! in CA associated with establishment
612	woodmass_ind(i,j) = &
613	& (biomass(i,j,isapabove,icarbon) + biomass(i,j,isapbelow,icarbon) &
614	& +biomass(i,j,iheartabove,icarbon) + biomass(i,j,iheartbelow,icarbon))*veget_cov_max(i,j) &
615	& /(ind(i,j) + d_ind(i,j))
616
617	ENDIF
618
619	! new diameter of PFT
620	dia(i) = (woodmass_ind(i,j)/(pipe_density(j)pi/4.pipe_tune2(j))) &
621	& **(1./(2.+pipe_tune3(j)))
622
623	!+++HACK+++
624	! Calculation of crown area has been changed
625	! Use the function wood_to_diahor
626	! ORIGINAL line of code : vn(i) = (ind(i,j) + d_ind(i,j))pipe_tune1(j)MIN(dia(i),maxdia(j))**pipe_tune_exp_coeff(j)
627	! Quick fix to make the model compile
628	vn(i) = (ind(i,j) + d_ind(i,j))crown_vertohor_dia(j)crown_to_height(j)pipe_tune2(j)(dia(i)**pipe_tune3(j))
629	!+++++++++++
630
631	ENDIF
632	ENDDO ! Loop over # pixels - domain size
633	ELSE ! for grasses, cnd=1, so the above calculation cancels
634	vn(:) = ind(:,j) + d_ind(:,j)
635	ENDIF
636
637	!! 4.3.2 without DGVM (static)\n
638	ELSE
639	DO i=1,npts ! Loop over # pixels - domain size
640	IF(is_tree(j).AND.(d_ind(i,j)+ind(i,j)).GT.min_stomate) THEN
641	IF(total_bm(i,icarbon).LE.min_stomate) THEN
642
643	! new wood mass of PFT
644	woodmass_ind(i,j) = &
645	& (((biomass(i,j,isapabove,icarbon) + biomass(i,j,isapbelow,icarbon) &
646	& + biomass(i,j,iheartabove,icarbon) + biomass(i,j,iheartbelow,icarbon))) &
647	& + (bm_sapl_2D(i,j,isapabove,icarbon) + bm_sapl_2D(i,j,isapbelow,icarbon) &
648	& + bm_sapl_2D(i,j,iheartabove,icarbon) + &
649	& bm_sapl_2D(i,j,iheartbelow,icarbon))*d_ind(i,j))/(ind(i,j)+d_ind(i,j))
650
651	ELSE
652
653	! new biomass is added to the labile pool, hence there is no change
654	! in CA associated with establishment
655	woodmass_ind(i,j) = &
656	& (biomass(i,j,isapabove,icarbon) + biomass(i,j,isapbelow,icarbon) &
657	& + biomass(i,j,iheartabove,icarbon) + biomass(i,j,iheartbelow,icarbon)) &
658	& /(ind(i,j) + d_ind(i,j))
659
660	ENDIF
661	ENDIF
662	ENDDO ! Loop over # pixels - domain size
663
664	vn(:) = un ! cannot change in static!, and veget_cov_max implicit in d_ind
665
666	ENDIF
667
668	!! 4.4 total biomass of PFT added by establishment defined over veget_cov_max ...
669
670	total_bm_sapl(:,:) = zero
671	total_bm_sapl_non(:,:) = zero
672	biomass_old(:,j,:,:)=biomass(:,j,:,:)
673	DO k = 1, nparts ! Loop over # litter tissues (nparts=8); see 'stomate_constants.f90'
674	WHERE(d_ind(:,j).GT.min_stomate.AND.total_bm(:,icarbon).GT.min_stomate.AND.vn(:).GT.min_stomate)
675
676	total_bm_sapl(:,icarbon) = total_bm_sapl(:,icarbon) + bm_sapl_2D(:,j,k,icarbon) * d_ind(:,j) / vn(:)
677	total_bm_sapl(:,initrogen) = total_bm_sapl(:,initrogen) + bm_sapl_2D(:,j,k,initrogen) * d_ind(:,j) / vn(:)
678
679	! non-effective establishment
680	total_bm_sapl_non(:,icarbon) = total_bm_sapl_non(:,icarbon) + &
681	bm_sapl_2D(:,j,k,icarbon) * (ind(:,j)+d_ind(:,j))*mortality(:,j) / vn(:)
682
683	total_bm_sapl_non(:,initrogen) = total_bm_sapl_non(:,initrogen) + &
684	bm_sapl_2D(:,j,k,initrogen) * (ind(:,j)+d_ind(:,j))*mortality(:,j) / vn(:)
685
686
687	ENDWHERE
688	ENDDO ! Loop over # litter tissues
689
690
691	total_bm_sapl(:,initrogen)=MIN(total_bm_sapl(:,initrogen),soil_n_constraint(:))
692	ntake(:)=zero
693	!Dan Zhu modification: there is a problem here, if DGVM is activated, co2_to_bm will never reach
694	!0 due to establishment, where d_ind is still large at equilibrium (=ind*mortality). So we
695	!need to subtract it from litter (not biomass, because the
696	!corresponding biomass has been lost in lpj_gap).
697
698	!! 4.5 Update biomass at each component
699	DO k = 1, nparts ! Loop over # litter tissues
700
701	bm_new(:,:) = zero
702	bm_non(:,:) = zero
703	bm_eff(:,:) = zero
704
705	! first ever establishment, C flows
706	WHERE( d_ind(:,j).GT.min_stomate .AND. &
707	total_bm(:,icarbon).LE.min_stomate .AND. &
708	veget_cov_max(:,j).GT.min_stomate)
709
710	bm_new(:,icarbon) = d_ind(:,j) * bm_sapl_2D(:,j,k,icarbon) / vn(:)
711	biomass(:,j,k,icarbon) = biomass(:,j,k,icarbon) + bm_new(:,icarbon)
712
713	bm_new(:,initrogen) = d_ind(:,j) * (bm_sapl_2D(:,j,k,initrogen)) / vn(:)
714	biomass(:,j,k,initrogen) = biomass(:,j,k,initrogen) + bm_new(:,initrogen)
715
716	! bm_to_litter minus the 'non-effective' establishment (mortality), but cannot be less than 0
717	WHERE((veget_cov_max_tree(:) .GT. 0.1) .AND. (veget_cov_max(:,j) .LT. veget_cov_max_tree(:)/nbtree) )
718
719	bm_non(:,icarbon) = MIN( biomass(:,j,k,icarbon)+bm_to_litter(:,j,k,icarbon), &
720	(ind(:,j)+d_ind(:,j))mortality(:,j) bm_sapl_2D(:,j,k,icarbon)/vn(:) )
721	bm_eff(:,icarbon) = MIN( npp_longterm(:,j)/one_year, bm_new(:,icarbon)-bm_non(:,icarbon) )
722	bm_non(:,icarbon) = MIN( biomass(:,j,k,icarbon)+bm_to_litter(:,j,k,icarbon), &
723	bm_new(:,icarbon)-bm_eff(:,icarbon) )
724
725	bm_non(:,initrogen) = MIN( biomass(:,j,k,initrogen)+bm_to_litter(:,j,k,initrogen), &
726	(ind(:,j)+d_ind(:,j))mortality(:,j) bm_sapl_2D(:,j,k,initrogen)/vn(:) )
727	bm_eff(:,initrogen) = bm_new(:,initrogen)-bm_non(:,initrogen)
728	bm_non(:,initrogen) = MIN( biomass(:,j,k,initrogen)+bm_to_litter(:,j,k,initrogen), &
729	bm_new(:,initrogen)-bm_eff(:,initrogen) )
730	! I think there was an error
731	! co2_to_bm(:,j)=co2_to_bm(:,j) + bm_new(:) - bm_non(:)
732	co2_to_bm(:,j)=co2_to_bm(:,j) + bm_new(:,icarbon)/dt - bm_non(:,icarbon)/dt
733	ntake(:) = ntake(:) + bm_new(:,initrogen) - bm_non(:,initrogen)
734
735	WHERE( bm_to_litter(:,j,k,icarbon) .LT. bm_non(:,icarbon) )
736	biomass(:,j,k,icarbon) = biomass(:,j,k,icarbon) - &
737	( bm_non(:,icarbon) - bm_to_litter(:,j,k,icarbon) )
738	ENDWHERE
739	bm_to_litter(:,j,k,icarbon) = bm_to_litter(:,j,k,icarbon) - &
740	MIN(bm_to_litter(:,j,k,icarbon), bm_non(:,icarbon) )
741
742	WHERE( bm_to_litter(:,j,k,initrogen) .LT. bm_non(:,initrogen) )
743	biomass(:,j,k,initrogen) = biomass(:,j,k,initrogen) - &
744	( bm_non(:,initrogen) - bm_to_litter(:,j,k,initrogen) )
745	ENDWHERE
746	bm_to_litter(:,j,k,initrogen) = bm_to_litter(:,j,k,initrogen) - &
747	MIN(bm_to_litter(:,j,k,initrogen), bm_non(:,initrogen) )
748
749	ELSEWHERE
750
751	bm_non(:,icarbon) = MIN( bm_to_litter(:,j,k,icarbon), &
752	(ind(:,j)+d_ind(:,j))mortality(:,j) bm_sapl_2D(:,j,k,icarbon)/vn(:) )
753	co2_to_bm(:,j)=co2_to_bm(:,j) + bm_new(:,icarbon)/dt - bm_non(:,icarbon)/dt
754	bm_to_litter(:,j,k,icarbon)=bm_to_litter(:,j,k,icarbon)- bm_non(:,icarbon)
755
756	bm_non(:,initrogen) = MIN( bm_to_litter(:,j,k,initrogen), &
757	(ind(:,j)+d_ind(:,j))mortality(:,j) bm_sapl_2D(:,j,k,initrogen)/vn(:) )
758	ntake(:)= ntake(:) + bm_new(:,initrogen) - bm_non(:,initrogen)
759	bm_to_litter(:,j,k,initrogen)=bm_to_litter(:,j,k,initrogen)- bm_non(:,initrogen)
760
761	ENDWHERE
762
763	ENDWHERE
764
765	! establishment into existing population, C flows
766	WHERE(d_ind(:,j).GT.min_stomate.AND.total_bm(:,icarbon).GT.min_stomate)
767
768	bm_new(:,icarbon) = total_bm_sapl(:,icarbon) * biomass_old(:,j,k,icarbon) / total_bm(:,icarbon)
769	biomass(:,j,k,icarbon) = biomass(:,j,k,icarbon) + bm_new(:,icarbon)
770
771	bm_new(:,initrogen) = total_bm_sapl(:,initrogen) * biomass_old(:,j,k,initrogen) / total_bm(:,initrogen)
772	biomass(:,j,k,initrogen) = biomass(:,j,k,initrogen) + bm_new(:,initrogen)
773
774	WHERE((veget_cov_max_tree(:) .GT. 0.1) .AND. (veget_cov_max(:,j) .LT. veget_cov_max_tree(:)/nbtree) )
775
776	bm_non(:,icarbon) = MIN( biomass(:,j,k,icarbon)+bm_to_litter(:,j,k,icarbon), &
777	total_bm_sapl_non(:,icarbon) *biomass_old(:,j,k,icarbon)/total_bm(:,icarbon) )
778	bm_eff(:,icarbon) = MIN( npp_longterm(:,j)/one_year, bm_new(:,icarbon)-bm_non(:,icarbon) )
779	bm_non(:,icarbon) = MAX( zero, MIN( biomass(:,j,k,icarbon)+bm_to_litter(:,j,k,icarbon)-min_stomate, &
780	bm_new(:,icarbon)-bm_eff(:,icarbon) ) )
781
782	bm_non(:,initrogen) = MIN( biomass(:,j,k,initrogen)+bm_to_litter(:,j,k,initrogen), &
783	total_bm_sapl_non(:,initrogen) *biomass_old(:,j,k,initrogen)/total_bm(:,initrogen) )
784	bm_eff(:,initrogen) = bm_new(:,initrogen)-bm_non(:,initrogen)
785	bm_non(:,initrogen) = MAX( zero, MIN( biomass(:,j,k,initrogen)+bm_to_litter(:,j,k,initrogen)-min_stomate, &
786	bm_new(:,initrogen)-bm_eff(:,initrogen) ) )
787
788	co2_to_bm(:,j)=co2_to_bm(:,j) + bm_new(:,icarbon) - bm_non(:,icarbon)
789	ntake(:) = ntake(:) + bm_new(:,initrogen) - bm_non(:,initrogen)
790
791	WHERE( bm_to_litter(:,j,k,icarbon) .LT. bm_non(:,icarbon) )
792	biomass(:,j,k,icarbon) = biomass(:,j,k,icarbon) - ( bm_non(:,icarbon) - bm_to_litter(:,j,k,icarbon) )
793	ENDWHERE
794	bm_to_litter(:,j,k,icarbon) = bm_to_litter(:,j,k,icarbon) - MIN(bm_to_litter(:,j,k,icarbon), bm_non(:,icarbon) )
795
796	WHERE( bm_to_litter(:,j,k,initrogen) .LT. bm_non(:,initrogen) )
797	biomass(:,j,k,initrogen) = biomass(:,j,k,initrogen) - ( bm_non(:,initrogen) - bm_to_litter(:,j,k,initrogen) )
798	ENDWHERE
799	bm_to_litter(:,j,k,initrogen) = bm_to_litter(:,j,k,initrogen) - &
800	MIN(bm_to_litter(:,j,k,initrogen), bm_non(:,initrogen) )
801
802	ELSEWHERE
803
804	bm_non(:,icarbon) = MIN( bm_to_litter(:,j,k,icarbon), &
805	total_bm_sapl_non(:,icarbon) *biomass_old(:,j,k,icarbon)/total_bm(:,icarbon) )
806	co2_to_bm(:,j) = co2_to_bm(:,j) + bm_new(:,icarbon)/dt - bm_non(:,icarbon)/dt
807	bm_to_litter(:,j,k,icarbon)=bm_to_litter(:,j,k,icarbon)- bm_non(:,icarbon)
808
809	bm_non(:,initrogen) = MIN( bm_to_litter(:,j,k,initrogen), &
810	total_bm_sapl_non(:,initrogen) *biomass_old(:,j,k,initrogen)/total_bm(:,initrogen) )
811	ntake(:) = ntake(:) + bm_new(:,initrogen) - bm_non(:,initrogen)
812	bm_to_litter(:,j,k,initrogen)=bm_to_litter(:,j,k,initrogen)- bm_non(:,initrogen)
813	ENDWHERE
814
815	ENDWHERE
816
817	ENDDO ! Loop over # litter tissues
818
819	ntake_nitrate(:)=MIN(ntake(:),MAX(soil_n_min(:,j,iammonium)-min_stomate,0.0))
820	soil_n_min(:,j,iammonium)=soil_n_min(:,j,iammonium)-ntake_nitrate(:)
821	n_uptake(:,j,iammonium)=n_uptake(:,j,iammonium) + ntake_nitrate(:)/dt
822
823	ntake(:)=MIN(ntake(:)-ntake_nitrate(:),MAX(soil_n_min(:,j,initrate)-min_stomate,0.0))
824	soil_n_min(:,j,initrate)=soil_n_min(:,j,initrate)-ntake(:)
825	n_uptake(:,j,initrate)=n_uptake(:,j,initrate) + ntake(:)/dt
826
827	IF (ANY( bm_to_litter(:,j,:,icarbon) .LT. 0.0 ) .OR. ANY( biomass(:,j,:,icarbon) .LT. 0.0 ) ) THEN
828	CALL ipslerr_p(3,'establish','something wrong in establish/gap.','','')
829	ENDIF
830
831	!! 4.6 Decrease leaf age in youngest class if new leaf biomass is higher than old one.
832	WHERE ( d_ind(:,j) * bm_sapl_2D(:,j,ileaf,icarbon) .GT. min_stomate )
833
834	! reset leaf ages. Should do a real calculation like in the npp routine,
835	! but this case is rare and not worth messing around.
836	! S. Zaehle 080806, added real calculation now, because otherwise leaf_age/leaf_frac
837	! are not initialised for the calculation of vmax, and hence no growth at all.
838	! logic follows that of stomate_npp.f90, just that it's been adjusted for the code here
839	leaf_age(:,j,1) = leaf_age(:,j,1) * leaf_mass_young(:) / &
840	( leaf_mass_young(:) + d_ind(:,j) * bm_sapl_2D(:,j,ileaf,icarbon) )
841
842	ENDWHERE
843
844	leaf_mass_young(:) = leaf_mass_young(:) + d_ind(:,j) * bm_sapl_2D(:,j,ileaf,icarbon)
845
846	!! 4.7 Youngest class: new mass in youngest class divided by total new mass
847	WHERE ( biomass(:,j,ileaf,icarbon) .GT. min_stomate )
848	! new age class fractions (fraction in youngest class increases)
849	leaf_frac(:,j,1) = leaf_mass_young(:) / biomass(:,j,ileaf,icarbon)
850
851	ENDWHERE
852
853	!! 4.8 Other classes: old mass in leaf age class divided by new mass
854	DO m = 2, nleafages
855
856	WHERE ( biomass(:,j,ileaf,icarbon) .GT. min_stomate )
857
858	leaf_frac(:,j,m) = leaf_frac(:,j,m) * &
859	( biomass(:,j,ileaf,icarbon) + d_ind(:,j) * bm_sapl_2D(:,j,ileaf,icarbon) ) / biomass(:,j,ileaf,icarbon)
860
861	ENDWHERE
862
863	ENDDO
864
865	!! 4.9 Update age and number of individuals
866	WHERE ( d_ind(:,j) .GT. min_stomate )
867
868	age(:,j) = age(:,j) * ind(:,j) / ( ind(:,j) + d_ind(:,j) )
869
870	ind(:,j) = ind(:,j) + d_ind(:,j)
871
872	ENDWHERE
873
874	!! 4.10 Convert excess sapwood to heartwood
875	!! No longer done : supressed by S. Zaehle given that the LPJ logic of carbon allocation was
876	!! contradictory to SLAVE allocation. See CVS tag 1_5 for initial formulation.
877
878
879	ENDIF ! natural
880
881	ENDDO ! Loop over # PFTs
882
883	!! 5. history
884
885	d_ind = d_ind / dt
886
887	CALL xios_orchidee_send_field("IND_ESTAB",d_ind)
888	CALL xios_orchidee_send_field("ESTABTREE",estab_rate_max_tree)
889	CALL xios_orchidee_send_field("ESTABGRASS",estab_rate_max_grass)
890
891	CALL histwrite_p (hist_id_stomate, 'IND_ESTAB', itime, d_ind, npts*nvm, horipft_index)
892	CALL histwrite_p (hist_id_stomate, 'ESTABTREE', itime, estab_rate_max_tree, npts, hori_index)
893	CALL histwrite_p (hist_id_stomate, 'ESTABGRASS', itime, estab_rate_max_grass, npts, hori_index)
894
895	IF (printlev>=4) WRITE(numout,*) 'Leaving establish'
896
897	END SUBROUTINE establish
898
899	END MODULE lpj_establish

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

Download in other formats: