Context Navigation

← Previous Revision
Latest Revision
Next Revision →
Normal
Revision Log

stomate_lcchange.f90 @ 7326

Last change on this file since 7326 was 7326, checked in by josefine.ghattas, 3 years ago
Corrected bug on carbon balance closure. See ticket #785 Integration in branch 2_2 done by P. Cadule
Property svn:keywords set to `HeadURL Date Author Revision`
File size: 20.1 KB

Rev	Line
[947]	1	! =================================================================================================================================
	2	! MODULE : stomate_lcchange
[8]	3	!
[4470]	4	! CONTACT : orchidee-help _at_ listes.ipsl.fr
[8]	5	!
[947]	6	! LICENCE : IPSL (2006)
	7	! This software is governed by the CeCILL licence see ORCHIDEE/ORCHIDEE_CeCILL.LIC
[8]	8	!
[947]	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$
	19	!! $Date$
	20	!! $Revision$
	21	!! \n
	22	!_ ================================================================================================================================
	23
	24
[8]	25	MODULE stomate_lcchange
	26
	27	! modules used:
[947]	28
[1078]	29	USE ioipsl_para
[8]	30	USE stomate_data
[511]	31	USE pft_parameters
	32	USE constantes
[947]	33
[8]	34	IMPLICIT NONE
[947]	35
[8]	36	PRIVATE
	37	PUBLIC lcchange_main
[947]	38
[8]	39	CONTAINS
	40
[947]	41
	42	!! ================================================================================================================================
	43	!! SUBROUTINE : lcchange_main
	44	!!
	45	!>\BRIEF Impact of land cover change on carbon stocks
	46	!!
[2718]	47	!! DESCRIPTION : This subroutine is always activate if VEGET_UPDATE>0Y in the configuration file, which means that the
	48	!! vegetation map is updated regulary. lcchange_main is called from stomateLpj the first time step after the vegetation
	49	!! map has been changed.
	50	!! The impact of land cover change on carbon stocks is computed in this subroutine. The land cover change is written
	51	!! by the difference of current and previous "maximal" coverage fraction of a PFT.
[947]	52	!! On the basis of this difference, the amount of 'new establishment'/'biomass export',
	53	!! and increase/decrease of each component, are estimated.\n
	54	!!
	55	!! Main structure of lpj_establish.f90 is:
	56	!! 1. Initialization
	57	!! 2. Calculation of changes in carbon stocks and biomass by land cover change
	58	!! 3. Update 10 year- and 100 year-turnover pool contents
	59	!! 4. History
	60	!!
	61	!! RECENT CHANGE(S) : None
	62	!!
	63	!! MAIN OUTPUT VARIABLE(S) : ::prod10, ::prod100, ::flux10, ::flux100,
	64	!! :: cflux_prod10 and :: cflux_prod100
	65	!!
	66	!! REFERENCES : None
	67	!!
	68	!! FLOWCHART :
	69	!! \latexonly
	70	!! \includegraphics[scale=0.5]{lcchange.png}
	71	!! \endlatexonly
	72	!! \n
	73	!_ ================================================================================================================================
	74
	75
[4647]	76	SUBROUTINE lcchange_main ( npts, dt_days, veget_cov_max_old, veget_cov_max_new, &
[947]	77	biomass, ind, age, PFTpresent, senescence, when_growthinit, everywhere, &
[1091]	78	co2_to_bm, bm_to_litter, turnover_daily, bm_sapl, cn_ind,flux10,flux100, &
[4657]	79	prod10,prod100,convflux,cflux_prod10,cflux_prod100,leaf_frac,&
[4872]	80	npp_longterm, lm_lastyearmax, litter, carbon,&
	81	convfluxpft, fDeforestToProduct, fLulccResidue)
	82
[947]	83
[8]	84	IMPLICIT NONE
[947]	85
	86	!! 0. Variable and parameter declaration
	87
	88	!! 0.1 Input variables
	89
	90	INTEGER, INTENT(in) :: npts !! Domain size - number of pixels (unitless)
	91	REAL(r_std), INTENT(in) :: dt_days !! Time step of vegetation dynamics for stomate
	92	!! (days)
[1170]	93	REAL(r_std), DIMENSION(nvm, nparts,nelements), INTENT(in) :: bm_sapl !! biomass of sapling
[947]	94	!! @tex ($gC individual^{-1}$) @endtex
[4647]	95	REAL(r_std), DIMENSION(npts,nvm), INTENT(in) :: veget_cov_max_old!! Current "maximal" coverage fraction of a PFT (LAI
[3094]	96	!! -> infinity) on ground
[4647]	97	REAL(r_std), DIMENSION(npts,nvm), INTENT(in) :: veget_cov_max_new!! New "maximal" coverage fraction of a PFT (LAI ->
[2718]	98	!! infinity) on ground (unitless)
[1091]	99
[947]	100	!! 0.2 Output variables
[8]	101
[947]	102	REAL(r_std), DIMENSION(npts), INTENT(out) :: convflux !! release during first year following land cover
	103	!! change
	104	REAL(r_std), DIMENSION(npts), INTENT(out) :: cflux_prod10 !! total annual release from the 10 year-turnover
	105	!! pool @tex ($gC m^{-2}$) @endtex
	106	REAL(r_std), DIMENSION(npts), INTENT(out) :: cflux_prod100 !! total annual release from the 100 year-
	107	!! turnover pool @tex ($gC m^{-2}$) @endtex
[1170]	108	REAL(r_std), DIMENSION(npts,nvm,nparts,nelements), INTENT(inout):: turnover_daily !! Turnover rates
[4872]	109	REAL(r_std), DIMENSION(npts,nvm), INTENT(out) :: convfluxpft !! release during first year following land cover
	110	!! change
	111	REAL(r_std), DIMENSION(npts,nvm), INTENT(out) :: fDeforestToProduct !! Deforested biomass into product pool due to anthorpogenic
	112	!! land use change
	113	REAL(r_std), DIMENSION(npts,nvm), INTENT(out) :: fLulccResidue !! carbon mass flux into soil and litter due to anthropogenic land use or land cover change
	114
[947]	115	!! 0.3 Modified variables
	116
[1170]	117	REAL(r_std), DIMENSION(npts,nvm,nparts,nelements), INTENT(inout):: biomass !! biomass @tex ($gC m^{-2}$) @endtex
[947]	118	REAL(r_std), DIMENSION(npts,nvm), INTENT(inout) :: ind !! Number of individuals @tex ($m^{-2}$) @endtex
	119	REAL(r_std), DIMENSION(npts,nvm), INTENT(inout) :: age !! mean age (years)
	120	LOGICAL, DIMENSION(npts,nvm), INTENT(inout) :: senescence !! plant senescent (only for deciduous trees) Set
	121	!! to .FALSE. if PFT is introduced or killed
	122	LOGICAL, DIMENSION(npts,nvm), INTENT(inout) :: PFTpresent !! Is pft there (unitless)
	123	REAL(r_std), DIMENSION(npts,nvm), INTENT(inout) :: everywhere !! is the PFT everywhere in the grid box or very
	124	!! localized (unitless)
	125	REAL(r_std), DIMENSION(npts,nvm), INTENT(inout) :: when_growthinit !! how many days ago was the beginning of the
	126	!! growing season (days)
	127	REAL(r_std), DIMENSION(npts,nvm), INTENT(inout) :: co2_to_bm !! biomass uptaken
	128	!! @tex ($gC m^{-2} day^{-1}$) @endtex
[1170]	129	REAL(r_std), DIMENSION(npts,nvm,nparts,nelements), INTENT(inout) :: bm_to_litter !! conversion of biomass to litter
[947]	130	!! @tex ($gC m^{-2} day^{-1}$) @endtex
	131	REAL(r_std), DIMENSION(npts,nvm), INTENT(inout) :: cn_ind !! crown area of individuals
	132	!! @tex ($m^{2}$) @endtex
	133	REAL(r_std), DIMENSION(npts,0:10), INTENT(inout) :: prod10 !! products remaining in the 10 year-turnover
	134	!! pool after the annual release for each
	135	!! compartment (10 + 1 : input from year of land
	136	!! cover change)
	137	REAL(r_std), DIMENSION(npts,0:100), INTENT(inout) :: prod100 !! products remaining in the 100 year-turnover
	138	!! pool after the annual release for each
	139	!! compartment (100 + 1 : input from year of land
	140	!! cover change)
	141	REAL(r_std), DIMENSION(npts,10), INTENT(inout) :: flux10 !! annual release from the 10/100 year-turnover
	142	!! pool compartments
	143	REAL(r_std), DIMENSION(npts,100), INTENT(inout) :: flux100 !! annual release from the 10/100 year-turnover
	144	!! pool compartments
	145	REAL(r_std), DIMENSION(npts,nvm,nleafages), INTENT(inout) :: leaf_frac !! fraction of leaves in leaf age class
	146	!! (unitless)
	147	REAL(r_std), DIMENSION(npts,nvm), INTENT(inout) :: lm_lastyearmax !! last year's maximum leaf mass for each PFT
	148	!! @tex ($gC m^{-2}$) @endtex
	149	REAL(r_std), DIMENSION(npts,nvm), INTENT(inout) :: npp_longterm !! "long term" net primary productivity
	150	!! @tex ($gC m^{-2} year^{-1}$) @endtex
[1170]	151	REAL(r_std),DIMENSION(npts,nlitt,nvm,nlevs,nelements), INTENT(inout):: litter !! metabolic and structural litter, above and
[947]	152	!! below ground @tex ($gC m^{-2}$) @endtex
	153	REAL(r_std),DIMENSION(npts,ncarb,nvm), INTENT(inout) :: carbon !! carbon pool: active, slow, or passive
[8]	154
[4657]	155	!! @tex ($gC m^{-2}$) @endtex
	156
[947]	157	!! 0.4 Local variables
[8]	158
[947]	159	INTEGER(i_std) :: i, j, k, l, m !! indices (unitless)
[1170]	160	REAL(r_std),DIMENSION(npts,nelements) :: bm_new !! biomass increase @tex ($gC m^{-2}$) @endtex
	161	REAL(r_std),DIMENSION(npts,nparts,nelements) :: biomass_loss !! biomass loss @tex ($gC m^{-2}$) @endtex
[947]	162	REAL(r_std) :: above !! aboveground biomass @tex ($gC m^{-2}$) @endtex
[1170]	163	REAL(r_std),DIMENSION(npts,nlitt,nlevs,nelements) :: dilu_lit !! Litter dilution @tex ($gC m^{-2}$) @endtex
[947]	164	REAL(r_std),DIMENSION(npts,ncarb) :: dilu_soil_carbon !! Soil Carbondilution @tex ($gC m^{-2}$) @endtex
	165	REAL(r_std),DIMENSION(nvm) :: delta_veg !! changes in "maximal" coverage fraction of PFT
	166	REAL(r_std) :: delta_veg_sum !! sum of delta_veg
[4872]	167	REAL(r_std),DIMENSION(npts,nvm) :: delta_ind !! change in number of individuals
	168
[947]	169	!_ ================================================================================================================================
[8]	170
[2348]	171	IF (printlev>=3) WRITE(numout,*) 'Entering lcchange_main'
[947]	172
	173	!! 1. initialization
	174
	175	prod10(:,0) = zero
	176	prod100(:,0) = zero
	177	above = zero
	178	convflux(:) = zero
[4872]	179	convfluxpft(:,:) = zero
[947]	180	cflux_prod10(:) = zero
	181	cflux_prod100(:) = zero
	182	delta_ind(:,:) = zero
	183	delta_veg(:) = zero
[4872]	184	fDeforestToProduct(:,:) = zero
	185	fLulccResidue(:,:) = zero
[4657]	186
	187
[947]	188
[4657]	189	!! 3. calculation of changes in carbon stocks and biomass by land cover change\n
	190
[947]	191	DO i = 1, npts ! Loop over # pixels - domain size
	192
[4657]	193	!! 3.1 initialization of carbon stocks\n
[4647]	194	delta_veg(:) = veget_cov_max_new(i,:)-veget_cov_max_old(i,:)
[8]	195	delta_veg_sum = SUM(delta_veg,MASK=delta_veg.LT.0.)
[947]	196
[1170]	197	dilu_lit(i,:,:,:) = zero
[8]	198	dilu_soil_carbon(i,:) = zero
[1170]	199	biomass_loss(i,:,:) = zero
[947]	200
[4657]	201	!! 3.2 if vegetation coverage decreases, compute dilution of litter, soil carbon, and biomass.\n
[7326]	202	DO j=1, nvm
[8]	203	IF ( delta_veg(j) < -min_stomate ) THEN
[1170]	204	dilu_lit(i,:,:,:) = dilu_lit(i,:,:,:) + delta_veg(j)*litter(i,:,j,:,:) / delta_veg_sum
[947]	205	dilu_soil_carbon(i,:) = dilu_soil_carbon(i,:) + delta_veg(j) * carbon(i,:,j) / delta_veg_sum
[1170]	206	biomass_loss(i,:,:) = biomass_loss(i,:,:) + biomass(i,j,:,:)*delta_veg(j) / delta_veg_sum
[8]	207	ENDIF
	208	ENDDO
[947]	209
[4657]	210	!! 3.3
[7326]	211	DO j=1, nvm ! Loop over # PFTs
[8]	212
[4657]	213	!! 3.3.1 The case that vegetation coverage of PFTj increases
[8]	214	IF ( delta_veg(j) > min_stomate) THEN
[947]	215
[4657]	216	!! 3.3.1.1 Initial setting of new establishment
[4647]	217	IF (veget_cov_max_old(i,j) .LT. min_stomate) THEN
[1091]	218	IF (is_tree(j)) THEN
[947]	219
	220	! cn_sapl(j)=0.5; stomate_data.f90
	221	cn_ind(i,j) = cn_sapl(j)
[8]	222	ELSE
[511]	223	cn_ind(i,j) = un
[8]	224	ENDIF
	225	ind(i,j)= delta_veg(j) / cn_ind(i,j)
	226	PFTpresent(i,j) = .TRUE.
[947]	227	everywhere(i,j) = 1.
[8]	228	senescence(i,j) = .FALSE.
[720]	229	age(i,j) = zero
[8]	230
[947]	231	! large_value = 1.E33_r_std
	232	when_growthinit(i,j) = large_value
	233	leaf_frac(i,j,1) = 1.0
[720]	234	npp_longterm(i,j) = npp_longterm_init
[1170]	235	lm_lastyearmax(i,j) = bm_sapl(j,ileaf,icarbon) * ind(i,j)
[8]	236	ENDIF
	237	IF ( cn_ind(i,j) > min_stomate ) THEN
	238	delta_ind(i,j) = delta_veg(j) / cn_ind(i,j)
	239	ENDIF
[947]	240
[4657]	241	!! 3.3.1.2 Update of biomass in each each carbon stock component
[947]	242	!! Update of biomass in each each carbon stock component (leaf, sapabove, sapbelow,
	243	!> heartabove, heartbelow, root, fruit, and carbres)\n
	244	DO k = 1, nparts ! loop over # carbon stock components, nparts = 8; stomate_constant.f90
[1170]	245	DO l = 1,nelements ! loop over # elements
[947]	246
[1170]	247	bm_new(i,l) = delta_ind(i,j) * bm_sapl(j,k,l)
[4647]	248	IF (veget_cov_max_old(i,j) .GT. min_stomate) THEN
[1170]	249
	250	! in the case that bm_new is overestimated compared with biomass?
	251	IF ((bm_new(i,l)/delta_veg(j)) > biomass(i,j,k,l)) THEN
	252	bm_new(i,l) = biomass(i,j,k,l)*delta_veg(j)
	253	ENDIF
[8]	254	ENDIF
[4647]	255	biomass(i,j,k,l) = ( biomass(i,j,k,l) * veget_cov_max_old(i,j) + bm_new(i,l) ) / veget_cov_max_new(i,j)
	256	co2_to_bm(i,j) = co2_to_bm(i,j) + (bm_new(i,icarbon)* dt_days) / (one_year * veget_cov_max_new(i,j))
[1170]	257	END DO ! loop over # elements
[947]	258	ENDDO ! loop over # carbon stock components
[1170]	259
[4657]	260	!! 3.3.1.3 Calculation of dilution in litter, soil carbon, and input of litter
[947]	261	!! In this 'IF statement', dilu_* is zero. Formulas for litter and soil carbon
	262	!! could be shortend?? Are the following formulas correct?
[8]	263
	264	! Litter
[4647]	265	litter(i,:,j,:,:)=(litter(i,:,j,:,:) * veget_cov_max_old(i,j) + &
	266	dilu_lit(i,:,:,:) * delta_veg(j)) / veget_cov_max_new(i,j)
[947]	267
[8]	268	! Soil carbon
[4647]	269	carbon(i,:,j)=(carbon(i,:,j) * veget_cov_max_old(i,j) + dilu_soil_carbon(i,:) * delta_veg(j)) / veget_cov_max_new(i,j)
[8]	270
[1170]	271	DO l = 1,nelements
	272
	273	! Litter input
[5536]	274	bm_to_litter(i,j,isapbelow,l) = (bm_to_litter(i,j,isapbelow,l) * veget_cov_max_old(i,j) + &
	275	biomass_loss(i,isapbelow,l)*delta_veg(j)) / veget_cov_max_new(i,j)
	276	bm_to_litter(i,j,iheartbelow,l) = (bm_to_litter(i,j,iheartbelow,l) * veget_cov_max_old(i,j) + &
	277	biomass_loss(i,iheartbelow,l) *delta_veg(j)) / veget_cov_max_new(i,j)
	278	bm_to_litter(i,j,iroot,l) = (bm_to_litter(i,j,iroot,l) * veget_cov_max_old(i,j) + &
	279	biomass_loss(i,iroot,l)*delta_veg(j)) / veget_cov_max_new(i,j)
	280	bm_to_litter(i,j,ifruit,l) = (bm_to_litter(i,j,ifruit,l) * veget_cov_max_old(i,j) + &
	281	biomass_loss(i,ifruit,l)*delta_veg(j)) / veget_cov_max_new(i,j)
	282	bm_to_litter(i,j,icarbres,l) = (bm_to_litter(i,j,icarbres,l) * veget_cov_max_old(i,j) + &
	283	biomass_loss(i,icarbres,l) *delta_veg(j)) / veget_cov_max_new(i,j)
	284	bm_to_litter(i,j,ileaf,l) = (bm_to_litter(i,j,ileaf,l) * veget_cov_max_old(i,j) + &
	285	biomass_loss(i,ileaf,l)*delta_veg(j)) / veget_cov_max_new(i,j)
[1170]	286
	287	END DO
	288
[4647]	289	age(i,j)=age(i,j)*veget_cov_max_old(i,j)/veget_cov_max_new(i,j)
[947]	290
[4657]	291	!! 3.3.2 The case that vegetation coverage of PFTj is no change or decreases
[947]	292	ELSE
	293
[4657]	294	!! 3.3.2.1 Biomass export
[947]	295	! coeff_lcchange_*: Coeff of biomass export for the year, decade, and century
[1170]	296	above = biomass(i,j,isapabove,icarbon) + biomass(i,j,iheartabove,icarbon)
[947]	297	convflux(i) = convflux(i) - ( coeff_lcchange_1(j) * above * delta_veg(j) )
[4872]	298	convfluxpft(i,j)= convfluxpft(i,j) - (coeff_lcchange_1(j) * above * delta_veg(j) )
[8]	299	prod10(i,0) = prod10(i,0) - ( coeff_lcchange_10(j) * above * delta_veg(j) )
	300	prod100(i,0) = prod100(i,0) - ( coeff_lcchange_100(j) * above * delta_veg(j) )
