Context Navigation

source: branches/ORCHIDEE_EXT/ORCHIDEE/src_stomate/stomate_season.f90 @ 135

Last change on this file since 135 was 135, checked in by didier.solyga, 13 years ago
Change herbivory model : herbivores is pft-dependent now (with N.Vuichard & N.Viovy)
File size: 31.8 KB

Line
1	! Calculate long-term meteorological parameters from daily temperatures
2	! and precipitations (essentially for phenology)
3	!
4	! $Header: /home/ssipsl/CVSREP/ORCHIDEE/src_stomate/stomate_season.f90,v 1.15 2010/04/06 16:06:34 ssipsl Exp $
5	! IPSL (2006)
6	! This software is governed by the CeCILL licence see ORCHIDEE/ORCHIDEE_CeCILL.LIC
7	!
8	MODULE stomate_season
9
10	! modules used:
11
12	USE ioipsl
13	USE stomate_data
14	USE constantes
15	USE pft_parameters
16
17	IMPLICIT NONE
18
19	! private & public routines
20
21	PRIVATE
22	PUBLIC season,season_clear
23
24	! first call
25	LOGICAL, SAVE :: firstcall = .TRUE.
26
27	CONTAINS
28
29
30	SUBROUTINE season_clear
31	firstcall=.TRUE.
32	END SUBROUTINE season_clear
33
34	SUBROUTINE season (npts, dt, EndOfYear, veget, veget_max, &
35	moiavail_daily, t2m_daily, tsoil_daily, soilhum_daily, &
36	precip_daily, npp_daily, biomass, turnover_daily, gpp_daily, when_growthinit, &
37	maxmoiavail_lastyear, maxmoiavail_thisyear, &
38	minmoiavail_lastyear, minmoiavail_thisyear, &
39	maxgppweek_lastyear, maxgppweek_thisyear, &
40	gdd0_lastyear, gdd0_thisyear, &
41	precip_lastyear, precip_thisyear, &
42	lm_lastyearmax, lm_thisyearmax, &
43	maxfpc_lastyear, maxfpc_thisyear, &
44	moiavail_month, moiavail_week, t2m_longterm, tlong_ref, t2m_month, t2m_week, &
45	tsoil_month, soilhum_month, &
46	npp_longterm, turnover_longterm, gpp_week, &
47	gdd_m5_dormance, gdd_midwinter, ncd_dormance, ngd_minus5, time_lowgpp, &
48	time_hum_min, hum_min_dormance, herbivores)
49
50	!
51	! 0 declarations
52	!
53
54	! 0.1 input
55
56	! Domain size
57	INTEGER(i_std), INTENT(in) :: npts
58	! time step in days
59	REAL(r_std), INTENT(in) :: dt
60	! update yearly variables?
61	LOGICAL, INTENT(in) :: EndOfYear
62	! coverage fraction of a PFT. Here: fraction of total ground.
63	REAL(r_std), DIMENSION(npts,nvm), INTENT(in) :: veget
64	! "maximal" coverage fraction of a PFT (for LAI -> infinity)
65	! Here: fraction of total ground.
66	REAL(r_std), DIMENSION(npts,nvm), INTENT(in) :: veget_max
67	! Daily moisture availability
68	REAL(r_std), DIMENSION(npts,nvm), INTENT(in) :: moiavail_daily
69	! Daily 2 meter temperature (K)
70	REAL(r_std), DIMENSION(npts), INTENT(in) :: t2m_daily
71	! Daily soil temperature (K)
72	REAL(r_std), DIMENSION(npts,nbdl), INTENT(in) :: tsoil_daily
73	! Daily soil humidity
74	REAL(r_std), DIMENSION(npts,nbdl), INTENT(in) :: soilhum_daily
75	! Daily mean precipitation (mm/day)
76	REAL(r_std), DIMENSION(npts), INTENT(in) :: precip_daily
77	! daily net primary productivity (gC/m**2/day)
78	REAL(r_std), DIMENSION(npts,nvm), INTENT(in) :: npp_daily
79	! biomass (gC/(m**2 of ground))
80	REAL(r_std), DIMENSION(npts,nvm,nparts), INTENT(in) :: biomass
81	! Turnover rates (gC/m**2/day)
82	REAL(r_std), DIMENSION(npts,nvm,nparts), INTENT(in) :: turnover_daily
83	! daily gross primary productivity (Here: gC/(m**2 of total ground)/day)
84	REAL(r_std), DIMENSION(npts,nvm), INTENT(in) :: gpp_daily
85	! how many days ago was the beginning of the growing season
86	REAL(r_std), DIMENSION(npts,nvm), INTENT(in) :: when_growthinit
87
88	! 0.2 modified fields
89
90	! last year's maximum moisture availability
91	REAL(r_std), DIMENSION(npts,nvm), INTENT(inout) :: maxmoiavail_lastyear
92	! this year's maximum moisture availability
93	REAL(r_std), DIMENSION(npts,nvm), INTENT(inout) :: maxmoiavail_thisyear
94	! last year's minimum moisture availability
95	REAL(r_std), DIMENSION(npts,nvm), INTENT(inout) :: minmoiavail_lastyear
96	! this year's minimum moisture availability
97	REAL(r_std), DIMENSION(npts,nvm), INTENT(inout) :: minmoiavail_thisyear
98	! last year's maximum weekly GPP
99	REAL(r_std), DIMENSION(npts,nvm), INTENT(inout) :: maxgppweek_lastyear
100	! this year's maximum weekly GPP
101	REAL(r_std), DIMENSION(npts,nvm), INTENT(inout) :: maxgppweek_thisyear
102	! last year's annual GDD0
103	REAL(r_std), DIMENSION(npts), INTENT(inout) :: gdd0_lastyear
104	! this year's annual GDD0
105	REAL(r_std), DIMENSION(npts), INTENT(inout) :: gdd0_thisyear
106	! last year's annual precipitation (mm/year)
107	REAL(r_std), DIMENSION(npts), INTENT(inout) :: precip_lastyear
108	! this year's annual precipitation (mm/year)
109	REAL(r_std), DIMENSION(npts), INTENT(inout) :: precip_thisyear
110	! last year's maximum leaf mass, for each PFT (gC/m**2)
111	REAL(r_std), DIMENSION(npts,nvm), INTENT(inout) :: lm_lastyearmax
112	! this year's maximum leaf mass, for each PFT (gC/m**2)
113	REAL(r_std), DIMENSION(npts,nvm), INTENT(inout) :: lm_thisyearmax
114	! last year's maximum fpc for each PFT, on ground
115	REAL(r_std), DIMENSION(npts,nvm), INTENT(inout) :: maxfpc_lastyear
116	! this year's maximum fpc for each PFT, on ground
117	REAL(r_std), DIMENSION(npts,nvm), INTENT(inout) :: maxfpc_thisyear
118	! "monthly" moisture availability
119	REAL(r_std), DIMENSION(npts,nvm), INTENT(inout) :: moiavail_month
120	! "weekly" moisture availability
121	REAL(r_std), DIMENSION(npts,nvm), INTENT(inout) :: moiavail_week
122	! "long term" 2-meter temperatures (K)
123	REAL(r_std), DIMENSION(npts), INTENT(inout) :: t2m_longterm
124	! "long term" refernce 2-meter temperatures (K)
125	REAL(r_std), DIMENSION(npts), INTENT(inout) :: tlong_ref
126	! "monthly" 2-meter temperatures (K)
127	REAL(r_std), DIMENSION(npts), INTENT(inout) :: t2m_month
128	! "weekly" 2-meter temperatures (K)
129	REAL(r_std), DIMENSION(npts), INTENT(inout) :: t2m_week
130	! "monthly" soil temperatures (K)
131	REAL(r_std), DIMENSION(npts,nbdl), INTENT(inout) :: tsoil_month
132	! "monthly" soil humidity
133	REAL(r_std), DIMENSION(npts,nbdl), INTENT(inout) :: soilhum_month
134	! "long term" net primary productivity (gC/m**2/year)
135	REAL(r_std), DIMENSION(npts,nvm), INTENT(inout) :: npp_longterm
136	! "long term" turnover rate (gC/m**2/year)
137	REAL(r_std), DIMENSION(npts,nvm,nparts), INTENT(inout) :: turnover_longterm
138	! "weekly" GPP (gC/day/(m**2 covered)
139	REAL(r_std), DIMENSION(npts,nvm), INTENT(inout) :: gpp_week
140	! growing degree days, threshold -5 deg. C
141	REAL(r_std), DIMENSION(npts,nvm), INTENT(inout) :: gdd_m5_dormance
142	! growing degree days since midwinter
143	REAL(r_std), DIMENSION(npts,nvm), INTENT(inout) :: gdd_midwinter
144	! number of chilling days since leaves were lost
145	REAL(r_std), DIMENSION(npts,nvm), INTENT(inout) :: ncd_dormance
146	! number of growing days
147	REAL(r_std), DIMENSION(npts,nvm), INTENT(inout) :: ngd_minus5
148	! duration of dormance (d)
149	REAL(r_std), DIMENSION(npts,nvm), INTENT(inout) :: time_lowgpp
150	! time elapsed since strongest moisture availability (d)
151	REAL(r_std), DIMENSION(npts,nvm), INTENT(inout) :: time_hum_min
152	! minimum moisture during dormance
153	REAL(r_std), DIMENSION(npts,nvm), INTENT(inout) :: hum_min_dormance
154
155	! 0.3 output (diagnostic)
156
157	! time constant of probability of a leaf to be eaten by a herbivore (days)
158	REAL(r_std), DIMENSION(npts,nvm), INTENT(out) :: herbivores
159
160	! 0.4 local
161
162	! indices
163	INTEGER(i_std) :: j
164	! maximum ncd (d) (to avoid floating point underflows)
165	REAL(r_std) :: ncd_max
166	! sum of natural fpcs
167	REAL(r_std), DIMENSION(npts) :: sumfpc_nat
168	!!$ 01/03/2011 NVuichard, NViovy and DS : correcting herbivores : now weighttot
169	!!$ nlflong_nat, green_age are pft-dependants
170	! weights
171	REAL(r_std), DIMENSION(npts,nvm) :: weighttot
172	! natural long-term leaf NPP ( gC/m**2/year)
173	REAL(r_std), DIMENSION(npts,nvm) :: nlflong_nat
174	! residence time of green tissue (years)
175	REAL(r_std), DIMENSION(npts,nvm) :: green_age
176
177	! ds 04/03 Old formulation for herbivores
178	!!$ ! weights
179	!!$ REAL(r_std), DIMENSION(npts) :: weighttot
180	!!$ ! natural long-term leaf NPP ( gC/m**2/year)
181	!!$ REAL(r_std), DIMENSION(npts) :: nlflong_nat
182	!!$ ! residence time of green tissue (years)
183	!!$ REAL(r_std), DIMENSION(npts) :: green_age
184
185	! herbivore consumption (gC/m**2/day)
186	REAL(r_std), DIMENSION(npts) :: consumption
187
188	! =========================================================================
189
190	IF (bavard.GE.3) WRITE(numout,*) 'Entering season'
191
192	!
193	! 1 Initializations
194	!
195	ncd_max = ncd_max_year * one_year
196
197	IF ( firstcall ) THEN
198
199	!
200	! 1.1 messages
201	!
202
203	IF ( bavard .GE. 1 ) THEN
204
205	WRITE(numout,*) 'season: '
206
207	WRITE(numout,*) ' > rapport maximal GPP/GGP_max pour dormance: ',gppfrac_dormance
208
209	WRITE(numout,*) ' > maximum possible ncd (d): ',ncd_max
210
211	WRITE(numout,*) ' > herbivore consumption C (gC/m2/day) as a function of NPP (gC/m2/d):'
212	WRITE(numout,) ' C=',hvc1,' NPP^',hvc2
213	WRITE(numout,) ' > for herbivores, suppose that ',leaf_frac_hvc100., &
214	'% of NPP is allocated to leaves'
215
216
217	ENDIF
218
219	!
220	! 1.2 Check whether longer-term meteorological parameters are initialized
221	! to zero
222	!
223
224	! 1.2.1 moisture availabilities
225
226	! 1.2.1.1 "monthly"
227	!MM PAS PARALLELISE!!
228	IF ( ABS( SUM( moiavail_month(:,2:nvm) ) ) .LT. min_stomate ) THEN
229
230	! in this case, set them it today's moisture availability
231	WRITE(numout,*) 'Warning! We have to initialize the ''monthly'' moisture availabilities.'
232	moiavail_month(:,:) = moiavail_daily(:,:)
233
234	ENDIF
235
236	! 1.2.1.2 "weekly"
237
238	IF ( ABS( SUM( moiavail_week(:,2:nvm) ) ) .LT. min_stomate ) THEN
239
240	! in this case, set them it today's moisture availability
241	WRITE(numout,*) 'Warning! We have to initialize the ''weekly'' moisture availabilities.'
242	moiavail_week(:,:) = moiavail_daily(:,:)
243
244	ENDIF
245
246	! 1.2.2 2-meter temperatures
247
248	! 1.2.2.1 "long term"
249
250	IF ( ABS( SUM( t2m_longterm(:) ) ) .LT. min_stomate ) THEN
251
252	! in this case, set them to today's temperature
253	WRITE(numout,*) 'Warning! We have to initialize the ''long term'' 2m temperatures.'
254	t2m_longterm(:) = t2m_daily(:)
255
256	ENDIF
257
258	! 1.2.2.2 "monthly"
259
260	IF ( ABS( SUM( t2m_month(:) ) ) .LT. min_stomate ) THEN
261
262	! in this case, set them to today's temperature
263	WRITE(numout,*) 'Warning! We have to initialize the ''monthly'' 2m temperatures.'
264	t2m_month(:) = t2m_daily(:)
265
266	ENDIF
267
268	! 1.2.2.3 "weekly"
269
270	IF ( ABS( SUM( t2m_week(:) ) ) .LT. min_stomate ) THEN
271
272	! in this case, set them to today's temperature
273	WRITE(numout,*) 'Warning! We have to initialize the ''weekly'' 2m temperatures.'
274	t2m_week(:) = t2m_daily(:)
275
276	ENDIF
277
278	! 1.2.3 "monthly" soil temperatures
279	!MM PAS PARALLELISE!!
280	IF ( ABS( SUM( tsoil_month(:,:) ) ) .LT. min_stomate ) THEN
281
282	! in this case, set them to today's temperature
283	WRITE(numout,*) 'Warning!'// &
284	' We have to initialize the ''monthly'' soil temperatures.'
285	tsoil_month(:,:) = tsoil_daily(:,:)
286
287	ENDIF
288
289	! 1.2.4 "monthly" soil humidity
290	IF ( ABS( SUM( soilhum_month(:,:) ) ) .LT. min_stomate ) THEN
291
292	! in this case, set them to today's humidity
293	WRITE(numout,*) 'Warning!'// &
294	' We have to initialize the ''monthly'' soil humidity.'
295	soilhum_month(:,:) = soilhum_daily(:,:)
296
297	ENDIF
298
299	! 1.2.5 growing degree days, threshold -5 deg C
300	IF ( ABS( SUM( gdd_m5_dormance(:,2:nvm) ) ) .LT. min_stomate ) THEN
301	WRITE(numout,*) 'Warning! Growing degree days (-5 deg) are initialized to ''undef''.'
302	gdd_m5_dormance(:,:) = undef
303	ENDIF
304
305	! 1.2.6 growing degree days since midwinter
306	IF ( ABS( SUM( gdd_midwinter(:,2:nvm) ) ) .LT. min_stomate ) THEN
307	WRITE(numout,*) 'Warning! Growing degree days since midwinter' // &
308	' are initialized to ''undef''.'
309	gdd_midwinter(:,:) = undef
310	ENDIF
311
312	! 1.2.7 number of chilling days since leaves were lost
313	IF ( ABS( SUM( ncd_dormance(:,2:nvm) ) ) .LT. min_stomate ) THEN
314	WRITE(numout,*) 'Warning! Number of chilling days is initialized to ''undef''.'
315	ncd_dormance(:,:) = undef
316	ENDIF
317
318	! 1.2.8 number of growing days, threshold -5 deg C
319	IF ( ABS( SUM( ngd_minus5(:,2:nvm) ) ) .LT. min_stomate ) THEN
320	WRITE(numout,*) 'Warning! Number of growing days (-5 deg) is initialized to 0.'
321	ENDIF
322
323	! 1.2.9 "long term" npp
324	IF ( ABS( SUM( npp_longterm(:,2:nvm) ) ) .LT. min_stomate ) THEN
325	WRITE(numout,*) 'Warning! Long term NPP is initialized to 0.'
326	ENDIF
327
328	! 1.2.10 "long term" turnover
329	IF ( ABS( SUM( turnover_longterm(:,2:nvm,:) ) ) .LT. min_stomate ) THEN
330	WRITE(numout,*) 'Warning! Long term turnover is initialized to 0.'
331	ENDIF
332
333	! 1.2.11 "weekly" GPP
334	IF ( ABS( SUM( gpp_week(:,2:nvm) ) ) .LT. min_stomate ) THEN
335	WRITE(numout,*) 'Warning! Weekly GPP is initialized to 0.'
336	ENDIF
337
338	! 1.2.12 minimum moisture availabilities
339	IF ( ABS( SUM( minmoiavail_thisyear(:,2:nvm) ) ) .LT. min_stomate ) THEN
340
341	! in this case, set them to a very high value
342	WRITE(numout,*) 'Warning! We have to initialize this year''s minimum '// &
343	'moisture availabilities.'
344	minmoiavail_thisyear(:,:) = large_value
345
346	ENDIF
347
348	!
349	! 1.3 reset flag
350	!
351
352	firstcall = .FALSE.
353
354	ENDIF
355
356	!
357	! 2 moisture availabilities
358	!
359
360	!
361	! 2.1 "monthly"
362	!
363
364	moiavail_month = ( moiavail_month * ( tau_hum_month - dt ) + &
365	moiavail_daily * dt ) / tau_hum_month
366
367	DO j = 2,nvm
368	WHERE ( ABS(moiavail_month(:,j)) .LT. EPSILON(zero) )
369	moiavail_month(:,j) = zero
370	ENDWHERE
371	ENDDO
372
373	!
374	! 2.2 "weekly"
375	!
376
377	moiavail_week = ( moiavail_week * ( tau_hum_week - dt ) + &
378	moiavail_daily * dt ) / tau_hum_week
379
380	DO j = 2,nvm
381	WHERE ( ABS(moiavail_week(:,j)) .LT. EPSILON(zero) )
382	moiavail_week(:,j) = zero
383	ENDWHERE
384	ENDDO
385
386	!
387	! 3 2-meter temperatures
388	!
389
390	!
391	! 3.1 "long term"
392	!
393
394	t2m_longterm = ( t2m_longterm * ( tau_longterm - dt ) + &
395	t2m_daily * dt ) / tau_longterm
396
397	WHERE ( ABS(t2m_longterm(:)) .LT. EPSILON(zero) )
398	t2m_longterm(:) = zero
399	ENDWHERE
400
401	CALL histwrite (hist_id_stomate, 'T2M_LONGTERM', itime, &
402	t2m_longterm, npts, hori_index)
403	!
404	! 3.2 "long term reference"
405	! This temperature is used for recalculating PFT-specific parameters such as
406	! critical photosynthesis temperatures of critical GDDs for phenology. This
407	! means that if the reference temperature varies, the PFTs adapt to them.
408	! Therefore the reference temperature can vary only if the vegetation is not
409	! static.
410	!
411
412	!!$ IF (control%ok_dgvm) THEN
413	tlong_ref(:) = MAX( tlong_ref_min, MIN( tlong_ref_max, t2m_longterm(:) ) )
414	!!$ ENDIF
415	!
416	! 3.3 "monthly"
417	!
418
419	t2m_month = ( t2m_month * ( tau_t2m_month - dt ) + &
420	t2m_daily * dt ) / tau_t2m_month
421
422	WHERE ( ABS(t2m_month(:)) .LT. EPSILON(zero) )
423	t2m_month(:) = zero
424	ENDWHERE
425
426	!
427	! 3.4 "weekly"
428	!
429
430	t2m_week = ( t2m_week * ( tau_t2m_week - dt ) + &
431	t2m_daily * dt ) / tau_t2m_week
432
433	WHERE ( ABS(t2m_week(:)) .LT. EPSILON(zero) )
434	t2m_week(:) = zero
435	ENDWHERE
436
437	!
438	! 4 ''monthly'' soil temperatures
439	!
440
441	tsoil_month = ( tsoil_month * ( tau_tsoil_month - dt ) + &
442	tsoil_daily(:,:) * dt ) / tau_tsoil_month
443
444	WHERE ( ABS(tsoil_month(:,:)) .LT. EPSILON(zero) )
445	tsoil_month(:,:) = zero
446	ENDWHERE
447
448	!
449	! 5 ''monthly'' soil humidity
450	!
451
452	soilhum_month = ( soilhum_month * ( tau_soilhum_month - dt ) + &
453	soilhum_daily * dt ) / tau_soilhum_month
454
455	WHERE ( ABS(soilhum_month(:,:)) .LT. EPSILON(zero) )
456	soilhum_month(:,:) = zero
457	ENDWHERE
458
459	!
460	! 6 dormance (d)
461	! when gpp is low, increase dormance time. Otherwise, set it to zero.
462	! NV: special case (3rd condition): plant is accumulating carbohydrates
463	! and does never use them. In this case, we allow the plant to
464	! detect a beginning of the growing season by declaring it dormant
465	!
466	!NVMODIF
467	DO j = 2,nvm
468	WHERE ( ( gpp_week(:,j) .LT. min_gpp_allowed ) .OR. &
469	( gpp_week(:,j) .LT. gppfrac_dormance * maxgppweek_lastyear(:,j) ) .OR. &
470	( ( when_growthinit(:,j) .GT. 2.*one_year ) .AND. &
471	( biomass(:,j,icarbres) .GT. biomass(:,j,ileaf)*4. ) ) )
472	! WHERE ( ( gpp_week(:,j) .EQ. zero ) .OR. &
473	! ( gpp_week(:,j) .LT. gppfrac_dormance * maxgppweek_lastyear(:,j) ) .OR. &
474	! ( ( when_growthinit(:,j) .GT. 2.*one_year ) .AND. &
475	! ( biomass(:,j,icarbres) .GT. biomass(:,j,ileaf)*4. ) ) )
476
477	time_lowgpp(:,j) = time_lowgpp(:,j) + dt
478
479	ELSEWHERE
480
481	time_lowgpp(:,j) = zero
482
483	ENDWHERE
484	ENDDO
485
486	!
487	! 7 growing degree days, threshold -5 deg C
488	!
489
490	DO j = 2,nvm
491
492	! only for PFTs for which critical gdd is defined
493	! gdd_m5_dormance is set to 0 at the end of the growing season. It is set to undef
494	! at the beginning of the growing season.
495
496	IF ( ALL(pheno_gdd_crit(j,:) .NE. undef) ) THEN
497
498	!
499	! 7.1 set to zero if undef and no gpp
500	!
501
502	WHERE ( ( time_lowgpp(:,j) .GT. zero ) .AND. ( gdd_m5_dormance(:,j) .EQ. undef ) )
503
504	gdd_m5_dormance(:,j) = zero
505
506	ENDWHERE
507
508	!
509	! 7.2 set to undef if there is gpp
510	!
511
512	WHERE ( time_lowgpp(:,j) .EQ. zero )
513
514	gdd_m5_dormance(:,j) = undef
515
516	ENDWHERE
517
518	!
519	! 7.3 normal update where gdd_m5_dormance is defined
520	!
521
522	WHERE ( ( t2m_daily(:) .GT. (ZeroCelsius-5.) ) .AND. &
523	( gdd_m5_dormance(:,j) .NE. undef ) )
524	gdd_m5_dormance(:,j) = gdd_m5_dormance(:,j) + &
525	dt * ( t2m_daily(:) - (ZeroCelsius-5.) )
526	ENDWHERE
527
528	WHERE ( gdd_m5_dormance(:,j) .NE. undef )
529	gdd_m5_dormance(:,j) = gdd_m5_dormance(:,j) * &
530	( tau_gdd - dt ) / tau_gdd
531	ENDWHERE
532
533	ENDIF
534
535	ENDDO
536
537	DO j = 2,nvm
538	WHERE ( ABS(gdd_m5_dormance(:,j)) .LT. EPSILON(zero) )
539	gdd_m5_dormance(:,j) = zero
540	ENDWHERE
541	ENDDO
542
543	!
544	! 8 growing degree days since midwinter
545	!
546
547	DO j = 2,nvm
548
549	! only for PFTs for which ncdgdd_crittemp is defined
550
551	IF ( ncdgdd_temp(j) .NE. undef ) THEN
552
553	!
554	! 8.1 set to 0 if undef and if we detect "midwinter"
555	!
556
557	WHERE ( ( gdd_midwinter(:,j) .EQ. undef ) .AND. &
558	( t2m_month(:) .LT. t2m_week(:) ) .AND. &
559	( t2m_month(:) .LT. t2m_longterm(:) ) )
560
561	gdd_midwinter(:,j) = zero
562
563	ENDWHERE
564
565	!
566	! 8.2 set to undef if we detect "midsummer"
567	!
568
569	WHERE ( ( t2m_month(:) .GT. t2m_week(:) ) .AND. &
570	( t2m_month(:) .GT. t2m_longterm(:) ) )
571
572	gdd_midwinter(:,j) = undef
573
574	ENDWHERE
575
576	!
577	! 8.3 normal update
578	!
579
580	WHERE ( ( gdd_midwinter(:,j) .NE. undef ) .AND. &
581	( t2m_daily(:) .GT. ncdgdd_temp(j)+ZeroCelsius ) )
582
583	gdd_midwinter(:,j) = &
584	gdd_midwinter(:,j) + &
585	dt * ( t2m_daily(:) - ( ncdgdd_temp(j)+ZeroCelsius ) )
586
587	ENDWHERE
588
589	ENDIF
590
591	ENDDO
592
593	!
594	! 9 number of chilling days since leaves were lost
595	!
596
597	DO j = 2,nvm
598
599	IF ( ncdgdd_temp(j) .NE. undef ) THEN
600
601	!
602	! 9.1 set to zero if undef and no gpp
603	!
604
605	WHERE ( ( time_lowgpp(:,j) .GT. zero ) .AND. ( ncd_dormance(:,j) .EQ. undef ) )
606
607	ncd_dormance(:,j) = zero
608
609	ENDWHERE
610
611	!
612	! 9.2 set to undef if there is gpp
613	!
614
615	WHERE ( time_lowgpp(:,j) .EQ. zero )
616
617	ncd_dormance(:,j) = undef
618
619	ENDWHERE
620
621	!
622	! 9.3 normal update where ncd_dormance is defined
623	!
624
625	WHERE ( ( ncd_dormance(:,j) .NE. undef ) .AND. &
626	( t2m_daily(:) .LE. ncdgdd_temp(j)+ZeroCelsius ) )
627
628	ncd_dormance(:,j) = MIN( ncd_dormance(:,j) + dt, ncd_max )
629
630	ENDWHERE
631
632	ENDIF
633
634	ENDDO
635
636	!
637	! 10 number of growing days, threshold -5 deg C
638	!
639
640	DO j = 2,nvm
641
642	!
643	! 10.1 Where there is GPP, set ngd to 0
644	! This means that we only take into account ngds when the leaves are off
645	!
646
647	WHERE ( time_lowgpp(:,j) .LT. min_stomate )
648	ngd_minus5(:,j) = zero
649	ENDWHERE
650
651	!
652	! 10.2 normal update
653	!
654
655	WHERE ( t2m_daily(:) .GT. (ZeroCelsius-5.) )
656	ngd_minus5(:,j) = ngd_minus5(:,j) + dt
657	ENDWHERE
658
659	ngd_minus5(:,j) = ngd_minus5(:,j) * ( tau_ngd - dt ) / tau_ngd
660
661	ENDDO
662
663	DO j = 2,nvm
664	WHERE ( ABS(ngd_minus5(:,j)) .LT. EPSILON(zero) )
665	ngd_minus5(:,j) = zero
666	ENDWHERE
667	ENDDO
668	!
669	! 11 minimum humidity since dormance began and time elapsed since this minimum
670	!
671
672	DO j = 2,nvm
673
674	IF ( hum_min_time(j) .NE. undef ) THEN
675
676	!
677	! 11.1 initialize if undef and no gpp
678	!
679
680	WHERE ( ( time_lowgpp(:,j) .GT. zero ) .AND. &
681	( ( time_hum_min(:,j) .EQ. undef ) .OR. ( hum_min_dormance(:,j) .EQ. undef ) ) )
682
683	time_hum_min(:,j) = zero
684	hum_min_dormance(:,j) = moiavail_month(:,j)
685
686	ENDWHERE
687
688	!
689	! 11.2 set to undef where there is gpp
690	!
691
692	WHERE ( time_lowgpp(:,j) .EQ. zero )
693
694	time_hum_min(:,j) = undef
695	hum_min_dormance(:,j) = undef
696
697	ENDWHERE
698
699	!
700	! 11.3 normal update where time_hum_min and hum_min_dormance are defined
701	!
702
703	! 11.3.1 increase time counter
704
705	WHERE ( ( time_hum_min(:,j) .NE. undef ) .AND. &
706	( hum_min_dormance(:,j) .NE. undef ) )
707
708	time_hum_min(:,j) = time_hum_min(:,j) + dt
709
710	ENDWHERE
711
712	! 11.3.2 set time to zero if minimum is reached
713
714	WHERE ( ( time_hum_min(:,j) .NE. undef ) .AND. &
715	( hum_min_dormance(:,j) .NE. undef ) .AND. &
716	( moiavail_month(:,j) .LE. hum_min_dormance(:,j) ) )
717
718	hum_min_dormance(:,j) = moiavail_month(:,j)
719	time_hum_min(:,j) = zero
720
721	ENDWHERE
722
723	ENDIF
724
725	ENDDO
726
727	!
728	! 12 "long term" NPP. npp_daily in gC/m2/day, npp_longterm in gC/m2/year.
729	!
730
731	npp_longterm = ( npp_longterm * ( tau_longterm - dt ) + &
732	(npp_dailyone_year) dt ) / &
733	tau_longterm
734
735	DO j = 2,nvm
736	WHERE ( ABS(npp_longterm(:,j)) .LT. EPSILON(zero) )
737	npp_longterm(:,j) = zero
738	ENDWHERE
739	ENDDO
740
741	!
742	! 13 "long term" turnover rates, in gC/m**2/year.
743	!
744
745	turnover_longterm = ( turnover_longterm * ( tau_longterm - dt ) + &
746	(turnover_dailyone_year) dt ) / &
747	tau_longterm
748
749	DO j = 2,nvm
750	WHERE ( ABS(turnover_longterm(:,j,:)) .LT. EPSILON(zero) )
751	turnover_longterm(:,j,:) = zero
752	ENDWHERE
753	ENDDO
754
755	!
756	! 14 "weekly" GPP, in gC/(m**2 covered)/day (!)
757	! i.e. divide daily gpp (in gC/m**2 of total ground/day) by vegetation fraction
758	! (m2 covered/m2 of total ground)
759	!
760
761	WHERE ( veget_max .GT. zero )
762
763	gpp_week = ( gpp_week * ( tau_gpp_week - dt ) + &
764	gpp_daily * dt ) / tau_gpp_week
765
766	ELSEWHERE
767
768	gpp_week = zero
769
770	ENDWHERE
771
772	DO j = 2,nvm
773	WHERE ( ABS(gpp_week(:,j)) .LT. EPSILON(zero) )
774	gpp_week(:,j) = zero
775	ENDWHERE
776	ENDDO
777
778	!
779	! 15 maximum and minimum moisture availabilities
780	!
781
782	WHERE ( moiavail_daily .GT. maxmoiavail_thisyear )
783	maxmoiavail_thisyear = moiavail_daily
784	ENDWHERE
785
786	WHERE ( moiavail_daily .LT. minmoiavail_thisyear )
787	minmoiavail_thisyear = moiavail_daily
788	ENDWHERE
789
790	!
791	! 16 annual maximum weekly GPP
792	!
793
794	WHERE ( gpp_week .GT. maxgppweek_thisyear )
795	maxgppweek_thisyear = gpp_week
796	ENDWHERE
797
798	!
799	! 17 annual GDD0
800	!
801
802	WHERE ( t2m_daily .GT. ZeroCelsius )
803	gdd0_thisyear = gdd0_thisyear + dt * ( t2m_daily - ZeroCelsius )
804	ENDWHERE
805
806	!
807	! 18 annual precipitation
808	!
809
810	precip_thisyear = precip_thisyear + dt * precip_daily
811
812	!
813	! 19 annual maximum leaf mass
814	! If STOMATE is not activated, this corresponds to the maximum possible
815	! LAI of the PFT
816	!
817
818	IF ( control%ok_stomate ) THEN
819
820	DO j = 2,nvm
821	WHERE ( biomass(:,j,ileaf) .GT. lm_thisyearmax(:,j) )
822	lm_thisyearmax(:,j) = biomass(:,j,ileaf)
823	ENDWHERE
824	ENDDO
825
826	ELSE
827
828	DO j = 2,nvm
829	lm_thisyearmax(:,j) = lai_max(j) / sla(j)
830	ENDDO
831
832	ENDIF
833
834	!
835	! 20 annual maximum fpc for each PFT
836	! "veget" is defined as fraction of total ground. Therefore, maxfpc_thisyear has
837	! the same unit.
838	!
839
840	WHERE ( veget(:,:) .GT. maxfpc_thisyear(:,:) )
841	maxfpc_thisyear(:,:) = veget(:,:)
842	ENDWHERE
843
844	!
845	! 21 Every year, replace last year's maximum and minimum moisture availability,
846	! annual GDD0, annual precipitation, annual max weekly GPP, and maximum leaf mass
847
848	IF ( EndOfYear ) THEN
849
850	!
851	! 21.1 replace old values
852	!
853	!NVMODIF
854	maxmoiavail_lastyear(:,:) = (maxmoiavail_lastyear(:,:)*(tau_climatology-1)+ maxmoiavail_thisyear(:,:))/tau_climatology
855	minmoiavail_lastyear(:,:) = (minmoiavail_lastyear(:,:)*(tau_climatology-1)+ minmoiavail_thisyear(:,:))/tau_climatology
856	maxgppweek_lastyear(:,:) =( maxgppweek_lastyear(:,:)*(tau_climatology-1)+ maxgppweek_thisyear(:,:))/tau_climatology
857	! maxmoiavail_lastyear(:,:) = maxmoiavail_thisyear(:,:)
858	! minmoiavail_lastyear(:,:) = minmoiavail_thisyear(:,:)
859	! maxgppweek_lastyear(:,:) = maxgppweek_thisyear(:,:)
860
861	gdd0_lastyear(:) = gdd0_thisyear(:)
862
863	precip_lastyear(:) = precip_thisyear(:)
864
865	lm_lastyearmax(:,:) = lm_thisyearmax(:,:)
866
867	maxfpc_lastyear(:,:) = maxfpc_thisyear(:,:)
868
869	!
870	! 21.2 reset new values
871	!
872
873	maxmoiavail_thisyear(:,:) = zero
874	minmoiavail_thisyear(:,:) = large_value
875
876	maxgppweek_thisyear(:,:) = zero
877
878	gdd0_thisyear(:) = zero
879
880	precip_thisyear(:) = zero
881
882	lm_thisyearmax(:,:) = zero
883
884	maxfpc_thisyear(:,:) = zero
885
886	!
887	! 21.3 Special treatment for maxfpc.
888	!
889
890	!
891	! 21.3.1 Only take into account natural PFTs
892	!
893
894	DO j = 2,nvm
895	IF ( .NOT. natural(j) ) THEN
896	maxfpc_lastyear(:,j) = zero
897	ENDIF
898	ENDDO
899
900	! 21.3.2 In Stomate, veget is defined as a fraction of ground, not as a fraction
901	! of total ground. maxfpc_lastyear will be compared to veget in lpj_light.
902	! Therefore, we have to transform maxfpc_lastyear.
903
904
905	! 21.3.3 The sum of the maxfpc_lastyear for natural PFT must not exceed fpc_crit (=.95).
906	! However, it can slightly exceed this value as not all PFTs reach their maximum
907	! fpc at the same time. Therefore, if sum(maxfpc_lastyear) for the natural PFTs
908	! exceeds fpc_crit, we scale the values of maxfpc_lastyear so that the sum is
909	! fpc_crit.
910
911	! calculate the sum of maxfpc_lastyear
912	sumfpc_nat(:) = zero
913	DO j = 2,nvm
914	sumfpc_nat(:) = sumfpc_nat(:) + maxfpc_lastyear(:,j)
915	ENDDO
916
917	! scale so that the new sum is fpc_crit
918	DO j = 2,nvm
919	WHERE ( sumfpc_nat(:) .GT. fpc_crit )
920	maxfpc_lastyear(:,j) = maxfpc_lastyear(:,j) * (fpc_crit/sumfpc_nat(:))
921	ENDWHERE
922	ENDDO
923
924	ENDIF ! EndOfYear
925
926	!
927	! 22 diagnose herbivore activity, determined through as probability for a leaf to be
928	! eaten in a day
929	! Follows McNaughton et al., Nat 341, 142-144, 1989.
930	!
931
932	!
933	! 22.1 first calculate mean long-term leaf NPP in grid box, mean residence
934	! time (years) of green tissue (i.e. tissue that will be eaten by
935	! herbivores) (crudely approximated: 6 months for seasonal and 2 years
936	! for evergreen) and mean length of growing season (6 months for
937	! seasonal and 1 year for evergreen).
938	!
939
940	!!$ 01/03/2011 NVuichard, NViovy and DS : correcting herbivory activity : now weighttot
941	!!$ nlflong_nat, green_age are pft-dependants
942
943	!!$ nlflong_nat(:) = zero
944	!!$ weighttot(:) = zero
945	!!$ green_age(:) = zero
946	!!$ !
947	!!$ DO j = 2,nvm
948	!!$ !
949	!!$ IF ( natural(j) ) THEN
950	!!$ !
951	!!$ weighttot(:) = weighttot(:) + lm_lastyearmax(:,j)
952	!!$ nlflong_nat(:) = nlflong_nat(:) + npp_longterm(:,j) * leaf_frac_hvc
953	!!$ !
954	!!$ IF ( pheno_model(j) .EQ. 'none' ) THEN
955	!!$ green_age(:) = green_age(:) + green_age_ever * lm_lastyearmax(:,j)
956	!!$ ELSE
957	!!$ green_age(:) = green_age(:) + green_age_dec * lm_lastyearmax(:,j)
958	!!$ ENDIF
959	!!$ !
960	!!$ ENDIF
961	!!$ !
962	!!$ ENDDO
963	!!$ !
964	!!$ WHERE ( weighttot(:) .GT. zero )
965	!!$ green_age(:) = green_age(:) / weighttot(:)
966	!!$ ELSEWHERE
967	!!$ green_age(:) = 1.
968	!!$ ENDWHERE
969	!!$
970	!!$ !
971	!!$ ! 22.2 McNaughton et al. give herbivore consumption as a function of annual leaf NPP.
972	!!$ ! The annual leaf NPP can give us an idea about the edible biomass:
973	!!$ !
974	!!$
975	!!$ DO j = 2,nvm
976	!!$ !
977	!!$ IF ( natural(j) ) THEN
978	!!$ !
979	!!$ WHERE ( nlflong_nat(:) .GT. zero )
980	!!$ consumption(:) = hvc1 * nlflong_nat(:) ** hvc2
981	!!$ herbivores(:,j) = one_year * green_age(:) * nlflong_nat(:) / consumption(:)
982	!!$ ELSEWHERE
983	!!$ herbivores(:,j) = 100000.
984	!!$ ENDWHERE
985	!!$ !
986	!!$ ELSE
987	!!$ !
988	!!$ herbivores(:,j) = 100000.
989	!!$ !
990	!!$ ENDIF
991	!!$ !
992	!!$ ENDDO
993	!!$ herbivores(:,ibare_sechiba) = zero
994
995	nlflong_nat(:,:) = zero
996	weighttot(:,:) = zero
997	green_age(:,:) = zero
998	!
999	DO j = 2,nvm
1000	!
1001	IF ( natural(j) ) THEN
1002	!
1003	weighttot(:,j) = lm_lastyearmax(:,j)
1004	nlflong_nat(:,j) = npp_longterm(:,j) * leaf_frac_hvc
1005	!
1006	IF ( pheno_model(j) .EQ. 'none' ) THEN
1007	green_age(:,j) = green_age_ever * lm_lastyearmax(:,j)
1008	ELSE
1009	green_age(:,j) = green_age_dec * lm_lastyearmax(:,j)
1010	ENDIF
1011	!
1012	ENDIF
1013	!
1014	ENDDO
1015	!
1016	WHERE ( weighttot(:,:) .GT. zero )
1017	green_age(:,:) = green_age(:,:) / weighttot(:,:)
1018	ELSEWHERE
1019	green_age(:,:) = un
1020	ENDWHERE
1021
1022	!
1023	! 22.2 McNaughton et al. give herbivore consumption as a function of annual leaf NPP.
1024	! The annual leaf NPP can give us an idea about the edible biomass:
1025	!
1026
1027	DO j = 2,nvm
1028	!
1029	IF ( natural(j) ) THEN
1030	!
1031	WHERE ( nlflong_nat(:,j) .GT. zero )
1032	consumption(:) = hvc1 * nlflong_nat(:,j) ** hvc2
1033	herbivores(:,j) = one_year * green_age(:,j) * nlflong_nat(:,j) / consumption(:)
1034	ELSEWHERE
1035	herbivores(:,j) = 100000.
1036	ENDWHERE
1037	!
1038	ELSE
1039	!
1040	herbivores(:,j) = 100000.
1041	!
1042	ENDIF
1043	!
1044	ENDDO
1045	herbivores(:,ibare_sechiba) = zero
1046
1047	IF (bavard.GE.4) WRITE(numout,*) 'Leaving season'
1048
1049	END SUBROUTINE season
1050
1051	END MODULE stomate_season

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

Download in other formats: