source: branches/publications/ORCHIDEE-GMv3.2/ORCHIDEE/src_stomate/stomate_glcc_Bioe1.f90 @ 6940

Last change on this file since 6940 was 6940, checked in by jinfeng.chang, 4 years ago

add missing files for ORCHIDEE-GMv3.2

File size: 49.0 KB
Line 
1! Remove definition of type_conversion
2
3! =================================================================================================================================
4! MODULE       : stomate_glcc_bioe1
5!
6! CONTACT      : orchidee-help _at_ ipsl.jussieu.fr
7!
8! LICENCE      : IPSL (2006)
9! This software is governed by the CeCILL licence see ORCHIDEE/ORCHIDEE_CeCILL.LIC
10!
11!>\BRIEF       This module implements gross land use change with age classes.
12!!
13!!\n DESCRIPTION: None
14!!
15!! RECENT CHANGE(S): Including permafrost carbon
16!!
17!! REFERENCE(S) : None
18!!
19!! SVN          :
20!! $HeadURL: svn://forge.ipsl.jussieu.fr/orchidee/perso/albert.jornet/ORCHIDEE-MICT/src_stomate/stomate_lcchange.f90 $
21!! $Date: 2015-07-30 15:38:45 +0200 (Thu, 30 Jul 2015) $
22!! $Revision: 2847 $
23!! \n
24!_ ================================================================================================================================
25
26
27MODULE stomate_glcc_bioe1
28
29  ! modules used:
30 
31  USE ioipsl_para
32  USE stomate_data
33  USE pft_parameters
34  USE constantes
35  USE constantes_soil_var
36  USE stomate_gluc_common
37  USE stomate_gluc_constants
38  USE stomate_glcchange_MulAgeC
39  USE xios_orchidee
40
41  IMPLICIT NONE
42 
43  PRIVATE
44  PUBLIC glcc_bioe1_firstday, glcc_bioe1 
45 
46CONTAINS
47
48! ================================================================================================================================
49!! SUBROUTINE   gross_lcchange
50!!
51!>\BRIEF       : Apply gross land cover change.
52!!
53!>\DESCRIPTION 
54!_ ================================================================================================================================
55  SUBROUTINE glcc_bioe1 (npts, dt_days, newvegfrac,  &
56               glccSecondShift,&
57               def_fuel_1hr_remain, def_fuel_10hr_remain,        &
58               def_fuel_100hr_remain, def_fuel_1000hr_remain,    &
59               deforest_litter_remain, deforest_biomass_remain,  &
60               convflux, cflux_prod10, cflux_prod100,                  &
61               glcc_pft, glcc_pftmtc,          &
62               veget_max, prod10, prod100, flux10, flux100,            &
63               PFTpresent, senescence, moiavail_month, moiavail_week,  &
64               gpp_week, ngd_minus5, resp_maint, resp_growth,          & 
65               resp_hetero, npp_daily, when_growthinit, npp_longterm,  &
66               ind, lm_lastyearmax, everywhere, age,                   &
67               co2_to_bm, gpp_daily, co2_fire,                         &
68               time_hum_min, gdd_midwinter, gdd_from_growthinit,       &
69               gdd_m5_dormance, ncd_dormance,                          &
70               lignin_struc, carbon, leaf_frac,                        &
71               deepC_a, deepC_s, deepC_p,                              &
72               leaf_age, bm_to_litter, biomass, litter,                &
73               fuel_1hr, fuel_10hr, fuel_100hr, fuel_1000hr)
74 
75    IMPLICIT NONE
76
77    !! 0.1 Input variables
78
79    INTEGER, INTENT(in)                                  :: npts             !! Domain size - number of pixels (unitless)
80    REAL(r_std), INTENT(in)                              :: dt_days          !! Time step of vegetation dynamics for stomate
81    REAL(r_std), DIMENSION (npts,12),INTENT(in)        :: glccSecondShift     !! the land-cover-change (LCC) matrix in case a gross LCC is
82                                                                              !! used.
83    REAL(r_std), DIMENSION (npts,nvmap),INTENT(in)       :: newvegfrac             !!
84                                                                             !!
85
86    REAL(r_std), DIMENSION(npts,nvm,nlitt,nelements), INTENT(in)                 :: def_fuel_1hr_remain
87    REAL(r_std), DIMENSION(npts,nvm,nlitt,nelements), INTENT(in)                 :: def_fuel_10hr_remain
88    REAL(r_std), DIMENSION(npts,nvm,nlitt,nelements), INTENT(in)                 :: def_fuel_100hr_remain
89    REAL(r_std), DIMENSION(npts,nvm,nlitt,nelements), INTENT(in)                 :: def_fuel_1000hr_remain
90    REAL(r_std), DIMENSION(npts,nlitt,nvm,nlevs,nelements), INTENT(in) :: deforest_litter_remain   !! Vegetmax-weighted remaining litter on the ground for
91                                                                                                      !! deforestation region.
92    REAL(r_std), DIMENSION(npts,nvm,nparts,nelements), INTENT(in)      :: deforest_biomass_remain  !! Vegetmax-weighted remaining biomass on the ground for
93                                                                                                      !! deforestation region.
94
95
96    !! 0.2 Output variables
97    REAL(r_std), DIMENSION(npts,nwp), INTENT(out)            :: convflux         !! release during first year following land cover
98                                                                             !! change
99    REAL(r_std), DIMENSION(npts,nwp), INTENT(out)            :: cflux_prod10     !! total annual release from the 10 year-turnover
100                                                                             !! pool @tex ($gC m^{-2}$) @endtex
101    REAL(r_std), DIMENSION(npts,nwp), INTENT(out)            :: cflux_prod100    !! total annual release from the 100 year-
102    REAL(r_std), DIMENSION(npts,nvm), INTENT(inout)      :: glcc_pft         !! Loss of fraction in each PFT
103    REAL(r_std), DIMENSION(npts,nvm,nvmap), INTENT(inout):: glcc_pftmtc      !! a temporary variable to hold the fractions each PFT is going to lose
104                                                                             !! i.e., the contribution of each PFT to the youngest age-class of MTC
105
106    !! 0.3 Modified variables
107    REAL(r_std), DIMENSION(npts,nvm), INTENT(inout)      :: veget_max        !! "maximal" coverage fraction of a PFT (LAI ->
108                                                                             !! infinity) on ground (unitless)
109    REAL(r_std), DIMENSION(npts,0:10,nwp), INTENT(inout)     :: prod10           !! products remaining in the 10 year-turnover
110                                                                             !! pool after the annual release for each
111                                                                             !! compartment (10 + 1 : input from year of land
112                                                                             !! cover change)
113    REAL(r_std), DIMENSION(npts,0:100,nwp), INTENT(inout)    :: prod100          !! products remaining in the 100 year-turnover
114                                                                             !! pool after the annual release for each
115                                                                             !! compartment (100 + 1 : input from year of land
116                                                                             !! cover change)
117    REAL(r_std), DIMENSION(npts,10,nwp), INTENT(inout)       :: flux10           !! annual release from the 10/100 year-turnover
118                                                                             !! pool compartments
119    REAL(r_std), DIMENSION(npts,100,nwp), INTENT(inout)      :: flux100          !! annual release from the 10/100 year-turnover
120                                                                             !! pool compartments
121    LOGICAL, DIMENSION(:,:), INTENT(inout)               :: PFTpresent       !! Tab indicating which PFTs are present in
122                                                                             !! each pixel
123    LOGICAL, DIMENSION(:,:), INTENT(inout)               :: senescence       !! Flag for setting senescence stage (only
124                                                                             !! for deciduous trees)
125    REAL(r_std), DIMENSION(:,:), INTENT(inout)           :: moiavail_month   !! "Monthly" moisture availability (0 to 1,
126                                                                             !! unitless)
127    REAL(r_std), DIMENSION(:,:), INTENT(inout)           :: moiavail_week    !! "Weekly" moisture availability
128                                                                             !! (0 to 1, unitless)
129    REAL(r_std), DIMENSION(:,:), INTENT(inout)           :: gpp_week         !! Mean weekly gross primary productivity
130                                                                             !! @tex $(gC m^{-2} day^{-1})$ @endtex
131    REAL(r_std), DIMENSION(:,:), INTENT(inout)           :: ngd_minus5       !! Number of growing days (days), threshold
132                                                                             !! -5 deg C (for phenology)   
133    REAL(r_std), DIMENSION(:,:), INTENT(inout)           :: resp_maint       !! Maintenance respiration 
134                                                                             !! @tex $(gC m^{-2} dtslow^{-1})$ @endtex
135    REAL(r_std), DIMENSION(:,:), INTENT(inout)           :: resp_growth      !! Growth respiration 
136                                                                             !! @tex $(gC m^{-2} dtslow^{-1})$ @endtex
137    REAL(r_std), DIMENSION(:,:), INTENT(inout)           :: resp_hetero      !! Heterotrophic respiration 
138                                                                             !! @tex $(gC m^{-2} dtslow^{-1})$ @endtex
139    REAL(r_std), DIMENSION(:,:), INTENT(inout)           :: npp_daily        !! Net primary productivity
140                                                                             !! @tex $(gC m^{-2} dtslow^{-1})$ @endtex
141    REAL(r_std), DIMENSION(:,:), INTENT(inout)           :: when_growthinit  !! How many days ago was the beginning of
142                                                                             !! the growing season (days)
143    REAL(r_std), DIMENSION(:,:), INTENT(inout)           :: npp_longterm     !! "Long term" mean yearly primary productivity
144    REAL(r_std), DIMENSION(:,:), INTENT(inout)           :: ind              !! Number of individuals at the stand level
145                                                                             !! @tex $(m^{-2})$ @endtex
146    REAL(r_std), DIMENSION(:,:), INTENT(inout)           :: lm_lastyearmax   !! last year's maximum leaf mass for each PFT
147                                                                             !! @tex ($gC m^{-2}$) @endtex
148    REAL(r_std), DIMENSION(:,:), INTENT(inout)           :: everywhere       !! is the PFT everywhere in the grid box or
149                                                                             !! very localized (after its introduction) (?)
150    REAL(r_std), DIMENSION(:,:), INTENT(inout)           :: age              !! mean age (years)
151    REAL(r_std), DIMENSION(:,:), INTENT(inout)           :: co2_to_bm        !! CO2 taken from the atmosphere to get C to create 
152                                                                             !! the seedlings @tex (gC.m^{-2}dt^{-1})$ @endtex
153    REAL(r_std), DIMENSION(:,:), INTENT(inout)           :: gpp_daily        !! Daily gross primary productivity 
154                                                                             !! @tex $(gC m^{-2} dtslow^{-1})$ @endtex
155    REAL(r_std), DIMENSION(:,:), INTENT(inout)           :: co2_fire         !! Fire carbon emissions
156                                                                             !! @tex $(gC m^{-2} dtslow^{-1})$ @endtex
157    REAL(r_std), DIMENSION(:,:), INTENT(inout)           :: time_hum_min     !! Time elapsed since strongest moisture
158                                                                             !! availability (days)
159    REAL(r_std), DIMENSION(:,:), INTENT(inout)           :: gdd_midwinter    !! Growing degree days (K), since midwinter
160                                                                             !! (for phenology) - this is written to the
161                                                                             !!  history files
162    REAL(r_std), DIMENSION(:,:), INTENT(inout)           :: gdd_from_growthinit !! growing degree days, since growthinit
163                                                                             !! for crops
164    REAL(r_std), DIMENSION(:,:), INTENT(inout)           :: gdd_m5_dormance  !! Growing degree days (K), threshold -5 deg
165                                                                             !! C (for phenology)
166    REAL(r_std), DIMENSION(:,:), INTENT(inout)           :: ncd_dormance     !! Number of chilling days (days), since
167                                                                             !! leaves were lost (for phenology)
168    REAL(r_std), DIMENSION(:,:,:), INTENT(inout)         :: lignin_struc     !! ratio Lignine/Carbon in structural litter,
169                                                                             !! above and below ground
170    REAL(r_std), DIMENSION(:,:,:), INTENT(inout)         :: carbon           !! carbon pool: active, slow, or passive
171                                                                             !! @tex ($gC m^{-2}$) @endtex
172    REAL(r_std), DIMENSION(:,:,:), INTENT(inout)         :: deepC_a          !! Permafrost soil carbon (g/m**3) active
173    REAL(r_std), DIMENSION(:,:,:), INTENT(inout)         :: deepC_s          !! Permafrost soil carbon (g/m**3) slow
174    REAL(r_std), DIMENSION(:,:,:), INTENT(inout)         :: deepC_p          !! Permafrost soil carbon (g/m**3) passive
175    REAL(r_std), DIMENSION(:,:,:), INTENT(inout)         :: leaf_frac        !! fraction of leaves in leaf age class (unitless;0-1)
176    REAL(r_std), DIMENSION(:,:,:), INTENT(inout)         :: leaf_age         !! Leaf age (days)
177    REAL(r_std), DIMENSION(:,:,:,:), INTENT(inout)       :: bm_to_litter     !! Transfer of biomass to litter
178                                                                             !! @tex $(gC m^{-2} dtslow^{-1})$ @endtex
179    REAL(r_std), DIMENSION(:,:,:,:), INTENT(inout)       :: biomass          !! Stand level biomass @tex $(gC.m^{-2})$ @endtex
180    REAL(r_std), DIMENSION(:,:,:,:,:), INTENT(inout)     :: litter           !! metabolic and structural litter, above and
181                                                                             !! below ground @tex ($gC m^{-2}$) @endtex
182    REAL(r_std), DIMENSION(:,:,:,:), INTENT(inout)       :: fuel_1hr
183    REAL(r_std), DIMENSION(:,:,:,:), INTENT(inout)       :: fuel_10hr
184    REAL(r_std), DIMENSION(:,:,:,:), INTENT(inout)       :: fuel_100hr
185    REAL(r_std), DIMENSION(:,:,:,:), INTENT(inout)       :: fuel_1000hr
186
187    !! 0.4 Local variables
188    REAL(r_std), DIMENSION(nvmap,nparts,nelements)             :: bm_to_litter_pro !! conversion of biomass to litter
189                                                                             !! @tex ($gC m^{-2} day^{-1}$) @endtex
190    REAL(r_std), DIMENSION(nvmap,nparts,nelements)             :: biomass_pro      !! biomass @tex ($gC m^{-2}$) @endtex
191    REAL(r_std), DIMENSION(nvmap)                        :: veget_max_pro    !! "maximal" coverage fraction of a PFT (LAI ->
192                                                                             !! infinity) on ground (unitless)
193    REAL(r_std), DIMENSION(nvmap,ncarb)                        :: carbon_pro       !! carbon pool: active, slow, or passive
194                                                                             !! @tex ($gC m^{-2}$) @endtex
195    REAL(r_std), DIMENSION(nvmap,ndeep)                        :: deepC_a_pro      !! Permafrost carbon pool: active, slow, or passive
196                                                                             !! @tex ($gC m^{-3}$) @endtex
197    REAL(r_std), DIMENSION(nvmap,ndeep)                        :: deepC_s_pro      !! Permafrost carbon pool: active, slow, or passive
198                                                                             !! @tex ($gC m^{-3}$) @endtex
199    REAL(r_std), DIMENSION(nvmap,ndeep)                        :: deepC_p_pro      !! Permafrost carbon pool: active, slow, or passive
200                                                                             !! @tex ($gC m^{-3}$) @endtex
201    REAL(r_std), DIMENSION(nvmap,nlitt,nlevs,nelements)        :: litter_pro       !! metabolic and structural litter, above and
202                                                                             !! below ground @tex ($gC m^{-2}$) @endtex
203    REAL(r_std), DIMENSION(nvmap,nlitt,nelements)              :: fuel_1hr_pro
204    REAL(r_std), DIMENSION(nvmap,nlitt,nelements)              :: fuel_10hr_pro
205    REAL(r_std), DIMENSION(nvmap,nlitt,nelements)              :: fuel_100hr_pro
206    REAL(r_std), DIMENSION(nvmap,nlitt,nelements)              :: fuel_1000hr_pro
207    REAL(r_std), DIMENSION(nvmap,nlevs)                        :: lignin_struc_pro !! ratio Lignine/Carbon in structural litter
208                                                                             !! above and below ground
209    REAL(r_std), DIMENSION(nvmap,nleafages)                    :: leaf_frac_pro    !! fraction of leaves in leaf age class
210    REAL(r_std), DIMENSION(nvmap,nleafages)                    :: leaf_age_pro     !! fraction of leaves in leaf age class
211    LOGICAL, DIMENSION(nvmap)                :: PFTpresent_pro, senescence_pro                 !! Is pft there (unitless)
212    REAL(r_std), DIMENSION(nvmap)            :: ind_pro, age_pro, lm_lastyearmax_pro, npp_longterm_pro
213    REAL(r_std), DIMENSION(nvmap)            :: everywhere_pro
214    REAL(r_std), DIMENSION(nvmap)            :: gpp_daily_pro, npp_daily_pro, co2_to_bm_pro
215    REAL(r_std), DIMENSION(nvmap)            :: resp_maint_pro, resp_growth_pro
216    REAL(r_std), DIMENSION(nvmap)            :: resp_hetero_pro, co2_fire_pro
217 
218    INTEGER                :: ipts,ivm,ivma,l,m,ipft_young_agec
219    CHARACTER(LEN=10)      :: part_str                               !! string suffix indicating an index
220
221    REAL(r_std), DIMENSION(npts,nvmap)       :: glcc_mtc             !! Increase in fraction of each MTC in its youngest age-class
222    REAL(r_std), DIMENSION(npts,nvm)         :: glccReal_tmp         !! A temporary variable
223
224    WRITE(numout,*) 'Entering glcc_MulAgeC'
225    glccReal_tmp(:,:) = zero
226
227    CALL glcc_bioe1_firstday(npts,veget_max,newvegfrac,   &
228                          glccSecondShift,glcc_pft,glcc_pftmtc)
229
230    glcc_mtc(:,:) = SUM(glcc_pftmtc,DIM=2)
231
232    DO ipts=1,npts
233
234      !! Initialize the _pro variables
235      bm_to_litter_pro(:,:,:)=zero                                               
236      biomass_pro(:,:,:)=zero
237      veget_max_pro(:)=zero                                                       
238      carbon_pro(:,:)=zero                                                       
239      deepC_a_pro(:,:)=zero                                                       
240      deepC_s_pro(:,:)=zero                                                       
241      deepC_p_pro(:,:)=zero                                                       
242      litter_pro(:,:,:,:)=zero                                                   
243      fuel_1hr_pro(:,:,:)=zero                                                   
244      fuel_10hr_pro(:,:,:)=zero                                                   
245      fuel_100hr_pro(:,:,:)=zero                                                 
246      fuel_1000hr_pro(:,:,:)=zero                                                 
247      lignin_struc_pro(:,:)=zero                                                 
248
249      leaf_frac_pro = zero
250      leaf_age_pro = zero
251      PFTpresent_pro(:) = .FALSE.
252      senescence_pro(:) = .TRUE.
253      ind_pro = zero
254      age_pro = zero
255      lm_lastyearmax_pro = zero
256      npp_longterm_pro = zero
257      everywhere_pro = zero
258     
259      gpp_daily_pro=zero                                                       
260      npp_daily_pro=zero                                                       
261      co2_to_bm_pro=zero                                                       
262      resp_maint_pro=zero                                                     
263      resp_growth_pro=zero                                                     
264      resp_hetero_pro=zero                                                     
265      co2_fire_pro=zero                                                       
266
267      ! Note that we assume people don't intentionally change baresoil to
268      ! vegetated land.
269      DO ivma = 2,nvmap
270
271        ! here we set glcc_mtc(ipts,ivma) > min_stomate as a condition,
272        ! this is necessary because later on in the subroutine of
273        ! add_incoming_proxy_pft we have to merge the newly established
274        ! youngest proxy with potentially exisiting youngest receiving MTC,
275        ! thus have to devide a new fraction of (frac_proxy + frac_exist),
276        ! but in case frac_exist = zero, we risk deviding by a very small value
277        ! of frac_proxy and thus we want it to be bigger than min_stomate.
278        IF ( glcc_mtc(ipts,ivma) .GT. min_stomate ) THEN
279
280          ! 1. we accumulate the scalar variables that will be inherited.
281          !    note that in the subroutine of `collect_legacy_pft`, all
282          !    zero transitions will be automatically skipped.
283          CALL collect_legacy_pft(npts, ipts, ivma, glcc_pftmtc,    &
284                  biomass, bm_to_litter, carbon, litter,            &
285                  deepC_a, deepC_s, deepC_p,                        &
286                  fuel_1hr, fuel_10hr, fuel_100hr, fuel_1000hr,     &
287                  lignin_struc, co2_to_bm, gpp_daily, npp_daily,    &
288                  resp_maint, resp_growth, resp_hetero, co2_fire,   &
289                  def_fuel_1hr_remain, def_fuel_10hr_remain,        &
290                  def_fuel_100hr_remain, def_fuel_1000hr_remain,    &
291                  deforest_litter_remain, deforest_biomass_remain,  &
292                  veget_max_pro(ivma), carbon_pro(ivma,:),          &
293                  lignin_struc_pro(ivma,:), litter_pro(ivma,:,:,:), &
294                  deepC_a_pro(ivma,:), deepC_s_pro(ivma,:), deepC_p_pro(ivma,:), &
295                  fuel_1hr_pro(ivma,:,:), fuel_10hr_pro(ivma,:,:),  &
296                  fuel_100hr_pro(ivma,:,:), fuel_1000hr_pro(ivma,:,:), &
297                  bm_to_litter_pro(ivma,:,:), co2_to_bm_pro(ivma),  &
298                  gpp_daily_pro(ivma), npp_daily_pro(ivma),         &
299                  resp_maint_pro(ivma), resp_growth_pro(ivma),      &
300                  resp_hetero_pro(ivma), co2_fire_pro(ivma),        &
301                  convflux,prod10,prod100)
302
303          !++TEMP++
304          ! Here we substract the outgoing fraction from the source PFT.
305          ! If a too small fraction remains in this source PFT, then it is
306          ! exhausted, we empty it. The subroutine 'empty_pft' might be
307          ! combined with 'collect_legacy_pft', but now we just put it here.
308          DO ivm = 1,nvm
309            IF( glcc_pftmtc(ipts,ivm,ivma)>zero ) THEN
310              veget_max(ipts,ivm) = veget_max(ipts,ivm)-glcc_pftmtc(ipts,ivm,ivma)
311              IF ( veget_max(ipts,ivm)<min_stomate ) THEN
312                CALL empty_pft(ipts, ivm, veget_max, biomass, ind,       &
313                       carbon, litter, lignin_struc, bm_to_litter,       &
314                       deepC_a, deepC_s, deepC_p,                        &
315                       fuel_1hr, fuel_10hr, fuel_100hr, fuel_1000hr,     &
316                       gpp_daily, npp_daily, gpp_week, npp_longterm,     &
317                       co2_to_bm, resp_maint, resp_growth, resp_hetero,  &
318                       lm_lastyearmax, leaf_frac, leaf_age, age,         &
319                       everywhere, PFTpresent, when_growthinit,          &
320                       senescence, gdd_from_growthinit, gdd_midwinter,   &
321                       time_hum_min, gdd_m5_dormance, ncd_dormance,      &
322                       moiavail_month, moiavail_week, ngd_minus5)
323              ENDIF
324            ENDIF
325          ENDDO
326
327        ENDIF !IF ( glcc_mtc(ipts,ivma) .GT. min_stomate )
328      ENDDO ! (DO ivma = 2,nvmap)
329
330      ! We can only establish new youngest proxy and add it to the
331      ! existing youngest-age PFT after all wood product extraction is done, to
332      ! avoid the dilution of extractable biomass by the young proxy
333      ! and ensure consistency. Therefore now we have to loop again
334      ! over nvmap.
335      DO ivma = 2,nvmap
336        !IF ( glcc_mtc(ipts,ivma) .GT. min_stomate ) THEN
337
338          ipft_young_agec = start_index(ivma)
339
340          ! 2. we establish a proxy PFT with the fraction of veget_max_pro,
341          !    which is going to be either merged with existing target
342          !    `ipft_young_agec` PFT, or fill the place if no existing target PFT
343          !    exits.
344          CALL initialize_proxy_pft(ipts,ipft_young_agec,veget_max_pro(ivma), &
345                 biomass_pro(ivma,:,:), co2_to_bm_pro(ivma), ind_pro(ivma),   &
346                 age_pro(ivma),                                               & 
347                 senescence_pro(ivma), PFTpresent_pro(ivma),                  &
348                 lm_lastyearmax_pro(ivma), everywhere_pro(ivma),              &
349                 npp_longterm_pro(ivma),                                      &
350                 leaf_frac_pro(ivma,:),leaf_age_pro(ivma,:))
351
352          CALL sap_take (ipts,ivma,veget_max,biomass_pro(ivma,:,:), &
353                         biomass,co2_to_bm_pro(ivma))
354
355          ! 3. we merge the newly initiazlized proxy PFT into existing one
356          !    or use it to fill an empty PFT slot.
357          CALL add_incoming_proxy_pft(npts, ipts, ipft_young_agec, veget_max_pro(ivma),&
358                 carbon_pro(ivma,:), litter_pro(ivma,:,:,:), lignin_struc_pro(ivma,:), &
359                 bm_to_litter_pro(ivma,:,:),    &
360                 deepC_a_pro(ivma,:), deepC_s_pro(ivma,:), deepC_p_pro(ivma,:), &
361                 fuel_1hr_pro(ivma,:,:), fuel_10hr_pro(ivma,:,:),               &
362                 fuel_100hr_pro(ivma,:,:), fuel_1000hr_pro(ivma,:,:),           &
363                 biomass_pro(ivma,:,:), co2_to_bm_pro(ivma),                    &
364                 npp_longterm_pro(ivma), ind_pro(ivma),                         &
365                 lm_lastyearmax_pro(ivma), age_pro(ivma), everywhere_pro(ivma), & 
366                 leaf_frac_pro(ivma,:), leaf_age_pro(ivma,:),                   &
367                 PFTpresent_pro(ivma), senescence_pro(ivma),                &
368                 gpp_daily_pro(ivma), npp_daily_pro(ivma),                      &
369                 resp_maint_pro(ivma), resp_growth_pro(ivma),                   &
370                 resp_hetero_pro(ivma), co2_fire_pro(ivma),                     &
371                 veget_max, carbon, litter, lignin_struc, bm_to_litter,         &
372                 deepC_a, deepC_s, deepC_p,                                     &
373                 fuel_1hr, fuel_10hr, fuel_100hr, fuel_1000hr,                  &
374                 biomass, co2_to_bm, npp_longterm, ind,                         &
375                 lm_lastyearmax, age, everywhere,                               &
376                 leaf_frac, leaf_age, PFTpresent, senescence,                   &
377                 gpp_daily, npp_daily, resp_maint, resp_growth,                 &
378                 resp_hetero, co2_fire)
379         
380        !ENDIF !IF ( glcc_mtc(ipts,ivma) .GT. min_stomate )
381      ENDDO !(DO ivma=1,nvmap)
382
383    ENDDO !(DO ipts=1,npts)
384
385    !! Update 10 year-turnover pool content following flux emission
386    !!     (linear decay (10%) of the initial carbon input)
387    DO  l = 0, 8
388      m = 10 - l
389      cflux_prod10(:,:) =  cflux_prod10(:,:) + flux10(:,m,:)
390      prod10(:,m,:)     =  prod10(:,m-1,:)   - flux10(:,m-1,:)
391      flux10(:,m,:)     =  flux10(:,m-1,:)
392    ENDDO
393   
394    cflux_prod10(:,:) = cflux_prod10(:,:) + flux10(:,1,:) 
395    flux10(:,1,:)     = 0.1 * prod10(:,0,:)
396    prod10(:,1,:)     = prod10(:,0,:)
397   
398    !! 2.4.3 update 100 year-turnover pool content following flux emission\n
399    DO l = 0, 98
400       m = 100 - l
401       cflux_prod100(:,:)  =  cflux_prod100(:,:) + flux100(:,m,:)
402       prod100(:,m,:)      =  prod100(:,m-1,:)   - flux100(:,m-1,:)
403       flux100(:,m,:)      =  flux100(:,m-1,:)
404    ENDDO
405   
406    cflux_prod100(:,:)  = cflux_prod100(:,:) + flux100(:,1,:) 
407    flux100(:,1,:)      = 0.01 * prod100(:,0,:)
408    prod100(:,1,:)      = prod100(:,0,:)
409    prod10(:,0,:)        = zero
410    prod100(:,0,:)       = zero 
411
412    ! Write out history files
413    CALL histwrite_p (hist_id_stomate, 'glcc_pft', itime, &
414         glcc_pft, npts*nvm, horipft_index)
415    CALL xios_orchidee_send_field ('glcc_pft', glcc_pft)
416
417    DO ivma = 1, nvmap
418      WRITE(part_str,'(I2)') ivma
419      IF (ivma < 10) part_str(1:1) = '0'
420      CALL histwrite_p (hist_id_stomate, 'glcc_pftmtc_'//part_str(1:LEN_TRIM(part_str)), &
421           itime, glcc_pftmtc(:,:,ivma), npts*nvm, horipft_index)
422    ENDDO
423    CALL xios_orchidee_send_field ('glcc_pftmtc', glcc_pftmtc) ! kjpindex,nvm,nvmap
424
425  END SUBROUTINE glcc_bioe1
426
427
428! ================================================================================================================================
429!! SUBROUTINE   : glcc_bioe1_firstday
430!!
431!>\BRIEF        : When necessary, adjust input glcc matrix, and allocate it
432!!                into different contributing age classes and receiving
433!!                youngest age classes.
434!! \n
435!_ ================================================================================================================================
436
437  ! Note: it has this name because this subroutine will also be called
438  ! the first day of each year to precalculate the forest loss for the
439  ! deforestation fire module.
440  SUBROUTINE glcc_bioe1_firstday(npts,veget_max_org,newvegfrac, &
441                          glccSecondShift, glcc_pft, glcc_pftmtc)
442
443    IMPLICIT NONE
444
445    !! 0.1 Input variables
446
447    INTEGER, INTENT(in)                                     :: npts           !! Domain size - number of pixels (unitless)
448    REAL(r_std), DIMENSION(npts,nvm), INTENT(in)            :: veget_max_org  !! "maximal" coverage fraction of a PFT on the ground
449                                                                              !! May sum to
450                                                                              !! less than unity if the pixel has
451                                                                              !! nobio area. (unitless, 0-1)
452    REAL(r_std), DIMENSION(npts,nvmap), INTENT(in)          :: newvegfrac     !! used to guid the allocation of new PFTs.
453                                                                              !!
454    REAL(r_std), DIMENSION (npts,12),INTENT(in)        :: glccSecondShift     !! the land-cover-change (LCC) matrix in case a gross LCC is
455                                                                              !! used.
456
457    !! 0.2 Output variables
458    REAL(r_std), DIMENSION(npts,nvm,nvmap), INTENT(inout)   :: glcc_pftmtc    !! a temporary variable to hold the fractions each PFT is going to lose
459    REAL(r_std), DIMENSION(npts,nvm), INTENT(inout)         :: glcc_pft       !! Loss of fraction in each PFT
460
461    !! 0.3 Modified variables
462   
463    !! 0.4 Local variables
464    REAL(r_std), DIMENSION(npts,nvmap)              :: veget_mtc           !! "maximal" coverage fraction of a PFT on the ground
465    REAL(r_std), DIMENSION(npts,nvmap)              :: veget_mtc_begin     !! "maximal" coverage fraction of a PFT on the ground
466    REAL(r_std), DIMENSION(npts,nagec_tree)         :: vegagec_tree        !! fraction of tree age-class groups, in sequence of old->young
467    REAL(r_std), DIMENSION(npts,nagec_herb)         :: vegagec_grass       !! fraction of grass age-class groups, in sequence of old->young
468    REAL(r_std), DIMENSION(npts,nagec_herb)         :: vegagec_pasture     !! fraction of pasture age-class groups, in sequence of old->young
469    REAL(r_std), DIMENSION(npts,nagec_herb)         :: vegagec_crop        !! fraction of crop age-class groups, in sequence of old->young
470    REAL(r_std), DIMENSION(npts,nagec_bioe1)        :: vegagec_bioe1       !! fraction of crop age-class groups, in sequence of old->young
471
472   
473    REAL(r_std), DIMENSION(npts,4)                  :: veget_4veg      !! "maximal" coverage fraction of a PFT on the ground
474    REAL(r_std), DIMENSION(npts)                    :: veget_tree      !! "maximal" coverage fraction of a PFT on the ground
475    REAL(r_std), DIMENSION(npts)                    :: veget_grass     !! "maximal" coverage fraction of a PFT on the ground
476    REAL(r_std), DIMENSION(npts)                    :: veget_pasture   !! "maximal" coverage fraction of a PFT on the ground
477    REAL(r_std), DIMENSION(npts)                    :: veget_crop      !! "maximal" coverage fraction of a PFT on the ground
478    REAL(r_std), DIMENSION(npts)                    :: veget_bioe1     !! "maximal" coverage fraction of a PFT on the ground
479
480    REAL(r_std), DIMENSION(npts,nvm)         :: veget_max              !! "maximal" coverage fraction of a PFT on the ground
481    REAL(r_std), DIMENSION(npts,nvm)         :: veget_max_old          !! "maximal" coverage fraction of a PFT on the ground
482    REAL(r_std), DIMENSION(npts,nvm)         :: glcc_pft_tmp           !! Loss of fraction in each PFT
483
484    REAL(r_std), DIMENSION(npts,nvm,nvmap)  :: glcc_pftmtc_SecShift   !! a temporary variable to hold the fractions each PFT is going to lose
485    REAL(r_std), DIMENSION(npts,12)      :: glccRealSecShift   !! real matrix applied for secondary shifting cultivation.
486    REAL(r_std), DIMENSION(npts,12)      :: glccDefSecShift    !! deficit for the glccSecondShift
487    REAL(r_std), DIMENSION(npts,12)         :: glccRemain      !!
488    REAL(r_std), DIMENSION(npts,12)         :: glccSecondShift_remain      !!
489
490    REAL(r_std), DIMENSION(npts,nvm)         :: glccReal_tmp         !! A temporary variable
491
492    LOGICAL, SAVE  :: glcc_bioe1_firstday_done = .FALSE.
493
494    ! Different indexes for convenient local uses
495    INTEGER :: f_to_bioe1=1, g_to_bioe1=2, p_to_bioe1=3, c_to_bioe1=4, & 
496               bioe1_to_f=5, bioe1_to_g=6, bioe1_to_p=7, bioe1_to_c=8
497
498    INTEGER :: ivma
499
500    INTEGER :: ipts,IndStart_f,IndEnd_f
501    CHARACTER(LEN=10)      :: part_str                               !! string suffix indicating an index
502
503    !Some more local configurations
504    LOGICAL                                 :: allow_youngest_forest_convert = .TRUE.
505   
506    ! Initialization
507    glcc_pftmtc = zero
508    glcc_pft = zero
509    glcc_pft_tmp = zero
510
511    !!! ** Land cover change processes start here ** !!!
512    ! we make copies of original input veget_max (which is veget_max_org
513    ! in the subroutine parameter list).
514    ! veget_max will be modified through different operations in order to
515    ! check various purposes, e.g., whether input glcc matrix
516    ! is compatible with existing veget_max and how to allocate it etc.
517    ! veget_max_old will not be modified
518    veget_max(:,:) = veget_max_org(:,:)
519    veget_max_old(:,:) = veget_max_org(:,:)
520
521    !! 3. Treat secondary-agriculture shifting cultivation transition matrix
522    !! [The primary-agriculture shifting cultivation will be treated together
523    !!  with the netLCC transitions, with the conversion sequence of oldest->
524    !!  youngest is applied.]
525    ! When we prepare the driving data, secondary-agriculture shifting cultivation
526    ! is intended to include the "constant transitions" over time. Ideally, we
527    ! should start applying this secondary-agriculture shifting cultivation with
528    ! the "secondary forest" in the model. Here we tentatively start with the 3rd
529    ! youngest age class and move to the 2ne youngest age class. But if the prescribed
530    ! transition fraction is not met, we then move further to 4th youngest age class
531    ! and then move to the oldest age class sequentially.
532
533    CALL calc_cover_bioe1(npts,veget_max,veget_mtc_begin,vegagec_tree,vegagec_grass, &
534           vegagec_pasture,vegagec_crop,vegagec_bioe1)
535    veget_mtc = veget_mtc_begin
536 
537    !! 3.1 We start treating secondary-agriculture cultivation from the 3rd youngest
538    !! age class and then move to the younger age class.
539    ! Because it's rather complicated to calculate which transtion fraction between
540    ! which vegetation types should occur in case there is deficit occuring
541    ! for the overall donation vegetation type, we will just start from some
542    ! priority and leave the unrealized parts into the latter section.
543
544    ! For this purpose, we should first make a copy of glccSecondShift into
545    ! glccRemain. glccRemain will tell us the transition fractions that have to
546    ! be treated starting from `IndStart_f+1` oldest age class and moving torward older
547    ! age class.
548    glccRemain(:,:) = glccSecondShift(:,:)
549
550    ! Now we will call type_conversion for each of the 12 transitions, starting
551    ! from `IndStart_f` age class moving to the 2nd youngest age class. We use glccRemain
552    ! to track the transtion fractions we should leave for the second case.
553    ! To make the code more flexible, we will store the start and end indecies
554    ! in variables.
555
556    !*[Note: we do above process only for forest now, as we assume the conversion
557    !  of crop/pasture/grass to other types will start always from the oldest
558    !  age class]
559
560    IndStart_f = nagec_tree-1  ! note the indecies and vegetfrac for tree age class
561                               ! is from old to young
562    IndEnd_f = nagec_tree    ! nagec_tree-2: The 3rd youngest age class
563                               ! nagec_tree-1: The 2nd youngest age class
564                               ! nagec_tree: The youngest age class
565
566    IF (IndStart_f .LE. 0 .OR. IndEnd_f .LE. 0) THEN
567      write(numout,*) 'glcc_bioe1_firstday: Age class index cannot be negative or zero!'
568      STOP
569    ENDIF
570
571    DO ipts=1,npts
572      !f_to_bioe1
573      CALL type_conversion(ipts,f_to_bioe1,glccSecondShift,veget_mtc,newvegfrac,       &
574                       indold_tree,indagec_tree,indagec_bioe1,num_bioe1_mulagec,     &
575                       IndEnd_f,nagec_bioe1,                    &
576                       vegagec_tree,veget_max,glcc_pft,glcc_pftmtc,glcc_pft_tmp, &
577                       glccRemain, &
578                       .TRUE., iagec_start=IndStart_f)
579      !g_to_bioe1
580      CALL type_conversion(ipts,g_to_bioe1,glccSecondShift,veget_mtc,newvegfrac,       &
581                       indold_grass,indagec_grass,indagec_bioe1,num_bioe1_mulagec,     &
582                       nagec_herb,nagec_bioe1,                    &
583                       vegagec_grass,veget_max,glcc_pft,glcc_pftmtc,glcc_pft_tmp,&
584                       glccRemain, &
585                       .TRUE.)
586      !p_to_bioe1
587      CALL type_conversion(ipts,p_to_bioe1,glccSecondShift,veget_mtc,newvegfrac,       &
588                       indold_pasture,indagec_pasture,indagec_bioe1,num_bioe1_mulagec,     &
589                       nagec_herb,nagec_bioe1,                    &
590                       vegagec_pasture,veget_max,glcc_pft,glcc_pftmtc,glcc_pft_tmp,&
591                       glccRemain, &
592                       .TRUE.)
593      !c_to_bioe1
594      CALL type_conversion(ipts,c_to_bioe1,glccSecondShift,veget_mtc,newvegfrac,       &
595                       indold_crop,indagec_crop,indagec_bioe1,num_bioe1_mulagec,     &
596                       nagec_herb,nagec_bioe1,                    &
597                       vegagec_crop,veget_max,glcc_pft,glcc_pftmtc,glcc_pft_tmp,&
598                       glccRemain, &
599                       .TRUE.)
600      !bioe1_to_f
601      CALL type_conversion(ipts,bioe1_to_f,glccSecondShift,veget_mtc,newvegfrac,       &
602                       indold_bioe1,indagec_bioe1,indagec_tree,num_tree_mulagec,     &
603                       nagec_bioe1,nagec_tree,                    &
604                       vegagec_bioe1,veget_max,glcc_pft,glcc_pftmtc,glcc_pft_tmp,&
605                       glccRemain, &
606                       .TRUE.)
607      !bioe1_to_g
608      CALL type_conversion(ipts,bioe1_to_g,glccSecondShift,veget_mtc,newvegfrac,       &
609                       indold_bioe1,indagec_bioe1,indagec_grass,num_grass_mulagec,     &
610                       nagec_bioe1,nagec_herb,                    &
611                       vegagec_bioe1,veget_max,glcc_pft,glcc_pftmtc,glcc_pft_tmp,&
612                       glccRemain, &
613                       .TRUE.)
614      !bioe1_to_p
615      CALL type_conversion(ipts,bioe1_to_p,glccSecondShift,veget_mtc,newvegfrac,       &
616                       indold_bioe1,indagec_bioe1,indagec_pasture,num_pasture_mulagec,     &
617                       nagec_bioe1,nagec_herb,                    &
618                       vegagec_bioe1,veget_max,glcc_pft,glcc_pftmtc,glcc_pft_tmp,&
619                       glccRemain, &
620                       .TRUE.)
621      !bioe1_to_c
622      CALL type_conversion(ipts,bioe1_to_c,glccSecondShift,veget_mtc,newvegfrac,       &
623                       indold_bioe1,indagec_bioe1,indagec_crop,num_crop_mulagec,     &
624                       nagec_bioe1,nagec_herb,                    &
625                       vegagec_bioe1,veget_max,glcc_pft,glcc_pftmtc,glcc_pft_tmp,&
626                       glccRemain, &
627                       .TRUE.)
628    ENDDO
629    glccSecondShift_remain(:,:) = glccRemain(:,:)
630
631    !! 3.2 We treat the remaing unrealized transtions from forest. Now we will
632    !! start with the 3rd oldest age class and then move to the oldest age class.
633
634    CALL calc_cover_bioe1(npts,veget_max,veget_mtc,vegagec_tree,vegagec_grass, &
635           vegagec_pasture,vegagec_crop,vegagec_bioe1)
636    veget_mtc = veget_mtc_begin
637 
638    IndStart_f = nagec_tree  ! note the indecies and vegetfrac for tree age class
639                               ! is from old to young.
640                               ! nagec_tree -3: The 4th youngest age class.
641
642    IndEnd_f = 1               ! oldest-age class forest.
643
644    IF (IndStart_f .LE. 0 .OR. IndEnd_f .LE. 0) THEN
645      write(numout,*) 'glcc_bioe1_firstday: Age class index cannot be negative or zero!'
646      STOP
647    ENDIF
648
649    ! we start with the 3rd youngest age class and move up to the oldest age
650    ! class in the sequence of young->old, as indicated by the .FALSE. parameter
651    ! when calling the subroutine type_conversion.
652    DO ipts=1,npts
653      !f_to_bioe1
654      CALL type_conversion(ipts,f_to_bioe1,glccSecondShift_remain,veget_mtc,newvegfrac,       &
655                       indold_tree,indagec_tree,indagec_bioe1,num_bioe1_mulagec,     &
656                       IndEnd_f,nagec_bioe1,                    &
657                       vegagec_tree,veget_max,glcc_pft,glcc_pftmtc,glcc_pft_tmp, &
658                       glccRemain, &
659                       .FALSE., iagec_start=IndStart_f)
660    ENDDO
661
662    IF (allow_youngest_forest_convert) THEN
663      !!++Temp++!!
664      !! this block of 3.3 could be commented to remove this process as desribed
665      !! below.
666
667      ! [2016-04-20] This is temporarily added: Normally we assume the youngest
668      ! forest age cohort will not be cut because in a shifting cultivation, they
669      ! are grown to let the land recover from agricultural process. (Or at least)
670      ! we can set the threshold of youngest age cohort to be very small. But there
671      ! are two reasons we allow the youngest forest cohort to be cut for shifting
672      ! cultivation purpose: a). Farmers may decide to harvest the wood of a forest
673      ! and then convert to crop. We don't simulate explicitly this process because
674      ! this will depend on input land change matrix and land use data assumptions.
675      ! However,we can implicitly account for this by assuming "farmers plant young
676      ! trees after harvesting the wood, and immediately convert this young trees
677      ! to crops. b). For the sake of conserving the total sum of veget_max before
678      ! and after the transition as one, we need to allow the youngest forest cohort
679      ! eligible for cutting.
680
681      !! 3.3 We treat the remaing unrealized transtions from forest, allowing
682      !! the youngest forest cohort to be cut. For this purpose, we will
683      !! start with the 2nd youngest age class and then move to the youngest one.
684
685      glccSecondShift_remain(:,:) = glccRemain(:,:)
686
687      CALL calc_cover_bioe1(npts,veget_max,veget_mtc,vegagec_tree,vegagec_grass, &
688           vegagec_pasture,vegagec_crop,vegagec_bioe1)
689      veget_mtc = veget_mtc_begin
690 
691      ! Note: the setting of index here must be consistent with those of 3.1 and 3.2
692      IndStart_f = nagec_tree-1  ! note the indecies and vegetfrac for tree age class
693                                 ! is from old to young.
694                                 ! nagec_tree -1: The 2nd youngest age class.
695
696      IndEnd_f = nagec_tree      ! youngest class forest.
697
698      IF (IndStart_f .LE. 0 .OR. IndEnd_f .LE. 0) THEN
699        write(numout,*) 'glcc_bioe1_firstday: Age class index cannot be negative or zero!'
700        STOP
701      ENDIF
702
703      ! we start with the 3rd youngest age class and move up to the oldest age
704      ! class in the sequence of young->old, as indicated by the .FALSE. parameter
705      ! when calling the subroutine type_conversion.
706      DO ipts=1,npts
707        !f_to_bioe1
708        CALL type_conversion(ipts,f_to_bioe1,glccSecondShift_remain,veget_mtc,newvegfrac,       &
709                         indold_tree,indagec_tree,indagec_bioe1,num_bioe1_mulagec,     &
710                         IndEnd_f,nagec_bioe1,                    &
711                         vegagec_tree,veget_max,glcc_pft,glcc_pftmtc,glcc_pft_tmp, &
712                         glccRemain, &
713                         .FALSE., iagec_start=IndStart_f)
714      ENDDO
715      !! End of ++Temp++ Section 3.3
716    ENDIF
717
718    ! Final handling of some output variables.
719    ! we separate the glcc_pftmtc_SecShift
720    glcc_pftmtc_SecShift = glcc_pftmtc
721
722    ! we put the remaining glccRemain into the deficit
723    glccDefSecShift = -1 * glccRemain
724    glccRealSecShift = glccSecondShift - glccRemain
725
726    !*****end block to handle glcc involving bioenergy vegtation type *******
727
728    IF (.NOT. glcc_bioe1_firstday_done) THEN
729
730      glccReal_tmp = zero
731
732      glccReal_tmp(:,1:12) = glccRealSecShift
733      CALL histwrite_p (hist_id_stomate, 'glccRealSecShift', itime, &
734           glccReal_tmp, npts*nvm, horipft_index)
735      CALL xios_orchidee_send_field ('glccRealSecShift', glccReal_tmp) ! kjpindex,nvm
736
737      glccReal_tmp(:,1:12) = glccDefSecShift
738      CALL histwrite_p (hist_id_stomate, 'glccDefSecShift', itime, &
739           glccReal_tmp, npts*nvm, horipft_index)
740      CALL xios_orchidee_send_field ('glccDefSecShift', glccReal_tmp) ! kjpindex,nvm
741
742      DO ivma = 1, nvmap
743        WRITE(part_str,'(I2)') ivma
744        IF (ivma < 10) part_str(1:1) = '0'
745        CALL histwrite_p (hist_id_stomate, 'glcc_pftmtc_SF_'//part_str(1:LEN_TRIM(part_str)), &
746             itime, glcc_pftmtc_SecShift(:,:,ivma), npts*nvm, horipft_index)
747      ENDDO
748      CALL xios_orchidee_send_field ('glcc_pftmtc_SecShift', glcc_pftmtc_SecShift) ! kjpindex,nvm,nvmap
749
750      glcc_bioe1_firstday_done = .TRUE.
751    ENDIF
752
753  END SUBROUTINE glcc_bioe1_firstday
754
755  SUBROUTINE calc_cover_bioe1(npts,veget_max,veget_mtc,vegagec_tree,vegagec_grass, &
756                 vegagec_pasture,vegagec_crop,vegagec_bioe1)
757
758   
759    IMPLICIT NONE
760
761    !! Input variables
762    INTEGER, INTENT(in)                                       :: npts             !! Domain size - number of pixels (unitless)
763    REAL(r_std), DIMENSION(npts,nvm), INTENT(in)         :: veget_max           !! "maximal" coverage fraction of a PFT on the ground
764
765    !! Output variables
766    REAL(r_std), DIMENSION(npts,nvmap), INTENT(inout)         :: veget_mtc        !! "maximal" coverage fraction of a PFT on the ground
767    REAL(r_std), DIMENSION(npts,nagec_tree), INTENT(inout)    :: vegagec_tree     !! fraction of tree age-class groups, in sequence of old->young
768    REAL(r_std), DIMENSION(npts,nagec_herb), INTENT(inout)    :: vegagec_grass    !! fraction of grass age-class groups, in sequence of old->young
769    REAL(r_std), DIMENSION(npts,nagec_herb), INTENT(inout)    :: vegagec_pasture  !! fraction of pasture age-class groups, in sequence of old->young
770    REAL(r_std), DIMENSION(npts,nagec_herb), INTENT(inout)    :: vegagec_crop     !! fraction of crop age-class groups, in sequence of old->young
771    REAL(r_std), DIMENSION(npts,nagec_bioe1), INTENT(inout)    :: vegagec_bioe1     !! fraction of bioenergy tree age-class groups, in sequence of old->young
772
773    !! Local variables
774    INTEGER(i_std)                                          :: ivma,staind,endind,j    !! indices (unitless)
775
776    veget_mtc(:,:) = 0.
777    vegagec_tree(:,:) = 0.
778    vegagec_grass(:,:) = 0.
779    vegagec_pasture(:,:) = 0.
780    vegagec_crop(:,:) = 0.
781    vegagec_bioe1(:,:) = 0.
782
783    ! Calculate veget_max for MTCs
784    DO ivma = 1,nvmap
785      staind = start_index(ivma)
786      IF (nagec_pft(ivma) == 1) THEN
787        veget_mtc(:,ivma) = veget_max(:,staind)
788      ELSE
789        veget_mtc(:,ivma) = \
790          SUM(veget_max(:,staind:staind+nagec_pft(ivma)-1),DIM=2)
791      ENDIF
792    ENDDO
793
794    ! Calculate veget_max for each age class
795    DO ivma = 2,nvmap  !here we start with 2 to exclude baresoil (always PFT1)
796      staind = start_index(ivma)
797      endind = staind+nagec_pft(ivma)-1
798
799      ! Single-age-class MTC goest to oldest age class.
800      IF (nagec_pft(ivma) == 1) THEN
801        WRITE(numout,*) "Error: metaclass has only a single age group: ",ivma
802        STOP
803      ELSE
804        IF (is_tree(staind)) THEN
805          DO j=1,nagec_tree
806            vegagec_tree(:,j) = vegagec_tree(:,j)+veget_max(:,endind-j+1)
807          ENDDO
808        ELSE IF (is_bioe1(staind)) THEN
809          DO j=1,nagec_bioe1
810            vegagec_bioe1(:,j) = vegagec_bioe1(:,j)+veget_max(:,endind-j+1)
811          ENDDO
812        ELSE IF (is_grassland_manag(staind)) THEN
813          DO j=1,nagec_herb
814            vegagec_pasture(:,j) = vegagec_pasture(:,j)+veget_max(:,endind-j+1)
815          ENDDO
816        ELSE IF (natural(staind)) THEN
817          DO j=1,nagec_herb
818            vegagec_grass(:,j) = vegagec_grass(:,j)+veget_max(:,endind-j+1)
819          ENDDO
820        ELSE
821          DO j=1,nagec_herb
822            vegagec_crop(:,j) = vegagec_crop(:,j)+veget_max(:,endind-j+1)
823          ENDDO
824        ENDIF
825      ENDIF
826    ENDDO
827
828  END SUBROUTINE calc_cover_bioe1
829
830END MODULE stomate_glcc_bioe1
Note: See TracBrowser for help on using the repository browser.