[2718]	301
[4872]	302	fDeforestToProduct(i,j) = - above * delta_veg(j)
	303	fLulccResidue(i,j) = - ( biomass(i,j,isapbelow,icarbon) &
	304	+ biomass(i,j,iheartbelow,icarbon) &
	305	+ biomass(i,j,iroot,icarbon) &
	306	+ biomass(i,j,ifruit,icarbon) &
	307	+ biomass(i,j,icarbres,icarbon) &
	308	+ biomass(i,j,ileaf,icarbon) ) * delta_veg(j)
[4657]	309	!! 3.3.2.2 Total reduction
[3094]	310	!! If the vegetation is to small, it has been set to 0.
[4647]	311	IF ( veget_cov_max_new(i,j) .LT. min_stomate ) THEN
[3094]	312
	313	ind(i,j) = zero
	314	biomass(i,j,:,:) = zero
	315	PFTpresent(i,j) = .FALSE.
	316	senescence(i,j) = .FALSE.
	317	age(i,j) = zero
	318	when_growthinit(i,j) = undef
	319	everywhere(i,j) = zero
	320	carbon(i,:,j) = zero
	321	litter(i,:,j,:,:) = zero
	322	bm_to_litter(i,j,:,:) = zero
	323	turnover_daily(i,j,:,:) = zero
	324
	325	ENDIF
	326
[947]	327	ENDIF ! End if PFT's coverage reduction
	328
	329	ENDDO ! Loop over # PFTs
	330
[4657]	331	!! 3.4 update 10 year-turnover pool content following flux emission
[947]	332	!! (linear decay (10%) of the initial carbon input)
	333	DO l = 0, 8
[8]	334	m = 10 - l
	335	cflux_prod10(i) = cflux_prod10(i) + flux10(i,m)
	336	prod10(i,m) = prod10(i,m-1) - flux10(i,m-1)
	337	flux10(i,m) = flux10(i,m-1)
[947]	338
[511]	339	IF (prod10(i,m) .LT. 1.0) prod10(i,m) = zero
[8]	340	ENDDO
[947]	341
[8]	342	cflux_prod10(i) = cflux_prod10(i) + flux10(i,1)
	343	flux10(i,1) = 0.1 * prod10(i,0)
	344	prod10(i,1) = prod10(i,0)
[947]	345
[4657]	346	!! 3.5 update 100 year-turnover pool content following flux emission\n
[8]	347	DO l = 0, 98
	348	m = 100 - l
	349	cflux_prod100(i) = cflux_prod100(i) + flux100(i,m)
	350	prod100(i,m) = prod100(i,m-1) - flux100(i,m-1)
	351	flux100(i,m) = flux100(i,m-1)
[947]	352
[511]	353	IF (prod100(i,m).LT.1.0) prod100(i,m) = zero
[8]	354	ENDDO
[947]	355
[8]	356	cflux_prod100(i) = cflux_prod100(i) + flux100(i,1)
	357	flux100(i,1) = 0.01 * prod100(i,0)
	358	prod100(i,1) = prod100(i,0)
[42]	359	prod10(i,0) = zero
[947]	360	prod100(i,0) = zero
[4657]	361
[947]	362
	363	ENDDO ! Loop over # pixels - domain size
	364
[4657]	365	!! 4. history
[8]	366	convflux = convflux/one_year*dt_days
[4872]	367	convfluxpft = convfluxpft/one_year*dt_days
	368	fDeforestToProduct= fDeforestToProduct/one_year*dt_days
	369	fLulccResidue = fLulccResidue/one_year*dt_days
[8]	370	cflux_prod10 = cflux_prod10/one_year*dt_days
	371	cflux_prod100 = cflux_prod100/one_year*dt_days
[4657]	372
	373
[2348]	374	IF (printlev>=4) WRITE(numout,*) 'Leaving lcchange_main'
[947]	375
[8]	376	END SUBROUTINE lcchange_main
[947]	377
[8]	378	END MODULE stomate_lcchange

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

Context Navigation

source: branches/ORCHIDEE_2_2/ORCHIDEE/src_stomate/stomate_lcchange.f90 @ 7326

Download in other formats: