Context Navigation

source: tags/ORCHIDEE_1_9_5_2/ORCHIDEE/src_stomate/lpj_light.f90 @ 2061

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

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

Download in other formats: