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

Last change on this file since 6393 was 4872, checked in by nicolas.vuichard, 6 years ago

new output variables for LUMIP

File size: 9.9 KB
Line 
1! =================================================================================================================================
2! MODULE       : stomate_woodharvest
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       Impact of land cover change on carbon stocks
10!!
11!!\n DESCRIPTION: None
12!!
13!! RECENT CHANGE(S): None
14!!
15!! REFERENCE(S) : None
16!!
17!! SVN          :
18!! $HeadURL: svn://forge.ipsl.jussieu.fr/orchidee/trunk/ORCHIDEE/src_stomate/stomate_woodharvest.f90 $
19!! $Date: 2016-01-06 13:15:53 +0100 (mer. 06 janv. 2016) $
20!! $Revision: 3094 $
21!! \n
22!_ ================================================================================================================================
23
24
25MODULE stomate_woodharvest
26
27  USE ioipsl_para
28  USE stomate_data
29  USE pft_parameters
30  USE constantes
31 
32  IMPLICIT NONE
33 
34  PRIVATE
35  PUBLIC woodharvest_main
36 
37CONTAINS
38
39
40!! ================================================================================================================================
41!! SUBROUTINE   : woodharvest_main
42!!
43!>\BRIEF        Impact of wood harvest on carbon stocks
44!!
45!! DESCRIPTION  : This subroutine is always activate if DO_WOOD_HARVEST=y in the configuration file.
46!! woodharvest_main is called from stomateLpj every year (at the beginning)
47!! The impact of wood harvest on carbon stocks is computed in this subroutine. The land cover change is written
48!! by the difference of current and previous "maximal" coverage fraction of a PFT.
49!! On the basis of this difference, the amount of 'new establishment'/'biomass export',
50!! and increase/decrease of each component, are estimated.\n
51!!
52!! Main structure of woodharvest is:
53!! 1. Initialization
54!! 2. Calculation of changes in carbon stocks and biomass by wood harvest
55!! 3. Update 10 year- and 100 year-turnover pool contents
56!! 4. History
57!!
58!! RECENT CHANGE(S) : None
59!!
60!! MAIN OUTPUT VARIABLE(S) : ::prod10, ::prod100, ::flux10, ::flux100,
61!!   :: cflux_prod10 and :: cflux_prod100
62!!
63!! REFERENCES   : None
64!!
65!! FLOWCHART    :
66!! \n
67!_ ================================================================================================================================
68
69 
70  SUBROUTINE woodharvest_main ( npts, dt_days, veget_max_old, &
71       biomass, &       
72       flux10_harvest,flux100_harvest, prod10_harvest,prod100_harvest,&
73       convflux_harvest,cflux_prod10_harvest,cflux_prod100_harvest,&
74       harvestwood,harvestpft,fHarvestToProduct)
75   
76    IMPLICIT NONE
77   
78  !! 0. Variable and parameter declaration
79   
80    !! 0.1 Input variables
81   
82    INTEGER, INTENT(in)                            :: npts             !! Domain size - number of pixels (unitless)
83    REAL(r_std), INTENT(in)                        :: dt_days          !! Time step of vegetation dynamics for stomate (days)
84    REAL(r_std), DIMENSION(npts,nvm), INTENT(in)   :: veget_max_old    !! Current "maximal" coverage fraction of a PFT (LAI
85                                                                       !! -> infinity) on ground
86    REAL(r_std), DIMENSION(npts), INTENT(in)       :: harvestwood      !! Harvested wood biomass (gC m-2 yr-1)
87 
88    !! 0.2 Output variables
89
90    REAL(r_std), DIMENSION(npts), INTENT(out)      :: convflux_harvest        !! flux release during first year following wood harvest
91    REAL(r_std), DIMENSION(npts), INTENT(out)      :: cflux_prod10_harvest    !! total annual release from the 10 year-turnover
92                                                                              !! pool @tex ($gC m^{-2}$) @endtex
93    REAL(r_std), DIMENSION(npts), INTENT(out)      :: cflux_prod100_harvest   !! total annual release from the 100 year-
94                                                                              !! turnover pool @tex ($gC m^{-2}$) @endtex
95    REAL(r_std), DIMENSION(npts,nvm), INTENT(out)  :: harvestpft              !! wood harvested for each PFT
96                                                                              !! @tex ($gC m^{-2} day^{-1}$) @endtex
97    REAL(r_std), DIMENSION(npts,nvm), INTENT(out)  :: fHarvestToProduct       !!  Harvested biomass into product pool due to anthorpogenic
98                                                                              !! land use (gC m-2 )
99!! 0.3 Modified variables   
100   
101    REAL(r_std), DIMENSION(npts,nvm,nparts,nelements), INTENT(inout) :: biomass    !! biomass @tex ($gC m^{-2}$) @endtex
102    REAL(r_std), DIMENSION(npts,0:10), INTENT(inout)  :: prod10_harvest  !! products remaining in the 10/100 year-turnover pool
103                                                                         !! after the annual release for each compartment
104                                                                         !! (10 or 100 + 1 : input from year of land cover change)
105    REAL(r_std), DIMENSION(npts,0:100), INTENT(inout) :: prod100_harvest !! see prod10_havest
106    REAL(r_std), DIMENSION(npts,10), INTENT(inout)    :: flux10_harvest  !! annual release from the 10 year-turnover pool compartments
107
108    REAL(r_std), DIMENSION(npts,100), INTENT(inout)   :: flux100_harvest !! annual release from the 100 year-turnover pool compartments
109                                                   
110    !! 0.4 Local variables
111
112    INTEGER(i_std)            :: i, j, l, m    !! indices (unitless)
113    REAL(r_std)               :: harvest       !! wood harvest renomalized by forest fraction
114    REAL(r_std)               :: relharvest    !! relative harvest fraction of each pft
115    REAL(r_std)               :: relbiomass    !! relative biomass of each PFT
116    REAL(r_std)               :: vegettree     !! vegfrac of trees
117    REAL(r_std)               :: abovemean     !! mean aboveground biomass @tex ($gC m^{-2}$) @endtex
118!_ ================================================================================================================================
119
120    IF (printlev>=3) WRITE(numout,*) 'Entering woodharvest_main'
121   
122  !! 1. initialization
123   
124    !! 2. calculation of changes in carbon stocks and biomass by wood harvest \n
125    !! 2.1 initialization of carbon stocks\n
126    prod10_harvest(:,0)           = zero
127    prod100_harvest(:,0)          = zero   
128    convflux_harvest(:)           = zero
129    cflux_prod10_harvest(:)       = zero
130    cflux_prod100_harvest(:)      = zero
131    harvestpft(:,:)               = zero
132
133    fHarvestToProduct(:,:)        = zero
134
135    DO i = 1, npts 
136       abovemean=0         
137       vegettree=0
138       !! 2.3 calculation of total forest fraction and mean biomass
139       DO j=1, nvm
140          if (is_tree(j)) THEN
141             abovemean=abovemean+((biomass(i,j,iheartabove,icarbon)+biomass(i,j,isapabove,icarbon)))*veget_max_old(i,j)
142             vegettree=vegettree+veget_max_old(i,j)
143          ENDIF
144       ENDDO
145       IF (vegettree > min_stomate) THEN
146         
147          harvest=harvestwood(i) 
148          abovemean=abovemean/vegettree
149
150          !! 2.4 wood harvest for each PFT is estimated proportionally to the relative biomass.
151          DO j=1,nvm
152             relharvest=0
153             IF ((is_tree(j)) .AND. (abovemean > min_stomate)  .AND. &
154                  (biomass(i,j,iheartabove,icarbon)+biomass(i,j,isapabove,icarbon)>min_stomate))  THEN
155                relbiomass=((biomass(i,j,iheartabove,icarbon)+biomass(i,j,isapabove,icarbon))/abovemean) 
156                relharvest=harvest*relbiomass
157
158                !! 2.5 the harvest could not be higher than 20% of the biomasse to avoid rapid killing of vegetation
159                IF ((biomass(i,j,iheartabove,icarbon)+biomass(i,j,isapabove,icarbon))*0.2<relharvest) &
160                     relharvest=(biomass(i,j,iheartabove,icarbon)+biomass(i,j,isapabove,icarbon))*0.2; 
161                biomass(i,j,iheartabove,:)= biomass(i,j,iheartabove,:) - &
162                     relharvest*(biomass(i,j,iheartabove,:)/(biomass(i,j,iheartabove,:)+biomass(i,j,isapabove,:)))
163                biomass(i,j,isapabove,:)= biomass(i,j,isapabove,:) - &
164                     relharvest*(biomass(i,j,isapabove,:)/(biomass(i,j,iheartabove,:)+biomass(i,j,isapabove,:)))
165
166                fHarvestToProduct(i,j)=relharvest*veget_max_old(i,j)
167                ! 2.6 calculation of different flux like for LU
168                convflux_harvest(i)  = convflux_harvest(i) +coeff_lcchange_1(j)*relharvest*veget_max_old(i,j)
169                harvestpft(i,j)=relharvest*dt_days/one_year
170                prod10_harvest(i,0)  = prod10_harvest(i,0)  +coeff_lcchange_10(j)*relharvest*veget_max_old(i,j)
171                prod100_harvest(i,0) = prod100_harvest(i,0) +coeff_lcchange_100(j)*relharvest*veget_max_old(i,j)
172             ENDIF
173          ENDDO
174       ENDIF
175
176       cflux_prod10_harvest(i) = cflux_prod10_harvest(i) + flux10_harvest(i,1) 
177       flux10_harvest(i,1)     = 0.1 * prod10_harvest(i,0)
178       prod10_harvest(i,1)     = prod10_harvest(i,0)
179       
180       !! 3.5 update 100 year-turnover pool content following flux emission\n
181       DO   l = 0, 98
182          m = 100 - l
183          cflux_prod100_harvest(i)  =  cflux_prod100_harvest(i) + flux100_harvest(i,m)
184          prod100_harvest(i,m)      =  prod100_harvest(i,m-1)   - flux100_harvest(i,m-1)
185          flux100_harvest(i,m)      =  flux100_harvest(i,m-1)
186          IF (prod100_harvest(i,m).LT.1.0) prod100_harvest(i,m) = zero
187       ENDDO
188       
189       cflux_prod100_harvest(i)  = cflux_prod100_harvest(i) + flux100_harvest(i,1) 
190       flux100_harvest(i,1)      = 0.01 * prod100_harvest(i,0)
191       prod100_harvest(i,1)      = prod100_harvest(i,0)
192       prod10_harvest(i,0)        = 0.0
193       prod100_harvest(i,0)       = 0.0 
194       
195    ENDDO ! Loop over # pixels - domain size
196   
197  !! 4. history
198    convflux_harvest        = convflux_harvest/one_year*dt_days
199    cflux_prod10_harvest    = cflux_prod10_harvest/one_year*dt_days
200    cflux_prod100_harvest   = cflux_prod100_harvest/one_year*dt_days
201    fHarvestToProduct       = fHarvestToProduct/one_year*dt_days   
202    IF (printlev>=4) WRITE(numout,*) 'Leaving woodharvest_main'
203   
204  END SUBROUTINE woodharvest_main
205 
206END MODULE stomate_woodharvest
Note: See TracBrowser for help on using the repository browser.