Context Navigation

stomate_prescribe.f90 @ 7108

Last change on this file since 7108 was 6849, checked in by yidi.xu, 4 years ago
ORCHIDEE-MICT-OP for oil palm growth modelling
Property svn:keywords set to `HeadURL Date Author Revision`
File size: 14.0 KB

Line
1	! =================================================================================================================================
2	! MODULE : stomate_prescribe
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 Initialize and update density, crown area.
10	!!
11	!!\n DESCRIPTION: None
12	!!
13	!! RECENT CHANGE(S): None
14	!!
15	!! REFERENCE(S) :
16	!!
17	!! SVN :
18	!! $HeadURL$
19	!! $Date$
20	!! $Revision$
21	!! \n
22	!_ ================================================================================================================================
23
24	MODULE stomate_prescribe
25
26	! modules used:
27
28	USE ioipsl_para
29	USE stomate_data
30	USE pft_parameters
31	USE constantes
32
33	IMPLICIT NONE
34
35	! private & public routines
36
37	PRIVATE
38	PUBLIC prescribe,prescribe_clear
39
40	! first call
41	LOGICAL, SAVE :: firstcall_prescribe = .TRUE.
42	!$OMP THREADPRIVATE(firstcall_prescribe)
43
44	CONTAINS
45
46	! =================================================================================================================================
47	!! SUBROUTINE : prescribe_clear
48	!!
49	!>\BRIEF : Set the firstcall_prescribe flag back to .TRUE. to prepare for the next simulation.
50	!_=================================================================================================================================
51
52	SUBROUTINE prescribe_clear
53	firstcall_prescribe=.TRUE.
54	END SUBROUTINE prescribe_clear
55
56	!! ================================================================================================================================
57	!! SUBROUTINE : prescribe
58	!!
59	!>\BRIEF Works only with static vegetation and agricultural PFT. Initialize biomass,
60	!! density, crown area in the first call and update them in the following.
61	!!
62	!! DESCRIPTION (functional, design, flags): \n
63	!! This module works only with static vegetation and agricultural PFT.
64	!! In the first call, initialize density of individuals, biomass, crown area,
65	!! and leaf age distribution to some reasonable value. In the following calls,
66	!! these variables are updated.
67	!!
68	!! To fulfill these purposes, pipe model are used:
69	!! \latexonly
70	!! \input{prescribe1.tex}
71	!! \input{prescribe2.tex}
72	!! \endlatexonly
73	!!
74	!! RECENT CHANGE(S): None
75	!!
76	!! MAIN OUTPUT VARIABLES(S): ::ind, ::cn_ind, ::leaf_frac
77	!!
78	!! REFERENCES :
79	!! - Krinner, G., N. Viovy, et al. (2005). "A dynamic global vegetation model
80	!! for studies of the coupled atmosphere-biosphere system." Global
81	!! Biogeochemical Cycles 19: GB1015, doi:1010.1029/2003GB002199.
82	!! - Sitch, S., B. Smith, et al. (2003), Evaluation of ecosystem dynamics,
83	!! plant geography and terrestrial carbon cycling in the LPJ dynamic
84	!! global vegetation model, Global Change Biology, 9, 161-185.
85	!!
86	!! FLOWCHART : None
87	!! \n
88	!_ ================================================================================================================================
89
90	SUBROUTINE prescribe (npts, &
91	veget_cov_max, dt, PFTpresent, everywhere, when_growthinit, &
92	biomass, leaf_frac, ind, cn_ind, co2_to_bm)
93
94	!! 0. Parameters and variables declaration
95
96	!! 0.1 Input variables
97
98	INTEGER(i_std), INTENT(in) :: npts !! Domain size (unitless)
99	REAL(r_std), DIMENSION(npts,nvm), INTENT(in) :: veget_cov_max !! "maximal" coverage fraction of a PFT (LAI -> infinity) on ground (unitless;0-1)
100	REAL(r_std), INTENT(in) :: dt !! time step (dt_days)
101	!! 0.2 Output variables
102
103	!! 0.3 Modified variables
104
105	LOGICAL, DIMENSION(npts,nvm), INTENT(inout) :: PFTpresent !! PFT present (0 or 1)
106	REAL(r_std), DIMENSION(npts,nvm), INTENT(inout) :: everywhere !! is the PFT everywhere in the grid box or very localized (after its introduction) (?)
107	REAL(r_std), DIMENSION(npts,nvm), INTENT(inout) :: when_growthinit !! how many days ago was the beginning of the growing season (days)
108	REAL(r_std), DIMENSION(npts,nvm,nparts,nelements), INTENT(inout) :: biomass !! biomass (gC/(m^2 of ground))
109	REAL(r_std), DIMENSION(npts,nvm,nleafages), INTENT(inout) :: leaf_frac !! fraction of leaves in leaf age class (unitless;0-1)
110	REAL(r_std), DIMENSION(npts,nvm), INTENT(inout) :: ind !! density of individuals (1/(m^2 of ground))
111	REAL(r_std), DIMENSION(npts,nvm), INTENT(inout) :: cn_ind !! crown area per individual (m^2)
112	REAL(r_std), DIMENSION(npts,nvm), INTENT(inout) :: co2_to_bm !! co2 taken up by carbohydrate
113	!! reserve at the beginning of the
114	!! growing season @tex ($gC m^{-2}
115	!! of total ground/day$) @endtex
116
117	!! 0.4 Local variables
118
119	REAL(r_std), DIMENSION(npts) :: dia !! stem diameter (m)
120	REAL(r_std), DIMENSION(npts) :: woodmass !! woodmass (gC/(m^2 of ground))
121	REAL(r_std), DIMENSION(npts) :: woodmass_ind !! woodmass of an individual (gC)
122	INTEGER(i_std) :: i,j !! index (unitless)
123
124	!_ ================================================================================================================================
125
126	DO j = 2,nvm
127
128	! only when the DGVM is not activated or agricultural PFT.
129
130	IF ( ( .NOT. ok_dgvm .AND. lpj_gap_const_mort ) .OR. ( .NOT. natural(j) .OR. pasture(j)) ) THEN
131
132	!
133	!! 1.Update crown area
134	!
135
136	cn_ind(:,j) = zero
137
138	IF ( is_tree(j) ) THEN
139
140	!
141	!! 1.1 treat for trees
142	!
143
144	dia(:) = zero
145
146	DO i = 1, npts ! loop over grid points
147
148	IF ( veget_cov_max(i,j) .GT. zero ) THEN
149
150	!! 1.1.1 calculate wood mass on an area basis, which include sapwood and heartwood aboveground and belowground.
151
152	woodmass(i) = (biomass(i,j,isapabove,icarbon) + biomass(i,j,isapbelow,icarbon) + &
153	biomass(i,j,iheartabove,icarbon) + biomass(i,j,iheartbelow,icarbon)) * veget_cov_max(i,j)
154
155	IF ( woodmass(i) .GT. min_stomate ) THEN
156
157	!! 1.1.2 calculate critical individual density
158	!?? the logic for 1.1.3 and 1.1.2 is strange, it should be the case that first to calculate critical woodmass per individual,
159	!?? then calculate critical density.
160
161
162	! how to derive the following equation:
163	! first, TreeHeight=pipe_tune2 * Diameter^{pipe_tune3}
164	! we assume the tree is an ideal cylinder, so it volume is: Volume = pi(Dia/2)^2Height = pi/4 * Dia * pipe_tune2*Dia^{pipe_tune3}
165	! = pi/4 * pipe_tune2 * Dia^{2+pipe_tune3}
166	! last, the woodmass_per_individual = pipe_density * Volume = pipe_densitypi/4.pipe_tune2 * Dia^{2+pipe_tune3}
167	ind(i,j) = woodmass(i) / &
168	( pipe_densitypi/4.pipe_tune2 * maxdia(j)**(2.+pipe_tune3) )
169
170	!! 1.1.3 individual biomass corresponding to this critical density of individuals
171
172	woodmass_ind(i) = woodmass(i) / ind(i,j)
173
174	!! 1.1.4 calculate stem diameter per individual tree
175
176	dia(i) = ( woodmass_ind(i) / ( pipe_density * pi/4. * pipe_tune2 ) ) ** &
177	( un / ( 2. + pipe_tune3 ) )
178
179	!! 1.1.5 calculate provisional tree crown area for per individual tree
180
181	! equation: CrownArea=pipe_tune1 * Diameter^{1.6}
182	cn_ind(i,j) = pipe_tune1 * MIN( maxdia(j), dia(i) ) ** pipe_tune_exp_coeff
183
184	!! 1.1.6 If total tree crown area for this tree PFT exceeds its veget_cov_max, tree density is recalculated.
185
186	IF ( cn_ind(i,j) * ind(i,j) .GT. 1.002* veget_cov_max(i,j) ) THEN
187
188	ind(i,j) = veget_cov_max(i,j) / cn_ind(i,j)
189
190	ELSE
191
192	ind(i,j) = ( veget_cov_max(i,j) / &
193	& ( pipe_tune1 * (woodmass(i)/(pipe_densitypi/4.pipe_tune2)) &
194	& **(pipe_tune_exp_coeff/(2.+pipe_tune3)) ) ) &
195	& ** (1./(1.-(pipe_tune_exp_coeff/(2.+pipe_tune3))))
196
197
198	woodmass_ind(i) = woodmass(i) / ind(i,j)
199
200	dia(i) = ( woodmass_ind(i) / ( pipe_density * pi/4. * pipe_tune2 ) ) ** &
201	( un / ( 2. + pipe_tune3 ) )
202
203	! final crown area
204	cn_ind(i,j) = pipe_tune1 * MIN( maxdia(j), dia(i) ) ** pipe_tune_exp_coeff
205
206	ENDIF
207
208	ELSE !woodmas=0 => impose some value
209
210	dia(:) = maxdia(j)
211
212	cn_ind(i,j) = pipe_tune1 * MIN( maxdia(j), dia(i) ) ** pipe_tune_exp_coeff
213
214	ENDIF ! IF ( woodmass(i) .GT. min_stomate )
215
216	ENDIF ! veget_cov_max .GT. 0.
217
218	ENDDO ! loop over grid points
219
220	ELSE !grass
221
222	!
223	!! 1.2 grasses: crown area always set to 1m**2
224	!
225
226	WHERE ( veget_cov_max(:,j) .GT. zero )
227	cn_ind(:,j) = un
228	ENDWHERE
229
230	ENDIF !IF ( is_tree(j) )
231
232	!
233	!! 2 density of individuals
234	!
235
236	WHERE ( veget_cov_max(:,j) .GT. zero )
237
238	ind(:,j) = veget_cov_max(:,j) / cn_ind(:,j)
239
240	ELSEWHERE
241
242	ind(:,j) = zero
243
244	ENDWHERE
245
246	ENDIF ! IF ( ( .NOT. ok_dgvm .AND. lpj_gap_const_mort ) .OR. ( .NOT. natural(j) ) )
247
248	ENDDO ! loop over PFTs
249
250	!
251	!!? it's better to move the code for first call at the beginning of the module.
252	!! 2 If it's the first call for this module,
253	!
254
255	IF (( firstcall_prescribe ) .AND. (TRIM(stom_restname_in) == 'NONE')) THEN
256
257	IF (printlev >= 2) THEN
258	WRITE(numout,*) 'prescribe:'
259	! impose some biomass if zero and PFT prescribed
260	WRITE(numout,*) ' > Imposing initial biomass for prescribed trees, '// &
261	'initial reserve mass for prescribed grasses.'
262	WRITE(numout,*) ' > Declaring prescribed PFTs present.'
263	END IF
264
265	DO j = 2,nvm ! loop over PFTs
266	DO i = 1, npts ! loop over grid points
267
268	! is vegetation static or PFT agricultural?
269	! Static vegetation or agricultural PFT
270	IF ( ( .NOT. ok_dgvm ) .OR. &
271	( ( .NOT. natural(j) .OR. pasture(j)) .AND. ( veget_cov_max(i,j) .GT. min_stomate ) ) ) THEN
272
273	!
274	!! 2.1 if tree biomass is extremely small, prescribe the biomass by assuming they have sapling biomass, which is a constant in the model.
275	!! then set all the leaf age as 1.
276	!
277	! if tree PFT and biomass too small, prescribe the biomass to a value.
278	IF ( is_tree(j) .AND. &
279	( veget_cov_max(i,j) .GT. min_stomate ) .AND. &
280	( SUM( biomass(i,j,:,icarbon) ) .LE. min_stomate ) ) THEN
281	!!? here the code is redundant, as "veget_cov_max(i,j) .GT. min_stomate" is already met in the above if condition.
282	IF (veget_cov_max(i,j) .GT. min_stomate) THEN
283	biomass(i,j,:,:) = (bm_sapl_rescale * bm_sapl(j,:,:) * ind(i,j)) / veget_cov_max(i,j)
284	ELSE
285	biomass(i,j,:,:) = zero
286	ENDIF
287
288	! set leaf age classes
289	!! yidi
290	!! no bm_sapl added to bm_phytomer PHYbm FFBbm, bm_sapl to bm_sapling first
291	!! IF (ok_oilpalm) THEN
292	!! IF (is_oilpalm(j)) THEN
293	!! IF (biomass(i,j,:,:) .GT. min_stomate) THEN
294	!! bm_phytomer(:,j,1)= bm_phytomer(:,j,1) + &
295	!! (biomass(j,isapabove,icarbon) *
296	!! PHYbm(:,j) / biomass(:,j,isapabove,icarbon)
297	!! PHYbm(:,j)=PHYbm(:,j)+ &
298	!! biomass(j,isapabove,icarbon)* PHYbm(:,j) / biomass(:,j,isapabove,icarbon)
299	!! ENDIF
300	!! 1ENDIF
301	!!yidi
302	leaf_frac(i,j,:) = zero
303	leaf_frac(i,j,1) = un
304
305	! set time since last beginning of growing season
306	when_growthinit(i,j) = large_value
307
308	! seasonal trees: no leaves at beginning
309
310	IF ( pheno_model(j) .NE. 'none' ) THEN
311
312	biomass(i,j,ileaf,icarbon) = zero
313	leaf_frac(i,j,1) = zero
314
315	ENDIF
316
317	co2_to_bm(i,j) = co2_to_bm(i,j) + ( SUM(biomass(i,j,:,icarbon)) / dt )
318	ENDIF
319
320	!
321	!! 2.2 for grasses, set only the carbon reserve pool to "sapling" carbon reserve pool.
322	!! and set all leaf age to 1.
323
324	IF ( ( .NOT. is_tree(j) ) .AND. &
325	( veget_cov_max(i,j) .GT. min_stomate ) .AND. &
326	( SUM( biomass(i,j,:,icarbon) ) .LE. min_stomate ) ) THEN
327
328	biomass(i,j,icarbres,:) = bm_sapl(j,icarbres,:) * ind(i,j) / veget_cov_max(i,j)
329
330	! set leaf age classes
331	leaf_frac(i,j,:) = zero
332	leaf_frac(i,j,1) = un
333
334	! set time since last beginning of growing season
335	when_growthinit(i,j) = large_value
336
337	co2_to_bm(i,j) = co2_to_bm(i,j) + ( biomass(i,j,icarbres,icarbon) / dt )
338	ENDIF
339
340	!
341	!! 2.3 declare all PFTs with positive veget_cov_max as present everywhere in that grid box
342	!
343
344	IF ( veget_cov_max(i,j) .GT. min_stomate ) THEN
345	PFTpresent(i,j) = .TRUE.
346	everywhere(i,j) = un
347	ENDIF
348
349	ENDIF ! not ok_dgvm or agricultural
350
351	ENDDO ! loop over grid points
352	ENDDO ! loop over PFTs
353
354
355	ENDIF
356
357	firstcall_prescribe = .FALSE.
358
359	END SUBROUTINE prescribe
360
361	END MODULE stomate_prescribe

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

Context Navigation

source: branches/publications/ORCHIDEE-MICT-OP-r6850/src_stomate/stomate_prescribe.f90 @ 7108

Download in other formats: