Context Navigation

source: tags/ORCHIDEE_1_9_6/ORCHIDEE/src_stomate/lpj_light.f90 @ 5623

Last change on this file since 5623 was 720, checked in by didier.solyga, 12 years ago
Add svn headers for all modules. Improve documentation of the parameters. Replace two values by the corresponding parameters.
Property svn:keywords set to `HeadURL Date Author Revision`
File size: 21.4 KB

Line
1	! Light competition
2	!
3	! If canopy is almost closed (fpc > fpc_crit), then trees outcompete grasses.
4	! fpc_crit is normally fpc_crit.
5	! Here, fpc ("foilage protected cover") also takes into account the minimum fraction
6	! of space covered by trees through branches etc. This is done to prevent strong relative
7	! changes of fpc from one day to another for deciduous trees at the beginning of their
8	! growing season, which would yield to strong cutbacks (see 3.2.1.1.2)
9	! No competition between woody pfts (height of individuals is not considered) !
10	! Exception: when one woody pft is overwhelming (i.e. fpc > fpc_crit). In that
11	! case, eliminate all other woody pfts and reduce dominant pft to fpc_crit.
12	! Age of individuals is not considered. In reality, light competition would more
13	! easily kill young individuals, thus increasing the mean age of the stand.
14	! Exclude agricultural pfts from competition
15	!
16	! SZ: added light competition for the static case if the mortality is not
17	! assumed to be constant.
18	! other modifs:
19	! -1 FPC is now always calculated from lm_lastyearmax*sla, since the aim of this DGVM is
20	! to represent community ecology effects; seasonal variations in establishment related to phenology
21	! may be relevant, but beyond the scope of a 1st generation DGVM
22	! -2 problem, if agriculture is present, fpc can never reach 1.0 since natural veget_max < 1.0. To
23	! correct for this, ind must be recalculated to correspond to the natural density...
24	! since ind is 1/m2 grid cell, this can be achived by dividing ind by the agricultural fraction
25
26	!
27	!< $HeadURL$
28	!< $Date$
29	!< $Author$
30	!< $Revision$
31	! IPSL (2006)
32	! This software is governed by the CeCILL licence see ORCHIDEE/ORCHIDEE_CeCILL.LIC
33	!
34	MODULE lpj_light
35
36	! modules used:
37
38	USE ioipsl
39	USE constantes
40	USE stomate_data
41
42	IMPLICIT NONE
43
44	! private & public routines
45
46	PRIVATE
47	PUBLIC light, light_clear
48
49	! first call
50	LOGICAL, SAVE :: firstcall = .TRUE.
51
52	CONTAINS
53
54	SUBROUTINE light_clear
55	firstcall=.TRUE.
56	END SUBROUTINE light_clear
57
58	SUBROUTINE light (npts, dt, &
59	veget_max, fpc_max, PFTpresent, cn_ind, lai, maxfpc_lastyear, &
60	lm_lastyearmax, ind, biomass, veget_lastlight, bm_to_litter, mortality)
61
62	!
63	! 0 declarations
64	!
65
66	! 0.1 input
67
68	! Domain size
69	INTEGER(i_std), INTENT(in) :: npts
70	! Time step (days)
71	REAL(r_std), INTENT(in) :: dt
72	! Is pft there
73	LOGICAL, DIMENSION(npts,nvm), INTENT(in) :: PFTpresent
74	! crown area of individuals (m**2)
75	REAL(r_std), DIMENSION(npts,nvm), INTENT(in) :: cn_ind
76	! leaf area index OF AN INDIVIDUAL PLANT
77	REAL(r_std), DIMENSION(npts,nvm), INTENT(in) :: lai
78	! last year's maximum fpc for each natural PFT, on ground
79	REAL(r_std), DIMENSION(npts,nvm), INTENT(in) :: maxfpc_lastyear
80	! last year's maximum leafmass for each natural PFT, on ground
81	REAL(r_std), DIMENSION(npts,nvm), INTENT(in) :: lm_lastyearmax
82	! last year's maximum fpc for each natural PFT, on ground
83	REAL(r_std), DIMENSION(npts,nvm), INTENT(in) :: veget_max
84	! last year's maximum fpc for each natural PFT, on ground
85	REAL(r_std), DIMENSION(npts,nvm), INTENT(in) :: fpc_max
86
87	! 0.2 modified fields
88
89	! Number of individuals / m2
90	REAL(r_std), DIMENSION(npts,nvm), INTENT(inout) :: ind
91	! biomass (gC/(m**2 of ground))
92	REAL(r_std), DIMENSION(npts,nvm,nparts), INTENT(inout) :: biomass
93	! Vegetation cover after last light competition
94	REAL(r_std), DIMENSION(npts,nvm), INTENT(inout) :: veget_lastlight
95	! biomass taken away (gC/m**2)
96	REAL(r_std), DIMENSION(npts,nvm,nparts), INTENT(inout) :: bm_to_litter
97	! fraction of individuals that died this time step
98	REAL(r_std), DIMENSION(npts,nvm), INTENT(inout) :: mortality
99
100	! 0.3 local
101
102	! index
103	INTEGER(i_std) :: i,j,k,m
104	! total natural fpc
105	REAL(r_std), DIMENSION(npts) :: sumfpc
106	! fraction of natural vegetation at grid cell level
107	REAL(r_std), DIMENSION(npts) :: fracnat
108	! total natural woody fpc
109	REAL(r_std) :: sumfpc_wood
110	! change in total woody fpc
111	REAL(r_std) :: sumdelta_fpc_wood
112	! maximum wood fpc
113	REAL(r_std) :: maxfpc_wood
114	! which woody pft is maximum
115	INTEGER(i_std) :: optpft_wood
116	! total natural grass fpc
117	REAL(r_std) :: sumfpc_grass
118	! this year's foliage protected cover on natural part of the grid cell
119	REAL(r_std), DIMENSION(npts,nvm) :: fpc_nat
120	! fpc change within last year
121	REAL(r_std), DIMENSION(nvm) :: deltafpc
122	! Relative change of number of individuals for trees
123	REAL(r_std) :: reduct
124	! Fraction of plants that survive
125	REAL(r_std), DIMENSION(nvm) :: survive
126	! FPC for static mode
127	REAL(r_std), DIMENSION(npts) :: fpc_real
128	! FPC mortality for static mode
129	REAL(r_std), DIMENSION(npts) :: lai_ind
130	! number of grass PFTs present in the grid box
131	! INTEGER(i_std) :: num_grass
132	! New total grass fpc
133	REAL(r_std) :: sumfpc_grass2
134	! fraction of plants that dies each day (1/day)
135	REAL(r_std), DIMENSION(npts,nvm) :: light_death
136	! Relative change of number of individuals for trees
137	REAL(r_std) :: fpc_dec
138
139	! =========================================================================
140
141	IF (bavard.GE.3) WRITE(numout,*) 'Entering light'
142
143	!
144	! 1 first call
145	!
146
147	IF ( firstcall ) THEN
148
149	WRITE(numout,*) 'light:'
150
151	WRITE(numout,*) ' > For trees, minimum fraction of crown area covered'
152	WRITE(numout,*) ' (due to its branches etc.)', min_cover
153
154	WRITE(numout,*) ' > for diagnosis of fpc increase, compare today''s fpc'
155	IF ( annual_increase ) THEN
156	WRITE(numout,*) ' to last year''s maximum.'
157	ELSE
158	WRITE(numout,*) ' to fpc of the last time step.'
159	ENDIF
160
161	firstcall = .FALSE.
162
163	ENDIF
164
165	IF (control%ok_dgvm) THEN
166	!
167	! 2 fpc characteristics
168	!
169
170	! 2.0 Only natural part of the grid cell:
171	! calculate fraction of natural and agricultural (1-fracnat) surface
172
173	fracnat(:) = un
174	DO j = 2,nvm
175	IF ( .NOT. natural(j) ) THEN
176	fracnat(:) = fracnat(:) - veget_max(:,j)
177	ENDIF
178	ENDDO
179	!
180	! 2.1 calculate fpc on natural part of grid cell.
181	!
182	fpc_nat(:,:) = zero
183	fpc_nat(:,ibare_sechiba) = un
184
185	DO j = 2, nvm
186
187	IF ( natural(j) ) THEN
188
189	! 2.1.1 natural PFTs
190
191	IF ( tree(j) ) THEN
192
193	! 2.1.1.1 trees: minimum cover due to stems, branches etc.
194
195	! DO i = 1, npts
196	! IF (lai(i,j) == val_exp) THEN
197	! fpc_nat(i,j) = cn_ind(i,j) * ind(i,j)
198	! ELSE
199	! fpc_nat(i,j) = cn_ind(i,j) * ind(i,j) * &
200	! MAX( ( 1._r_std - exp( -lai(i,j) * ext_coeff(j) ) ), min_cover )
201	! ENDIF
202	! ENDDO
203
204	!NV : modif from SZ version : fpc is based on veget_max, not veget.
205	WHERE(fracnat(:).GE.min_stomate)
206	! WHERE(LAI(:,j) == val_exp)
207	! fpc_nat(:,j) = cn_ind(:,j) * ind(:,j) / fracnat(:)
208	! ELSEWHERE
209	! fpc_nat(:,j) = cn_ind(:,j) * ind(:,j) / fracnat(:) * &
210	! MAX( ( 1._r_std - exp( - lm_lastyearmax(:,j) * sla(j) * ext_coeff(j) ) ), min_cover )
211	! ENDWHERE
212	fpc_nat(:,j) = cn_ind(:,j) * ind(:,j) / fracnat(:)
213	ENDWHERE
214
215	ELSE
216
217	!NV : modif from SZ version : fpc is based on veget_max, not veget.
218	WHERE(fracnat(:).GE.min_stomate)
219	! WHERE(LAI(:,j) == val_exp)
220	! fpc_nat(:,j) = cn_ind(:,j) * ind(:,j) / fracnat(:)
221	! ELSEWHERE
222	! fpc_nat(:,j) = cn_ind(:,j) * ind(:,j) / fracnat(:) * &
223	! ( 1._r_std - exp( - lm_lastyearmax(:,j) * sla(j) * ext_coeff(j) ) )
224	! ENDWHERE
225	fpc_nat(:,j) = cn_ind(:,j) * ind(:,j) / fracnat(:)
226	ENDWHERE
227
228	!!$ ! 2.1.1.2 bare ground
229	!!$ IF (j == ibare_sechiba) THEN
230	!!$ fpc_nat(:,j) = cn_ind(:,j) * ind(:,j)
231	!!$
232	!!$ ! 2.1.1.3 grasses
233	!!$ ELSE
234	!!$ DO i = 1, npts
235	!!$ IF (lai(i,j) == val_exp) THEN
236	!!$ fpc_nat(i,j) = cn_ind(i,j) * ind(i,j)
237	!!$ ELSE
238	!!$ fpc_nat(i,j) = cn_ind(i,j) * ind(i,j) * &
239	!!$ ( 1._r_std - exp( -lai(i,j) * ext_coeff(j) ) )
240	!!$ ENDIF
241	!!$ ENDDO
242	!!$ ENDIF
243
244	ENDIF ! tree/grass
245
246	ELSE
247
248	! 2.1.2 agricultural PFTs: not present on natural part
249
250	fpc_nat(:,j) = zero
251
252	ENDIF ! natural/agricultural
253
254	ENDDO
255
256	!
257	! 2.2 sum natural fpc for every grid point
258	!
259
260	sumfpc(:) = zero
261	DO j = 2,nvm
262	!SZ bug correction MERGE: need to subtract agricultural area!
263	sumfpc(:) = sumfpc(:) + fpc_nat(:,j)
264	ENDDO
265
266	!
267	! 3 Light competition
268	!
269
270	light_death(:,:) = zero
271
272	DO i = 1, npts ! SZ why this loop and not a vector statement ?
273
274	! Only if vegetation cover is dense
275
276	IF ( sumfpc(i) .GT. fpc_crit ) THEN
277
278	! fpc change for each pft
279	! There are two possibilities: either we compare today's fpc with the fpc after the last
280	! time step, or we compare it to last year's maximum fpc of that PFT. In the first case,
281	! the fpc increase will be strong for seasonal PFTs at the beginning of the growing season.
282	! As for trees, the cutback is proportional to this increase, this means that seasonal trees
283	! will be disadvantaged compared to evergreen trees. In the original LPJ model, with its
284	! annual time step, the second method was used (this corresponds to annual_increase=.TRUE.)
285
286	IF ( annual_increase ) THEN
287	deltafpc(:) = MAX( (fpc_nat(i,:)-maxfpc_lastyear(i,:)), zero )
288	ELSE
289	deltafpc(:) = MAX( (fpc_nat(i,:)-veget_lastlight(i,:)), zero )
290	ENDIF
291
292	! default: survive
293
294	survive(:) = un
295
296	!
297	! 3.1 determine some characteristics of the fpc distribution
298	!
299
300	sumfpc_wood = zero
301	sumdelta_fpc_wood = zero
302	maxfpc_wood = zero
303	optpft_wood = 0
304	sumfpc_grass = zero
305	! num_grass = 0
306
307	DO j = 2,nvm
308
309	! only natural pfts
310
311	IF ( natural(j) ) THEN
312
313	IF ( tree(j) ) THEN
314
315	! trees
316
317	! total woody fpc
318
319	sumfpc_wood = sumfpc_wood + fpc_nat(i,j)
320
321	! how much did the woody fpc increase
322
323	sumdelta_fpc_wood = sumdelta_fpc_wood + deltafpc(j)
324
325	! which woody pft is preponderant
326
327	IF ( fpc_nat(i,j) .GT. maxfpc_wood ) THEN
328
329	optpft_wood = j
330
331	maxfpc_wood = fpc_nat(i,j)
332
333	ENDIF
334
335	ELSE
336
337	! grasses
338
339	! total (natural) grass fpc
340
341	sumfpc_grass = sumfpc_grass + fpc_nat(i,j)
342
343	! number of grass PFTs present in the grid box
344
345	! IF ( PFTpresent(i,j) ) THEN
346	! num_grass = num_grass + 1
347	! ENDIF
348
349	ENDIF ! tree or grass
350
351	ENDIF ! natural
352
353	ENDDO ! loop over pfts
354
355	!
356	! 3.2 light competition: assume wood outcompetes grass
357	!
358	!SZ
359	!!$ IF (sumfpc_wood .GE. fpc_crit ) THEN
360
361	!
362	! 3.2.1 all allowed natural space is covered by wood:
363	! cut back trees to fpc_crit.
364	! Original DGVM: kill grasses. Modified: we let a very
365	! small fraction of grasses survive.
366	!
367
368	DO j = 2,nvm
369
370	! only present and natural pfts compete
371
372	IF ( PFTpresent(i,j) .AND. natural(j) ) THEN
373
374	IF ( tree(j) ) THEN
375
376	!
377	! 3.2.1.1 tree
378	!
379
380	! no single woody pft is overwhelming
381	! (original DGVM: tree_mercy = 0.0 )
382	! The reduction rate is proportional to the ratio deltafpc/fpc.
383
384	IF (sumfpc_wood .GE. fpc_crit .AND. fpc_nat(i,j) .GT. min_stomate .AND. &
385	sumdelta_fpc_wood .GT. min_stomate) THEN
386
387	! reduct = MIN( ( ( deltafpc(j)/sumdelta_fpc_wood * &
388	! (sumfpc_wood-fpc_crit) ) / fpc_nat(i,j) ), &
389	! ( 1._r_std - 0.01 ) )
390	reduct = un - MIN((fpc_nat(i,j)-(sumfpc_wood-fpc_crit) &
391	* deltafpc(j)/sumdelta_fpc_wood)/fpc_nat(i,j), un )
392
393	ELSE
394
395	! tree fpc didn't icrease or it started from nothing
396
397	reduct = zero
398
399	ENDIF
400
401	survive(j) = un - reduct
402
403	ELSE
404
405	!
406	! 3.2.1.2 grass: let a very small fraction survive (the sum of all
407	! grass individuals may make up a maximum cover of
408	! grass_mercy [for lai -> infinity]).
409	! In the original DGVM, grasses were killed in that case,
410	! corresponding to grass_mercy = 0.
411	!
412
413	! survive(j) = ( 0.01 / REAL( num_grass,r_std ) ) / ind(i,j)
414
415	! survive(j) = MIN( 1._r_std, survive(j) )
416
417	IF(sumfpc_grass .GE. un-MIN(fpc_crit,sumfpc_wood).AND. &
418	sumfpc_grass.GE.min_stomate) THEN
419
420	fpc_dec = (sumfpc_grass - un + MIN(fpc_crit,sumfpc_wood))*fpc_nat(i,j)/sumfpc_grass
421
422	reduct = fpc_dec
423	ELSE
424	reduct = zero
425	ENDIF
426	survive(j) = ( un - reduct )
427	ENDIF ! tree or grass
428
429	ENDIF ! pft there and natural
430
431	ENDDO ! loop over pfts
432
433	!SZ
434	!!$ ELSE
435	!!$
436	!!$ !
437	!!$ ! 3.2.2 not too much wood so that grasses can subsist
438	!!$ !
439	!!$
440	!!$ ! new total grass fpc
441	!!$ sumfpc_grass2 = fpc_crit - sumfpc_wood
442	!!$
443	!!$ DO j = 2,nvm
444	!!$
445	!!$ ! only present and natural PFTs compete
446	!!$
447	!!$ IF ( PFTpresent(i,j) .AND. natural(j) ) THEN
448	!!$
449	!!$ IF ( tree(j) ) THEN
450	!!$
451	!!$ ! no change for trees
452	!!$
453	!!$ survive(j) = 1.0
454	!!$
455	!!$ ELSE
456	!!$
457	!!$ ! grass: fractional loss is the same for all grasses
458	!!$
459	!!$ IF ( sumfpc_grass .GT. min_stomate ) THEN
460	!!$ survive(j) = sumfpc_grass2 / sumfpc_grass
461	!!$ ELSE
462	!!$ survive(j)= zero
463	!!$ ENDIF
464	!!$
465	!!$ ENDIF
466	!!$
467	!!$ ENDIF ! pft there and natural
468	!!$
469	!!$ ENDDO ! loop over pfts
470	!!$
471	!!$ ENDIF ! sumfpc_wood > fpc_crit
472
473	!
474	! 3.3 update output variables
475	!
476
477	DO j = 2,nvm
478
479	IF ( PFTpresent(i,j) .AND. natural(j) ) THEN
480
481	bm_to_litter(i,j,:) = bm_to_litter(i,j,:) + &
482	biomass(i,j,:) * ( un - survive(j) )
483
484	biomass(i,j,:) = biomass(i,j,:) * survive(j)
485
486	IF ( control%ok_dgvm ) THEN
487	ind(i,j) = ind(i,j) * survive(j)
488	ENDIF
489
490	! fraction of plants that dies each day.
491	! exact formulation: light_death(i,j) = un - survive(j) / dt
492	light_death(i,j) = ( un - survive(j) ) / dt
493
494	ENDIF ! pft there and natural
495
496	ENDDO ! loop over pfts
497
498	ENDIF ! sumfpc > fpc_crit
499
500	ENDDO ! loop over grid points
501
502	!
503	! 4 recalculate fpc on natural part of grid cell (for next light competition)
504	!
505
506	DO j = 2,nvm
507
508	IF ( natural(j) ) THEN
509
510	!
511	! 4.1 natural PFTs
512	!
513
514	IF ( tree(j) ) THEN
515
516	! 4.1.1 trees: minimum cover due to stems, branches etc.
517
518	DO i = 1, npts
519	!NVMODIF
520	! IF (lai(i,j) == val_exp) THEN
521	! veget_lastlight(i,j) = cn_ind(i,j) * ind(i,j)
522	! ELSE
523	! veget_lastlight(i,j) = &
524	! cn_ind(i,j) * ind(i,j) * &
525	! MAX( ( un - exp( -lai(i,j) * ext_coeff(j) ) ), min_cover )
526	! ENDIF
527	!! veget_lastlight(i,j) = cn_ind(i,j) * ind(i,j)
528	IF (lai(i,j) == val_exp) THEN
529	veget_lastlight(i,j) = cn_ind(i,j) * ind(i,j)
530	ELSE
531	veget_lastlight(i,j) = &
532	cn_ind(i,j) * ind(i,j) * &
533	MAX( ( un - EXP( - lm_lastyearmax(i,j) * sla(j) * ext_coeff(j) ) ), min_cover )
534	ENDIF
535	ENDDO
536
537	ELSE
538
539	! 4.1.2 grasses
540	DO i = 1, npts
541	!NVMODIF
542	! IF (lai(i,j) == val_exp) THEN
543	! veget_lastlight(i,j) = cn_ind(i,j) * ind(i,j)
544	! ELSE
545	! veget_lastlight(i,j) = cn_ind(i,j) * ind(i,j) * &
546	! ( un - exp( -lai(i,j) * ext_coeff(j) ) )
547	! ENDIF
548	!!veget_lastlight(i,j) = cn_ind(i,j) * ind(i,j)
549	IF (lai(i,j) == val_exp) THEN
550	veget_lastlight(i,j) = cn_ind(i,j) * ind(i,j)
551	ELSE
552	veget_lastlight(i,j) = cn_ind(i,j) * ind(i,j) * &
553	( un - exp( - lm_lastyearmax(i,j) * sla(j) * ext_coeff(j) ) )
554	ENDIF
555	ENDDO
556	ENDIF ! tree/grass
557
558	ELSE
559
560	!
561	! 4.2 agricultural PFTs: not present on natural part
562	!
563
564	veget_lastlight(:,j) = zero
565
566	ENDIF ! natural/agricultural
567
568	ENDDO
569
570	ELSE ! static
571
572	light_death(:,:) = zero
573
574	DO j = 2, nvm
575
576	IF ( natural(j) ) THEN
577
578	! 2.1.1 natural PFTs, in the one PFT only case there needs to be no special case for grasses,
579	! neither a redistribution of mortality (delta fpc)
580
581	WHERE( ind(:,j)*cn_ind(:,j) .GT. min_stomate )
582	lai_ind(:) = sla(j) * lm_lastyearmax(:,j) / ( ind(:,j) * cn_ind(:,j) )
583	ELSEWHERE
584	lai_ind(:) = zero
585	ENDWHERE
586
587	fpc_nat(:,j) = cn_ind(:,j) * ind(:,j) * &
588	MAX( ( un - exp( - ext_coeff(j) * lai_ind(:) ) ), min_cover )
589
590	WHERE(fpc_nat(:,j).GT.fpc_max(:,j))
591
592	light_death(:,j) = MIN(un, un - fpc_max(:,j)/fpc_nat(:,j))
593
594	ENDWHERE
595
596	DO k=1,nparts
597
598	bm_to_litter(:,j,k) = bm_to_litter(:,j,k)+light_death(:,j)*biomass(:,j,k)
599	biomass(:,j,k) = biomass(:,j,k)-light_death(:,j)*biomass(:,j,k)
600
601	ENDDO
602	ind(:,j) = ind(:,j)-light_death(:,j)*ind(:,j)
603	! if (j==10) print ,'ind10bis=',ind(:,j),light_death(:,j)ind(:,j)
604	ENDIF
605	ENDDO
606
607	light_death(:,:) = light_death(:,:)/dt
608
609	ENDIF
610
611	!
612	! 5 history
613	!
614
615	CALL histwrite (hist_id_stomate, 'LIGHT_DEATH', itime, &
616	light_death, npts*nvm, horipft_index)
617
618	IF (bavard.GE.4) WRITE(numout,*) 'Leaving light'
619
620	END SUBROUTINE light
621
622	END MODULE lpj_light

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

Download in other formats: