source: branches/publications/ORCHIDEE_GLUC_r6545/src_stomate/stomate_vmax.f90 @ 6737

Last change on this file since 6737 was 5082, checked in by albert.jornet, 6 years ago

Fix: make sure to initialize variables before those are used (valgrind)
Fix: N_limfert must be passed as assume-shaped array. Otherwise the subroutine assumes a wrong size when ok_LAIdev is not enabled(1,1 vs npts,nvm). So data is written were it's not supposed to due to missmatch.
  • Property svn:keywords set to HeadURL Date Author Revision
File size: 12.8 KB
Line 
1! =================================================================================================================================
2! MODULE        : stomate_vmax
3!
4! CONTACT       : orchidee-help _at_ listes.ipsl.fr
5!
6! LICENCE       : IPSL (2006). This software is governed by the CeCILL licence see ORCHIDEE/ORCHIDEE_CeCILL.LIC
7!
8!>\BRIEF        calculates the leaf efficiency.
9!!     
10!!\n DESCRIPTION: None
11!!
12!! RECENT CHANGE(S): None
13!!
14!! SVN          :
15!! $HeadURL$
16!! $Date$
17!! $Revision$
18!! \n
19!_ =================================================================================================================================
20
21MODULE stomate_vmax
22
23  ! modules used:
24
25  USE ioipsl_para
26  USE stomate_data
27  USE constantes
28  USE pft_parameters
29
30  IMPLICIT NONE
31
32  ! private & public routines
33
34  PRIVATE
35  PUBLIC vmax, vmax_clear
36
37  ! first call
38  LOGICAL, SAVE                                              :: firstcall_vmax = .TRUE.
39!$OMP THREADPRIVATE(firstcall_vmax)
40!gmjc
41  LOGICAL, SAVE                                           :: ok_Nlim = .FALSE.
42!end gmjc
43CONTAINS
44
45!! ================================================================================================================================
46!! SUBROUTINE   : vmax_clear
47!!
48!>\BRIEF          Flag setting
49!!
50!!\n DESCRIPTION: This subroutine sets flags ::firstcall_vmax, to .TRUE., and therefore activates   
51!!                section 1.1 of the ::vmax subroutine which writes messages to the output. \n
52!!                This subroutine is called at the end of the subroutine ::stomate_clear, in the
53!!                module ::stomate.
54!!
55!! RECENT CHANGE(S):None
56!!
57!! MAIN OUTPUT VARIABLE(S): ::firstcall_vmax
58!!
59!! REFERENCE(S)  : None
60!!
61!! FLOWCHART     : None
62!! \n             
63!_ =================================================================================================================================
64
65  SUBROUTINE vmax_clear
66    firstcall_vmax=.TRUE.
67  END SUBROUTINE vmax_clear
68
69
70
71!! ================================================================================================================================
72!! SUBROUTINE    : vmax
73!!
74!>\BRIEF         This subroutine computes vcmax photosynthesis parameters
75!! given optimal vcmax parameter values and a leaf age-related efficiency.
76!!
77!! DESCRIPTION (functional, design, flags):
78!! Leaf age classes are introduced to take into account the fact that photosynthetic activity depends on leaf age
79!! (Ishida et al., 1999). There are \f$nleafages\f$ classes (constant defined in stomate_constants.f90).
80!! This subroutine first calculates the new age of each leaf age-class based on fraction of leaf
81!! that goes from one to another class.                                             
82!! Then calculation of the new fraction of leaf in each class is performed.     
83!! Last, leaf efficiency is calculated for each PFT and for each leaf age class.
84!! vcmax is defined as vcmax25 and vjmax_opt weighted by a mean leaf
85!! efficiency. vcmax25 is PFT-dependent constants defined in constants_mtc.f90.
86!!
87!! This routine is called once at the beginning by stomate_var_init and then at each stomate time step by stomateLpj.
88!!
89!! RECENT CHANGE(S): None
90!!
91!! MAIN OUTPUT VARIABLE(S): vcmax
92!!
93!! REFERENCE(S) :
94!! - Ishida, A., A. Uemura, N. Koike, Y. Matsumoto, and A. Lai Hoe (1999),
95!! Interactive effects of leaf age and self-shading on leaf structure, photosynthetic
96!! capacity and chlorophyll fluorescence in the rain forest tree,
97!! dryobalanops aromatica, Tree Physiol., 19, 741-747
98!!
99!! FLOWCHART    : None
100!!
101!! REVISION(S)  : None
102!! \n
103!_ ================================================================================================================================
104
105  SUBROUTINE vmax (npts, dt, &
106       leaf_age, leaf_frac, &
107       vcmax, &!)
108!gmjc
109       N_limfert)
110!end gmjc
111    !
112    !! 0. Variable and parameter declaration
113    !
114
115    !
116    !! 0.1 Input variables
117    !
118    INTEGER(i_std), INTENT(in)                                 :: npts                    !! Domain size (unitless)
119    REAL(r_std), INTENT(in)                                    :: dt                      !! time step of stomate (days)
120
121    !
122    !! 0.2 Output variables
123    !
124    REAL(r_std), DIMENSION(npts,nvm), INTENT(out)              :: vcmax                   !! Maximum rate of carboxylation
125                                                                                          !! @tex ($\mu mol m^{-2} s^{-1}$) @endtex
126
127    !
128    !! 0.3 Modified variables
129    !
130    REAL(r_std), DIMENSION(npts,nvm,nleafages), INTENT(inout)  :: leaf_age                !! Leaf age (days)
131    REAL(r_std), DIMENSION(npts,nvm,nleafages), INTENT(inout)  :: leaf_frac               !! fraction of leaves in leaf age
132                                                                                          !! classes
133                                                                                          !! (unitless)
134!gmjc
135    ! N fertilization limitation factor for grassland Vcmax and SLA
136    REAL(r_std), DIMENSION(:,:), INTENT(in)                    :: N_limfert
137!end gmjc
138    !
139    !! 0.4 Local variables
140    !
141    REAL(r_std), DIMENSION(npts)                               :: leaf_efficiency         !! leaf efficiency (vcmax/vcmax25)
142                                                                                          !! (unitless)
143    REAL(r_std), DIMENSION(npts,nvm,nleafages)                 :: d_leaf_frac             !! turnover between age classes
144                                                                                          !! (unitless)
145    REAL(r_std), DIMENSION(npts,nleafages)                     :: leaf_age_new            !! new leaf age (days)
146    REAL(r_std), DIMENSION(npts)                               :: sumfrac                 !! sum of leaf age fractions,
147                                                                                          !! for normalization
148                                                                                          !! (unitless)
149    REAL(r_std), DIMENSION(npts)                               :: rel_age                 !! relative leaf age (age/critical age)
150                                                                                          !! (unitless)
151    INTEGER(i_std)                                             :: j,m                     !! indices (unitless)
152
153!_ ================================================================================================================================
154
155    IF (printlev>=3) WRITE(numout,*) 'Entering vmax'
156
157    !
158    !! 1 Initialization
159    !
160
161    !
162    !! 1.1 first call: info about flags and parameters.
163    !
164
165    IF ( firstcall_vmax ) THEN
166!gmjc
167      ok_Nlim=.false.
168      CALL GETIN ('GRM_N_LIMITATION',ok_Nlim)
169      IF (printlev >= 2) THEN
170          WRITE(numout,*) 'GRM_N_LIMITATION',ok_Nlim
171!end gmjc
172          WRITE(numout,*) '   > offset (minimum vcmax/vmax_opt):' , vmax_offset
173          WRITE(numout,*) '   > relative leaf age at which vmax reaches vcmax_opt:', leafage_firstmax 
174          WRITE(numout,*) '   > relative leaf age at which vmax falls below vcmax_opt:', leafage_lastmax
175          WRITE(numout,*) '   > relative leaf age at which vmax reaches its minimum:', leafage_old
176      END IF
177      firstcall_vmax = .FALSE.
178
179    ENDIF
180
181    !
182    !! 1.2 initialize output
183    !
184
185    vcmax(:,:) = zero
186
187    !
188    !! 2 leaf age: general increase and turnover between age classes.
189    !
190
191    !
192    !! 2.1 increase leaf age
193    !
194!
195!! The age of the leaves in each leaf-age-class increases by 1 time step.
196    DO m = 1, nleafages ! Loop over # leaf age classes
197       DO j = 2,nvm     ! Loop over # PFTs
198          WHERE ( leaf_frac(:,j,m) .GT. min_stomate )
199
200             leaf_age(:,j,m) = leaf_age(:,j,m) + dt
201             
202          ENDWHERE
203       ENDDO    ! Loop over # PFTs
204
205    ENDDO       ! Loop over # leaf age classes
206
207    !
208    !! 2.2 turnover between leaf age classes
209    !     d_leaf_frac(:,:,m) = what leaves m-1 and goes into m
210    !
211
212    DO j = 2,nvm        ! Loop over # PFTs
213
214       !! 2.2.1 fluxes
215
216       !! nothing goes into first age class
217       d_leaf_frac(:,j,1) = zero
218
219       !! for others age classes (what goes from m-1 to m)
220       DO m = 2, nleafages 
221!! leaf_timecst is defined in stomate_constants.f90 as the quotient of the critical leaf age per the number of age classes.
222!! The critical leaf age is a PFT-dependent constant defined in stomate_constants.f90, that represents the leaf life span.
223!! This time constant (leaf_timecst) determines the turnover between the nleafages different leaf age classes
224!! (see section [118] in Krinner et al. (2005)).
225          d_leaf_frac(:,j,m) = leaf_frac(:,j,m-1) * dt/leaf_timecst(j)
226
227       ENDDO
228
229       !! 2.2.2 new leaf age in class
230       !!       new age = ( old age * (old fraction - fractional loss) + fractional increase * age of the source class ) / new fraction
231       !!       The leaf age of the youngest class (m=1) is updated into stomate_alloc         
232       leaf_age_new(:,:) = zero
233
234       DO m = 2, nleafages-1       ! Loop over age classes
235        !! For all age classes except first and last
236          WHERE ( d_leaf_frac(:,j,m) .GT. min_stomate )
237
238             leaf_age_new(:,m) = ( ( (leaf_frac(:,j,m)- d_leaf_frac(:,j,m+1)) * leaf_age(:,j,m) )  + &
239                  ( d_leaf_frac(:,j,m) * leaf_age(:,j,m-1) ) ) / &
240                  ( leaf_frac(:,j,m) + d_leaf_frac(:,j,m)- d_leaf_frac(:,j,m+1) )
241
242          ENDWHERE
243
244       ENDDO       ! Loop over age classes
245
246        !! For last age class, there is no leaf fraction leaving the class.
247
248       WHERE ( d_leaf_frac(:,j,nleafages) .GT. min_stomate )
249
250          leaf_age_new(:,nleafages) = ( ( leaf_frac(:,j,nleafages) * leaf_age(:,j,nleafages) )  + &
251               ( d_leaf_frac(:,j,nleafages) * leaf_age(:,j,nleafages-1) ) ) / &
252               ( leaf_frac(:,j,nleafages) + d_leaf_frac(:,j,nleafages) )
253
254       ENDWHERE
255
256       DO m = 2, nleafages       ! Loop over age classes
257
258          WHERE ( d_leaf_frac(:,j,m) .GT. min_stomate )
259
260             leaf_age(:,j,m) = leaf_age_new(:,m)
261
262          ENDWHERE
263
264       ENDDO       ! Loop over age classes
265
266       !! 2.2.3 calculate new fraction
267
268       DO m = 2, nleafages       ! Loop over age classes
269
270          ! where the change comes from
271          leaf_frac(:,j,m-1) = leaf_frac(:,j,m-1) - d_leaf_frac(:,j,m)
272
273          ! where it goes to
274          leaf_frac(:,j,m) = leaf_frac(:,j,m) + d_leaf_frac(:,j,m)
275
276       ENDDO       ! Loop over age classes
277
278       !! 2.2.4 renormalize fractions in order to prevent accumulation
279       !       of numerical errors
280
281       ! correct small negative values
282
283       DO m = 1, nleafages
284          leaf_frac(:,j,m) = MAX( zero, leaf_frac(:,j,m) )
285       ENDDO
286
287       ! total of fractions, should be very close to one where there is leaf mass
288
289       sumfrac(:) = zero
290
291       DO m = 1, nleafages       ! Loop over age classes
292
293          sumfrac(:) = sumfrac(:) + leaf_frac(:,j,m)
294
295       ENDDO       ! Loop over age classes
296
297       ! normalize
298
299       DO m = 1, nleafages       ! Loop over age classes
300
301          WHERE ( sumfrac(:) .GT. min_stomate )
302
303             leaf_frac(:,j,m) = leaf_frac(:,j,m) / sumfrac(:) 
304
305          ELSEWHERE
306
307             leaf_frac(:,j,m) = zero
308
309          ENDWHERE
310
311       ENDDO       ! Loop over age classes
312
313    ENDDO         ! Loop over PFTs
314
315    !
316    !! 3 calculate vmax as a function of the age
317    !
318
319    DO j = 2,nvm
320
321       vcmax(:,j) = zero
322
323       ! sum up over the different age classes
324       IF (ok_dgvm .AND. pheno_type(j)==1 .AND. leaf_tab(j)==2) THEN
325          ! pheno_typ=evergreen and leaf_tab=needleleaf
326          vcmax(:,j) = Vcmax25(j)
327
328       ELSE 
329          ! for deciduous tree
330          DO m = 1, nleafages       ! Loop over age classes
331
332             !
333             !! 3.1 efficiency in each of the age classes
334             !!     it varies from vmax_offset to 1
335             !!     linearly increases from vmax_offset to 1 for 0 < rel_age < leafage_firstmax
336             !!     is 1 when leafage_firstmax < rel_age < leafage_lastmax
337             !!     linearly decreases from 1 to vmax_offset for leafage_lastmax < rel_age < leafage_firstmax
338             !!     vmax_offset for rel_age >= leafage_old
339             !!     (Ishida et al., 1999)
340             rel_age(:) = leaf_age(:,j,m) / leafagecrit(j)
341
342             leaf_efficiency(:) = MAX( vmax_offset, MIN( un, &
343                  vmax_offset + (un - vmax_offset) * rel_age(:) / leafage_firstmax, &
344                  un - (un - vmax_offset) * ( rel_age(:) - leafage_lastmax ) / &
345                  ( leafage_old - leafage_lastmax ) ) )
346
347             !
348             !! 3.2 add to mean vmax
349             !             
350             IF (ok_Nlim .OR. ok_LAIdev(j)) THEN
351                vcmax(:,j) = vcmax(:,j) + Vcmax25(j) * N_limfert(:,j)*leaf_efficiency(:) * leaf_frac(:,j,m)
352             ELSE
353                vcmax(:,j) = vcmax(:,j) + vcmax25(j) * leaf_efficiency(:) * leaf_frac(:,j,m)
354             ENDIF
355          ENDDO ! loop over age classes
356       ENDIF
357
358    ENDDO       ! loop over PFTs
359
360    IF (printlev>=4) WRITE(numout,*) 'Leaving vmax'
361
362  END SUBROUTINE vmax
363
364END MODULE stomate_vmax
Note: See TracBrowser for help on using the repository browser.