source: tags/ORCHIDEE_1_9_5/ORCHIDEE/src_stomate/lpj_light.f90 @ 8

Last change on this file since 8 was 8, checked in by orchidee, 14 years ago

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