source: branches/ORCHIDEE_2_2/ORCHIDEE/src_stomate/lpj_pftinout.f90 @ 6393

Last change on this file since 6393 was 4647, checked in by josefine.ghattas, 7 years ago

Ticket #392 : Change variable name veget_max into veget_cov_max to avoid confusion with sechiba variable veget_max. The output variable VEGET_MAX is also changed into VEGET_COV_MAX.

  • Property svn:keywords set to HeadURL Date Author Revision
File size: 22.1 KB
Line 
1! =================================================================================================================================
2! MODULE       : lpj_pftinout
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       Introduce and eliminate PFT's from pixel
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!! - Thonicke, K., S. Venevsky, et al. (2001), The role of fire disturbance
20!!        for global vegetation dynamics: coupling fire into a Dynamic Global
21!!        Vegetation Model, Global Ecology and Biogeography, 10, 661-677.\n
22!! - Haxeltine, A. and I. C. Prentice (1996), BIOME3: An equilibrium
23!!        terrestrial biosphere model based on ecophysiological constraints,
24!!        resource availability, and competition among plant functional types,
25!!        Global Biogeochemical Cycles, 10(4), 693-709.\n
26!! - Smith, B., I. C. Prentice, et al. (2001), Representation of vegetation
27!!        dynamics in the modelling of terrestrial ecosystems: comparing two
28!!        contrasting approaches within European climate space,
29!!        Global Ecology and Biogeography, 10, 621-637.\n
30!!
31!! SVN          :
32!! $HeadURL$
33!! $Date$
34!! $Revision$
35!! \n
36!_ ==============================================================================================================================
37
38MODULE lpj_pftinout
39
40  ! modules used:
41
42  USE ioipsl_para
43  USE stomate_data
44  USE pft_parameters
45  USE constantes
46  USE grid
47
48  IMPLICIT NONE
49
50  ! private & public routines
51
52  PRIVATE
53  PUBLIC pftinout,pftinout_clear
54
55  LOGICAL, SAVE                       :: firstcall_pftinout = .TRUE.                !! first call
56!$OMP THREADPRIVATE(firstcall_pftinout)
57
58CONTAINS
59
60
61!! ================================================================================================================================
62!! SUBROUTINE  : pftinout_clear
63!!
64!>\BRIEF       Set flag ::firstcall_pftinout to true and initialize the variables
65!_ ================================================================================================================================
66
67  SUBROUTINE pftinout_clear
68    firstcall_pftinout = .TRUE.
69  END SUBROUTINE pftinout_clear
70
71!! ================================================================================================================================
72!! SUBROUTINE  : pftinout
73!!
74!>\BRIEF       Introduce and eliminate PFT's from pixel
75!!
76!! DESCRIPTION**3 : Introduction and elimination of PFTs on the basis of climate condition.
77!! For natural and woody PFTs the foliage projected coverage is calculated as follows:
78!! \latexonly
79!!  \input{equation_lpj_pftinout.tex}
80!! \endlatexonly
81!! \n
82!! where FPC is foliage projective cover (::fpc_nat), CN crown area (::cn_ind,
83!! @tex $ m^{2} $ @endtex), IND number of individuals (::ind, @tex $ m^{-2} $ @endtex,
84!! FRAC total fraction occupied by natural vegetation (::fracnat),
85!! @tex $ LM_{rm max} $ @endtex maximum leaf mass in last year
86!! (::lm_lastyearmax, @tex $ g C m^{-2} $ @endtex), SLA specific leaf area (sla,
87!! @tex $ m^{2} (g C)^{-1} $ @endtex), and coff coefficient (::ext_coeff). ::ext_coeff
88!! describes a property of the canopy (i.e. law of Lambert-Beer) and is defined in **2
89!!
90!! The foliage projective cover feeds into the calculation of the space available for
91!! expansion of existing and dispersion of new PFTs within a gridbox. In turn, available
92!! space is use to calculate the number of individuals with a PFT.
93!!
94!! Saplings are introduced under the condition that winter temperature is
95!! moderate, plant age is older than 1.25, (and for some PFTs at least one adjacent grid
96!! box exists for expansion), new saplings are introduced for narural PFT. In the simulation of
97!! agricultural grassland, if target PFT does not exist in the gridbox, it is introduced
98!! regardless of climate condition. When a new PFT is introduced CO_2 is taken from the
99!! atmosphere to account for CO_2 present in the seed and required by the germinated seeds
100!! to establish a sapling. These initial phases in ontology are not accounted for. However,
101!! by taking this small amount of CO2 from the atmosphere, mass balance closure for C is
102!! preserved.
103!!
104!! PFTs are eliminated under the condition that they are no longer adapted to the critical
105!! temperatures in winter. When a PFT is eliminated its number of indiviuals is set to zero and
106!! the rest of the elimination process is taken care of in lpj_kill.f90.
107!!
108!! RECENT CHANGE(S) : None
109!!
110!! MAIN OUTPUT VARIABLE(S): :: avail_tree (space availability for trees, unitless),   
111!! :: avail_grass (space availability for grasses, unitless), :: biomass (biomass, \f$gC m^{-2}\f$)
112!! and :: ind (density of individuals, \f$m^{-2}\f$)   
113!!
114!! REFERENCE(S) : None
115!!
116!! FLOWCHART    :
117!! \latexonly
118!!   \includegraphics[scale = 0.6]{pftinout.png}
119!! \endlatexonly
120!! \n
121!_ ================================================================================================================================
122
123
124
125
126
127
128  SUBROUTINE pftinout (npts, dt, adapted, regenerate, &
129       neighbours, veget_cov_max, &
130       biomass, ind, cn_ind, age, leaf_frac, npp_longterm, lm_lastyearmax, senescence, &
131       PFTpresent, everywhere, when_growthinit, need_adjacent, RIP_time, &
132       co2_to_bm, &
133       avail_tree, avail_grass)
134   
135  !! 0. Variable and parameter declaration
136   
137    !! 0.1 Input variables
138   
139    INTEGER(i_std), INTENT(in)                                :: npts            !! Domain size - number of pixels (unitless)           
140    REAL(r_std), INTENT(in)                                   :: dt              !! Time step of vegetation dynamics for stomate
141                                                                                 !! (days)               
142    REAL(r_std), DIMENSION(npts,nvm), INTENT(in)              :: adapted         !! Winter not too cold (unitless)   
143    REAL(r_std), DIMENSION(npts,nvm), INTENT(in)              :: regenerate      !! Winter sufficiently cold (unitless) 
144    INTEGER(i_std), DIMENSION(npts,NbNeighb), INTENT(in)      :: neighbours      !! Indices of the 8 neighbours of each grid point
145                                                                                 !! (unitless); 1=North and then clockwise.
146    REAL(r_std), DIMENSION(npts,nvm), INTENT(inout)           :: veget_cov_max   !! "maximal" coverage fraction of a PFT (LAI ->
147                                                                                 !! infinity) on ground (unitless)       
148
149    !! 0.2 Output variables
150
151    REAL(r_std), DIMENSION(npts), INTENT(out)                 :: avail_tree      !! Space availability for trees (unitless)   
152    REAL(r_std), DIMENSION(npts), INTENT(out)                 :: avail_grass     !! Space availability for grasses (unitless)     
153
154
155    !! 0.3 Modified variables   
156
157    REAL(r_std), DIMENSION(npts,nvm,nparts,nelements), INTENT(inout) :: biomass  !! Biomass (gC m^{-2})     
158    REAL(r_std), DIMENSION(npts,nvm), INTENT(inout)           :: ind             !! Density of individuals (m^{-2})
159    REAL(r_std), DIMENSION(npts,nvm), INTENT(in)              :: cn_ind          !! Crown area of individuals (m^2)
160    REAL(r_std), DIMENSION(npts,nvm), INTENT(inout)           :: age             !! Mean age (years)           
161    REAL(r_std), DIMENSION(npts,nvm,nleafages), INTENT(inout) :: leaf_frac       !! Fraction of leaves in leaf age class (unitless) 
162    REAL(r_std), DIMENSION(npts,nvm), INTENT(inout)           :: npp_longterm    !! "long term" net primary productivity
163                                                                                 !! (gC m^{-2} year^{-1})   
164    REAL(r_std), DIMENSION(npts,nvm), INTENT(inout)           :: lm_lastyearmax  !! Last year's maximum leaf mass, for each PFT
165                                                                                 !! (gC m^{-2}) 
166    LOGICAL, DIMENSION(npts,nvm), INTENT(inout)               :: senescence      !! Plant senescent for deciduous trees; .FALSE.
167                                                                                 !! if PFT is introduced or killed     
168    LOGICAL, DIMENSION(npts,nvm), INTENT(inout)               :: PFTpresent      !! PFT exists   
169    REAL(r_std), DIMENSION(npts,nvm), INTENT(inout)           :: everywhere      !! is the PFT everywhere in the grid box or very
170                                                                                 !! localized (unitless) 
171    REAL(r_std), DIMENSION(npts,nvm), INTENT(inout)           :: when_growthinit !! how many days ago was the beginning of the
172                                                                                 !! growing season (days)
173    LOGICAL, DIMENSION(npts,nvm), INTENT(inout)               :: need_adjacent   !! in order for this PFT to be introduced, does it
174                                                                                 !! have to be present in an adjacent grid box?   
175    REAL(r_std), DIMENSION(npts,nvm), INTENT(inout)           :: RIP_time        !! How much time ago was the PFT eliminated for
176                                                                                 !! the last time (years)   
177    REAL(r_std), DIMENSION(npts,nvm), INTENT(inout)           :: co2_to_bm       !! biomass uptaken (gC m^{-2} day^{-1})       
178 
179    !! 0.4 Local variables
180
181    REAL(r_std), DIMENSION(npts)                              :: avail           !! availability           
182    INTEGER(i_std)                                            :: i,j,m           !! indices     
183    REAL(r_std), DIMENSION(npts)                              :: sumfrac_wood    !! total woody vegetation cover     
184    INTEGER(i_std), DIMENSION(npts)                           :: n_present       !! number of adjacent grid cells where PFT is
185                                                                                 !! ubiquitous       
186    LOGICAL, DIMENSION(npts)                                  :: can_introduce   !! we can introduce this PFT   
187    REAL(r_std), DIMENSION(npts)                              :: fracnat         !! no real need for dimension(ntps) except for
188                                                                                 !! vectorisation         
189!_ ================================================================================================================================
190
191    IF (printlev>=3) WRITE(numout,*) 'Entering pftinout'
192
193  !! 1. Messages   
194
195    IF ( firstcall_pftinout ) THEN
196
197       WRITE(numout,*) 'pftinout: Minimum space availability: ', min_avail
198
199       firstcall_pftinout = .FALSE.
200
201    ENDIF
202
203  !! 2. Total woody fpc and space avaibility on grid
204
205    ! Only natural part of the grid cell\n
206    ! S. Zaehle bug correction MERGE: need to subtract agricultural area!
207    ! fraction of agricultural surface
208
209    !! 2.1 only natural PFT
210    fracnat(:) = un
211    DO j = 2,nvm ! Loop over # PFTs
212       IF ( .NOT. natural(j) ) THEN
213          fracnat(:) = fracnat(:) - veget_cov_max(:,j)
214       ENDIF
215    ENDDO ! Loop over # PFTs
216
217    !! 2.2 Total woody fractional plant cover
218    sumfrac_wood(:) = zero
219    DO j = 2,nvm ! Loop over # PFTs
220       
221       ! S. Zaehle problem here: agriculture, not convinced that this representation of LPJ is correct
222       ! if agriculture is present, ind must be recalculated to correspond to the natural density...
223       ! since ind is per grid cell, can be achived by discounting for agricultura fraction
224       IF ( natural(j).AND.is_tree(j) ) THEN
225          WHERE(fracnat(:).GT.min_stomate) 
226                sumfrac_wood(:) = sumfrac_wood(:) + cn_ind(:,j) * ind(:,j) / fracnat(:) & 
227                     * ( un - exp( - lm_lastyearmax(:,j) * sla(j) * ext_coeff(j) ) )
228                !lai changed to lm_last
229          ENDWHERE
230       ENDIF
231    ENDDO ! Loop over # PFTs
232
233    !! 2.3 Space availability
234    avail_grass(:) = MAX( ( un - sumfrac_wood(:) ), min_avail )
235    avail_tree(:) = MAX( ( fpc_crit - sumfrac_wood(:) ), min_avail )
236
237  !! 3. Time since last elimination (y)
238
239    RIP_time = RIP_time + dt / one_year
240
241  !! 4. Agicultural PFTs
242
243    ! Agricultural PFTs are only present if they are prescribed
244    DO j = 2,nvm ! Loop over # PFTs
245       
246       IF ( .NOT. natural(j) ) THEN
247         
248          IF (printlev>=4) WRITE(numout,*) 'pftinout: Agricultural PFTs'
249         
250          !! 4.1 Agricultural trees
251          !      Agricultural trees are not treated for the moment
252          IF ( is_tree(j) ) THEN
253
254             CALL ipslerr_p(3,'pftinout','Agricultural trees not treated.','','')
255
256          !! 4.2 Initialization of agricultural grass lands
257          !      Initialize parameter values of prescribed agricultural PFTs
258          ELSE
259
260             DO i = 1, npts ! Loop over # pixels - domain size
261
262                IF ( ( veget_cov_max(i,j) .GT. min_stomate ) .AND. ( .NOT. PFTpresent(i,j) ) ) THEN
263
264                   ! prescribed, but not yet there.
265                   ind(i,j) = veget_cov_max(i,j)
266                   biomass(i,j,:,:) = bm_sapl(j,:,:) * ind(i,j) /veget_cov_max(i,j) 
267                   co2_to_bm(i,j) =  co2_to_bm(i,j) +SUM( biomass(i,j,:,icarbon) ) / dt
268                   PFTpresent(i,j) = .TRUE.
269                   everywhere(i,j) = un
270                   senescence(i,j) = .FALSE.
271                   age(i,j) = zero
272
273                ENDIF  ! prescribed, but PFT not yet present
274
275             ENDDO ! Loop over # pixels - domain size
276
277          ENDIF
278
279       ENDIF ! not natural
280
281    ENDDO ! Loop over # PFTs
282
283  !! 5 Eliminate PFTs
284
285    DO j = 2,nvm ! Loop over # PFTs
286
287       ! only for natural PFTs
288       IF ( natural(j) ) THEN
289
290          ! Number of individuals are set to zero in the condition of cold winter   
291          ! 'adapted_crit' critical value for being adapted = 1.-(1./euler); see 'stomate_constants.f90'
292          WHERE (  PFTpresent(:,j) .AND. ( adapted(:,j) .LT. adapted_crit ) )
293
294             ! PFT there, but not adapted any more (ex: winter too cold): kill
295             ! set number of individuals to zero - rest will be done in lpj_kill
296             ind(:,j) = zero
297
298          ENDWHERE
299
300       ENDIF ! natural
301
302    ENDDO ! Loop over # PFTs
303
304  !! 6. Introduce PFTs
305
306    DO j = 2,nvm ! Loop over # PFTs
307
308       IF ( natural(j) ) THEN
309
310          ! space availability for this PFT
311          IF ( is_tree(j) ) THEN
312             avail(:) = avail_tree(:)
313          ELSE
314             avail(:) = avail_grass(:)
315          ENDIF
316         
317          !! 6.1 Check if PFT not present but (adapted and regenerative)     
318          can_introduce(:) = .FALSE.
319         
320          IF ( NbNeighb /= 8 ) THEN
321             CALL ipslerr(3, "pftinout", "This routine needs to be adapted to non rectengular grids", &
322                  &           "Talk to Jan Polcher", " ")
323          ENDIF
324
325          DO i = 1, npts ! Loop over # pixels - domain size
326
327             IF ( .NOT. PFTpresent(i,j) .AND. &
328                  ( adapted(i,j) .GT. adapted_crit ) .AND. &
329                  ( regenerate(i,j) .GT. regenerate_crit )  ) THEN
330
331                ! Seed are available nearby
332                IF ( need_adjacent(i,j) ) THEN
333
334                   !! 6.1.1 Climate allows introduction of the PFT but dispersion requires the
335                   !        presence of seeds. Seed are considered available if at least one
336                   !        neighbouring pixel is entirely invaded by the PFT. If that condition is
337                   !        satisfied, the PFT can establish in the new pixel.
338                   ! Count number of totally invaded neighbours.
339                   ! no loop so that it can vectorize
340                   n_present(i) = 0
341                   IF ( neighbours(i,1) .GT. 0 ) THEN
342                      IF ( everywhere(neighbours(i,1),j) .GE. un-min_stomate ) THEN
343                         n_present(i) = n_present(i)+1
344                      ENDIF
345                   ENDIF
346                   IF ( neighbours(i,3) .GT. 0 ) THEN
347                      IF ( everywhere(neighbours(i,3),j) .GE. un-min_stomate ) THEN
348                         n_present(i) = n_present(i)+1
349                      ENDIF
350                   ENDIF
351                   IF ( neighbours(i,5) .GT. 0 ) THEN
352                      IF ( everywhere(neighbours(i,5),j) .GE. un-min_stomate ) THEN
353                         n_present(i) = n_present(i)+1
354                      ENDIF
355                   ENDIF
356                   IF ( neighbours(i,7) .GT. 0 ) THEN
357                      IF ( everywhere(neighbours(i,7),j) .GE. un-min_stomate ) THEN
358                         n_present(i) = n_present(i)+1
359                      ENDIF
360                   ENDIF
361
362                   IF ( n_present(i) .GT. 0 ) THEN
363
364                      ! PFT is ubiquitous in at least one adjacent grid box
365                      can_introduce(i) = .TRUE.
366
367                   ENDIF
368
369                ELSE
370
371                   !! 6.1.2 No seed (trees) required for dispersion
372                   !        The PFT can establish without the presence of seed trees in
373                   !        neighbouring pixels.
374                   can_introduce(i) = .TRUE.
375
376                ENDIF ! do we have to look at the neighbours?
377
378             ENDIF ! we'd like to introduce the PFT
379
380          ENDDO ! Loop over # pixels - domain size
381
382          !! 6.2 Has the PFT been eliminated lately?
383          !      Additional test whether the PFT has been eliminated lately, i.e.
384          !      less than 1.25 years ago. Do not only take full years as success of
385          !      introduction, as introduction might depend on season.
386          WHERE ( RIP_time(:,j) .LT. RIP_time_min )
387
388             ! PFT was eliminated lately - cannot reintroduce
389             can_introduce(:) = .FALSE.
390
391          ENDWHERE
392
393          !! 6.3 Introduce that PFT where possible
394          !      "can_introduce" means that it either exists in neighbouring grid boxes
395          !      or that we do not look at neighbours, that it has not been eliminated
396          !      lately, and, of course, that the climate is good for that PFT.
397          WHERE ( can_introduce(:) )
398             
399             PFTpresent(:,j) = .TRUE.
400             
401             senescence(:,j) = .FALSE.
402             
403             ! introduce at least a few saplings, even if canopy is closed
404             ! initial density of individuals (ind_0) = 0.02, see 'stomate_constant.f90'
405             ind(:,j) = ind_0 * (dt/one_year) * avail(:)
406             
407             WHERE(veget_cov_max(:,j) .GT. min_stomate)
408                biomass(:,j,ileaf,icarbon) = bm_sapl(j,ileaf,icarbon) * ind(:,j) /veget_cov_max(:,j)
409                biomass(:,j,isapabove,icarbon) = bm_sapl(j,isapabove,icarbon) * ind(:,j) /veget_cov_max(:,j)
410                biomass(:,j,isapbelow,icarbon) = bm_sapl(j,isapbelow,icarbon) * ind(:,j)/veget_cov_max(:,j)
411                biomass(:,j,iheartabove,icarbon) = bm_sapl(j,iheartabove,icarbon) * ind(:,j)/veget_cov_max(:,j)
412                biomass(:,j,iheartbelow,icarbon) = bm_sapl(j,iheartbelow,icarbon) * ind(:,j)/veget_cov_max(:,j)
413                biomass(:,j,iroot,icarbon) = bm_sapl(j,iroot,icarbon) * ind(:,j)/veget_cov_max(:,j)
414                biomass(:,j,ifruit,icarbon) = bm_sapl(j,ifruit,icarbon) * ind(:,j)/veget_cov_max(:,j)
415                biomass(:,j,icarbres,icarbon) = bm_sapl(j,icarbres,icarbon) * ind(:,j)/veget_cov_max(:,j)
416             ELSEWHERE             
417                biomass(:,j,ileaf,icarbon) = bm_sapl(j,ileaf,icarbon) * ind(:,j)
418                biomass(:,j,isapabove,icarbon) = bm_sapl(j,isapabove,icarbon) * ind(:,j)
419                biomass(:,j,isapbelow,icarbon) = bm_sapl(j,isapbelow,icarbon) * ind(:,j)
420                biomass(:,j,iheartabove,icarbon) = bm_sapl(j,iheartabove,icarbon) * ind(:,j)
421                biomass(:,j,iheartbelow,icarbon) = bm_sapl(j,iheartbelow,icarbon) * ind(:,j)
422                biomass(:,j,iroot,icarbon) = bm_sapl(j,iroot,icarbon) * ind(:,j)
423                biomass(:,j,ifruit,icarbon) = bm_sapl(j,ifruit,icarbon) * ind(:,j)
424                biomass(:,j,icarbres,icarbon) = bm_sapl(j,icarbres,icarbon) * ind(:,j)
425             END WHERE
426           
427             co2_to_bm(:,j) = &
428                  co2_to_bm(:,j) +  &
429                  ( biomass(:,j,ileaf,icarbon) + biomass(:,j,isapabove,icarbon) + &
430                  biomass(:,j,isapbelow,icarbon) + biomass(:,j,iheartabove,icarbon) + &
431                  biomass(:,j,iheartbelow,icarbon) + biomass(:,j,iroot,icarbon) + &
432                  biomass(:,j,ifruit,icarbon) + biomass(:,j,icarbres,icarbon) )/dt
433
434             when_growthinit(:,j) = large_value
435
436             age(:,j) = zero
437
438             ! all leaves are young
439             leaf_frac(:,j,1) = un
440
441             ! non-zero "long term" npp and last year's leaf mass for saplings -
442             ! so they won't be killed off by gap or kill
443             npp_longterm(:,j) = npp_longterm_init
444
445             lm_lastyearmax(:,j) = bm_sapl(j,ileaf,icarbon) * ind(:,j)
446
447          ENDWHERE    ! we can introduce the PFT
448
449          !! 6.4 Expansion of the PFT within the grid box
450          !      PFT expansion/dispersion to a new grid box should not be confused with
451          !      expansion in areal coverage
452          IF ( treat_expansion ) THEN
453
454             WHERE ( can_introduce(:) )
455
456                ! low value at the beginning
457                everywhere(:,j) = everywhere_init
458             ENDWHERE
459
460          ELSE
461
462             ! expansion is not treated
463             WHERE ( can_introduce(:) )
464                everywhere(:,j) = un
465             ENDWHERE
466
467          ENDIF ! treat expansion
468
469       ENDIF ! only natural PFTs
470
471    ENDDO ! Loop over # PFTs
472
473  !! 7. If a PFT has been present once in a grid box, we suppose that it will survive
474
475    !   If a PFT has been present once in a grid box, we suppose that it will survive
476    !   in isolated places (e.g., an oasis) within that grid box, even if it gets
477    !   officially eliminated from it later. That means that if climate becomes favorable
478    !   again, it will not need to get seeds from adjacent grid cells.
479    WHERE ( PFTpresent )
480       need_adjacent = .FALSE.
481    ENDWHERE
482
483    IF (printlev>=4) WRITE(numout,*) 'Leaving pftinout'
484
485  END SUBROUTINE pftinout
486
487END MODULE lpj_pftinout
Note: See TracBrowser for help on using the repository browser.