source: branches/ORCHIDEE_2_2/ORCHIDEE/src_parameters/constantes_var.f90 @ 8377

Last change on this file since 8377 was 8377, checked in by jan.polcher, 6 months ago

The modifications performed on the trunk for the coupling to WRF and the interpolation by XIOS on curviliean grids has been backported.

  • Property svn:keywords set to Date Revision
File size: 59.6 KB
Line 
1! =================================================================================================================================
2! MODULE       : constantes_var
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        constantes_var module contains most constantes like pi, Earth radius, etc...
10!!              and all externalized parameters except pft-dependent constants.
11!!
12!!\n DESCRIPTION: This module contains most constantes and the externalized parameters of ORCHIDEE which
13!!                are not pft-dependent.\n
14!!                In this module, you can set the flag diag_qsat in order to detect the pixel where the
15!!                temperature is out of range (look qsatcalc and dev_qsatcalc in qsat_moisture.f90).\n
16!!                The Earth radius is approximated by the Equatorial radius.The Earth's equatorial radius a,
17!!                or semi-major axis, is the distance from its center to the equator and equals 6,378.1370 km.
18!!                The equatorial radius is often used to compare Earth with other planets.\n
19!!                The meridional mean is well approximated by the semicubic mean of the two axe yielding
20!!                6367.4491 km or less accurately by the quadratic mean of the two axes about 6,367.454 km
21!!                or even just the mean of the two axes about 6,367.445 km.\n
22!!                This module is already USE in module constantes. Therefor no need to USE it seperatly except
23!!                if the subroutines in module constantes are not needed.\n
24!!               
25!! RECENT CHANGE(S):
26!!
27!! REFERENCE(S) :
28!! - Louis, Jean-Francois (1979), A parametric model of vertical eddy fluxes in the atmosphere.
29!! Boundary Layer Meteorology, 187-202.\n
30!!
31!! SVN          :
32!! $HeadURL: $
33!! $Date$
34!! $Revision$
35!! \n
36!_ ================================================================================================================================
37
38MODULE constantes_var
39
40  USE defprec
41
42  IMPLICIT NONE
43!-
44
45                         !-----------------------!
46                         !  ORCHIDEE CONSTANTS   !
47                         !-----------------------!
48
49  !
50  ! FLAGS
51  !
52  LOGICAL :: river_routing      !! activate river routing
53!$OMP THREADPRIVATE(river_routing)
54  LOGICAL, SAVE :: ok_nudge_mc  !! Activate nudging of soil moisture
55!$OMP THREADPRIVATE(ok_nudge_mc)
56  LOGICAL, SAVE :: ok_nudge_snow!! Activate nudging of snow variables
57!$OMP THREADPRIVATE(ok_nudge_snow)
58  LOGICAL, SAVE :: nudge_interpol_with_xios  !! Activate reading and interpolation with XIOS for nudging fields
59!$OMP THREADPRIVATE(nudge_interpol_with_xios)
60  LOGICAL :: do_floodplains     !! activate flood plains
61!$OMP THREADPRIVATE(do_floodplains)
62  LOGICAL :: do_irrigation      !! activate computation of irrigation flux
63!$OMP THREADPRIVATE(do_irrigation)
64  LOGICAL :: ok_sechiba         !! activate physic of the model
65!$OMP THREADPRIVATE(ok_sechiba)
66  LOGICAL :: ok_stomate         !! activate carbon cycle
67!$OMP THREADPRIVATE(ok_stomate)
68  LOGICAL :: ok_dgvm            !! activate dynamic vegetation
69!$OMP THREADPRIVATE(ok_dgvm)
70  LOGICAL :: do_wood_harvest    !! activate wood harvest
71!$OMP THREADPRIVATE(do_wood_harvest)
72  LOGICAL :: ok_pheno           !! activate the calculation of lai using stomate rather than a prescription
73!$OMP THREADPRIVATE(ok_pheno)
74  LOGICAL :: ok_bvoc            !! activate biogenic volatile organic coumpounds
75!$OMP THREADPRIVATE(ok_bvoc)
76  LOGICAL :: ok_leafage         !! activate leafage
77!$OMP THREADPRIVATE(ok_leafage)
78  LOGICAL :: ok_radcanopy       !! use canopy radiative transfer model
79!$OMP THREADPRIVATE(ok_radcanopy)
80  LOGICAL :: ok_multilayer      !! use canopy radiative transfer model with multi-layers
81!$OMP THREADPRIVATE(ok_multilayer)
82  LOGICAL :: ok_pulse_NOx       !! calculate NOx emissions with pulse
83!$OMP THREADPRIVATE(ok_pulse_NOx)
84  LOGICAL :: ok_bbgfertil_NOx   !! calculate NOx emissions with bbg fertilizing effect
85!$OMP THREADPRIVATE(ok_bbgfertil_NOx)
86  LOGICAL :: ok_cropsfertil_NOx !! calculate NOx emissions with fertilizers use
87!$OMP THREADPRIVATE(ok_cropsfertil_NOx)
88
89  LOGICAL :: ok_co2bvoc_poss    !! CO2 inhibition on isoprene activated following Possell et al. (2005) model
90!$OMP THREADPRIVATE(ok_co2bvoc_poss)
91  LOGICAL :: ok_co2bvoc_wilk    !! CO2 inhibition on isoprene activated following Wilkinson et al. (2006) model
92!$OMP THREADPRIVATE(ok_co2bvoc_wilk)
93 
94  LOGICAL, SAVE :: OFF_LINE_MODE = .FALSE.  !! ORCHIDEE detects if it is coupled with a GCM or
95                                            !! just use with one driver in OFF-LINE. (true/false)
96!$OMP THREADPRIVATE(OFF_LINE_MODE) 
97  LOGICAL, SAVE :: impose_param = .TRUE.    !! Flag impos_param : read all the parameters in the run.def file
98!$OMP THREADPRIVATE(impose_param)
99  CHARACTER(LEN=80), SAVE     :: restname_in       = 'NONE'                 !! Input Restart files name for Sechiba component 
100!$OMP THREADPRIVATE(restname_in)
101  CHARACTER(LEN=80), SAVE     :: restname_out      = 'sechiba_rest_out.nc'  !! Output Restart files name for Sechiba component
102!$OMP THREADPRIVATE(restname_out)
103  CHARACTER(LEN=80), SAVE     :: stom_restname_in  = 'NONE'                 !! Input Restart files name for Stomate component
104!$OMP THREADPRIVATE(stom_restname_in)
105  CHARACTER(LEN=80), SAVE     :: stom_restname_out = 'stomate_rest_out.nc'  !! Output Restart files name for Stomate component
106!$OMP THREADPRIVATE(stom_restname_out)
107  INTEGER, SAVE :: printlev=2       !! Standard level for text output [0, 1, 2, 3]
108!$OMP THREADPRIVATE(printlev)
109
110  !
111  ! TIME
112  !
113  INTEGER(i_std), PARAMETER  :: spring_days_max = 40  !! Maximum number of days during which we watch for possible spring frost damage
114  !
115  ! SPECIAL VALUES
116  !
117  INTEGER(i_std), PARAMETER :: undef_int = 999999999     !! undef integer for integer arrays (unitless)
118  !-
119  REAL(r_std), SAVE :: val_exp = 999999.                 !! Specific value if no restart value  (unitless)
120!$OMP THREADPRIVATE(val_exp)
121  REAL(r_std), PARAMETER :: undef = -9999.               !! Special value for stomate (unitless)
122 
123  REAL(r_std), PARAMETER :: min_sechiba = 1.E-8_r_std    !! Epsilon to detect a near zero floating point (unitless)
124  REAL(r_std), PARAMETER :: undef_sechiba = 1.E+20_r_std !! The undef value used in SECHIBA (unitless)
125 
126  REAL(r_std), PARAMETER :: min_stomate = 1.E-8_r_std    !! Epsilon to detect a near zero floating point (unitless)
127  REAL(r_std), PARAMETER :: large_value = 1.E33_r_std    !! some large value (for stomate) (unitless)
128
129  REAL(r_std), SAVE :: alpha_nudge_mc                    !! Nudging constant for soil moisture
130!$OMP THREADPRIVATE(alpha_nudge_mc)
131  REAL(r_std), SAVE :: alpha_nudge_snow                  !! Nudging constant for snow variables
132!$OMP THREADPRIVATE(alpha_nudge_snow)
133
134  !
135  !  DIMENSIONING AND INDICES PARAMETERS 
136  !
137  INTEGER(i_std), PARAMETER :: ibare_sechiba = 1 !! Index for bare soil in Sechiba (unitless)
138  INTEGER(i_std), PARAMETER :: ivis = 1          !! index for albedo in visible range (unitless)
139  INTEGER(i_std), PARAMETER :: inir = 2          !! index for albeod i near-infrared range (unitless)
140  INTEGER(i_std), PARAMETER :: nnobio = 1        !! Number of other surface types: land ice (lakes,cities, ...) (unitless)
141  INTEGER(i_std), PARAMETER :: iice = 1          !! Index for land ice (see nnobio) (unitless)
142  !-
143  !! Soil
144  INTEGER(i_std), PARAMETER :: classnb = 9       !! Levels of soil colour classification (unitless)
145  !-
146  INTEGER(i_std), PARAMETER :: nleafages = 4     !! leaf age discretisation ( 1 = no discretisation )(unitless)
147  !-
148  !! litter fractions: indices (unitless)
149  INTEGER(i_std), PARAMETER :: ileaf = 1         !! Index for leaf compartment (unitless)
150  INTEGER(i_std), PARAMETER :: isapabove = 2     !! Index for sapwood above compartment (unitless)
151  INTEGER(i_std), PARAMETER :: isapbelow = 3     !! Index for sapwood below compartment (unitless)
152  INTEGER(i_std), PARAMETER :: iheartabove = 4   !! Index for heartwood above compartment (unitless)
153  INTEGER(i_std), PARAMETER :: iheartbelow = 5   !! Index for heartwood below compartment (unitless)
154  INTEGER(i_std), PARAMETER :: iroot = 6         !! Index for roots compartment (unitless)
155  INTEGER(i_std), PARAMETER :: ifruit = 7        !! Index for fruits compartment (unitless)
156  INTEGER(i_std), PARAMETER :: icarbres = 8      !! Index for reserve compartment (unitless)
157  INTEGER(i_std), PARAMETER :: nparts = 8        !! Number of biomass compartments (unitless)
158  !-
159  !! indices for assimilation parameters
160  INTEGER(i_std), PARAMETER :: ivcmax = 1        !! Index for vcmax (assimilation parameters) (unitless)
161  INTEGER(i_std), PARAMETER :: npco2 = 1         !! Number of assimilation parameters (unitless)
162  !-
163  !! trees and litter: indices for the parts of heart-
164  !! and sapwood above and below the ground
165  INTEGER(i_std), PARAMETER :: iabove = 1       !! Index for above part (unitless)
166  INTEGER(i_std), PARAMETER :: ibelow = 2       !! Index for below part (unitless)
167  INTEGER(i_std), PARAMETER :: nlevs = 2        !! Number of levels for trees and litter (unitless)
168  !-
169  !! litter: indices for metabolic and structural part
170  INTEGER(i_std), PARAMETER :: imetabolic = 1   !! Index for metabolic litter (unitless)
171  INTEGER(i_std), PARAMETER :: istructural = 2  !! Index for structural litter (unitless)
172  INTEGER(i_std), PARAMETER :: nlitt = 2        !! Number of levels for litter compartments (unitless)
173  !-
174  !! carbon pools: indices
175  INTEGER(i_std), PARAMETER :: iactive = 1      !! Index for active carbon pool (unitless)
176  INTEGER(i_std), PARAMETER :: islow = 2        !! Index for slow carbon pool (unitless)
177  INTEGER(i_std), PARAMETER :: ipassive = 3     !! Index for passive carbon pool (unitless)
178  INTEGER(i_std), PARAMETER :: ncarb = 3        !! Number of soil carbon pools (unitless)
179  !-
180  !! For isotopes and nitrogen
181  INTEGER(i_std), PARAMETER :: nelements = 1    !! Number of isotopes considered
182  INTEGER(i_std), PARAMETER :: icarbon = 1      !! Index for carbon
183  !
184  !! Indices used for analytical spin-up
185  INTEGER(i_std), PARAMETER :: nbpools = 7              !! Total number of carbon pools (unitless)
186  INTEGER(i_std), PARAMETER :: istructural_above = 1    !! Index for structural litter above (unitless)
187  INTEGER(i_std), PARAMETER :: istructural_below = 2    !! Index for structural litter below (unitless)
188  INTEGER(i_std), PARAMETER :: imetabolic_above = 3     !! Index for metabolic litter above (unitless)
189  INTEGER(i_std), PARAMETER :: imetabolic_below = 4     !! Index for metabolic litter below (unitless)
190  INTEGER(i_std), PARAMETER :: iactive_pool = 5         !! Index for active carbon pool (unitless)
191  INTEGER(i_std), PARAMETER :: islow_pool   = 6         !! Index for slow carbon pool (unitless)
192  INTEGER(i_std), PARAMETER :: ipassive_pool = 7        !! Index for passive carbon pool (unitless)
193  !
194  !! Indicies used for output variables on Landuse tiles defined according to LUMIP project
195  !! Note that ORCHIDEE do not represent pasture and urban land. Therefor the variables will have
196  !! val_exp as missing value for these tiles.
197  INTEGER(i_std), PARAMETER :: nlut=4                   !! Total number of landuse tiles according to LUMIP
198  INTEGER(i_std), PARAMETER :: id_psl=1                 !! Index for primary and secondary land
199  INTEGER(i_std), PARAMETER :: id_crp=2                 !! Index for crop land
200  INTEGER(i_std), PARAMETER :: id_pst=3                 !! Index for pasture land
201  INTEGER(i_std), PARAMETER :: id_urb=4                 !! Index for urban land
202
203
204  !
205  ! NUMERICAL AND PHYSICS CONSTANTS
206  !
207  !
208
209  !-
210  ! 1. Mathematical and numerical constants
211  !-
212  REAL(r_std), PARAMETER :: pi = 3.141592653589793238   !! pi souce : http://mathworld.wolfram.com/Pi.html (unitless)
213  REAL(r_std), PARAMETER :: euler = 2.71828182845904523 !! e source : http://mathworld.wolfram.com/e.html (unitless)
214  REAL(r_std), PARAMETER :: zero = 0._r_std             !! Numerical constant set to 0 (unitless)
215  REAL(r_std), PARAMETER :: undemi = 0.5_r_std          !! Numerical constant set to 1/2 (unitless)
216  REAL(r_std), PARAMETER :: un = 1._r_std               !! Numerical constant set to 1 (unitless)
217  REAL(r_std), PARAMETER :: moins_un = -1._r_std        !! Numerical constant set to -1 (unitless)
218  REAL(r_std), PARAMETER :: deux = 2._r_std             !! Numerical constant set to 2 (unitless)
219  REAL(r_std), PARAMETER :: trois = 3._r_std            !! Numerical constant set to 3 (unitless)
220  REAL(r_std), PARAMETER :: quatre = 4._r_std           !! Numerical constant set to 4 (unitless)
221  REAL(r_std), PARAMETER :: cinq = 5._r_std             !![DISPENSABLE] Numerical constant set to 5 (unitless)
222  REAL(r_std), PARAMETER :: six = 6._r_std              !![DISPENSABLE] Numerical constant set to 6 (unitless)
223  REAL(r_std), PARAMETER :: huit = 8._r_std             !! Numerical constant set to 8 (unitless)
224  REAL(r_std), PARAMETER :: mille = 1000._r_std         !! Numerical constant set to 1000 (unitless)
225
226  !-
227  ! 2 . Physics
228  !-
229#ifndef OASIS
230  REAL(r_std), PARAMETER :: R_Earth = 6378000.              !! Radius of the Earth : Earth radius ~= Equatorial radius (m)
231#else
232  REAL(r_std), PARAMETER :: R_Earth = 6370000.              !! When coupled through OASIS it is with WRF. We have to use another Earth radius 
233#endif
234  REAL(r_std), PARAMETER :: mincos  = 0.0001                !! Minimum cosine value used for interpolation (unitless)
235  REAL(r_std), PARAMETER :: pb_std = 1013.                  !! standard pressure (hPa)
236  REAL(r_std), PARAMETER :: ZeroCelsius = 273.15            !! 0 degre Celsius in degre Kelvin (K)
237  REAL(r_std), PARAMETER :: tp_00 = 273.15                  !! 0 degre Celsius in degre Kelvin (K)
238  REAL(r_std), PARAMETER :: chalsu0 = 2.8345E06             !! Latent heat of sublimation (J.kg^{-1})
239  REAL(r_std), PARAMETER :: chalev0 = 2.5008E06             !! Latent heat of evaporation (J.kg^{-1})
240  REAL(r_std), PARAMETER :: chalfu0 = chalsu0-chalev0       !! Latent heat of fusion (J.kg^{-1})
241  REAL(r_std), PARAMETER :: c_stefan = 5.6697E-8            !! Stefan-Boltzman constant (W.m^{-2}.K^{-4})
242  REAL(r_std), PARAMETER :: cp_air = 1004.675               !! Specific heat of dry air (J.kg^{-1}.K^{-1})
243  REAL(r_std), PARAMETER :: cte_molr = 287.05               !! Specific constant of dry air (kg.mol^{-1})
244  REAL(r_std), PARAMETER :: kappa = cte_molr/cp_air         !! Kappa : ratio between specific constant and specific heat
245                                                            !! of dry air (unitless)
246  REAL(r_std), PARAMETER :: msmlr_air = 28.964E-03          !! Molecular weight of dry air (kg.mol^{-1})
247  REAL(r_std), PARAMETER :: msmlr_h2o = 18.02E-03           !! Molecular weight of water vapor (kg.mol^{-1})
248  REAL(r_std), PARAMETER :: cp_h2o = &                      !! Specific heat of water vapor (J.kg^{-1}.K^{-1})
249       & cp_air*(quatre*msmlr_air)/( 3.5_r_std*msmlr_h2o) 
250  REAL(r_std), PARAMETER :: cte_molr_h2o = cte_molr/quatre  !! Specific constant of water vapor (J.kg^{-1}.K^{-1})
251  REAL(r_std), PARAMETER :: retv = msmlr_air/msmlr_h2o-un   !! Ratio between molecular weight of dry air and water
252                                                            !! vapor minus 1(unitless) 
253  REAL(r_std), PARAMETER :: rvtmp2 = cp_h2o/cp_air-un       !! Ratio between specific heat of water vapor and dry air
254                                                            !! minus 1 (unitless)
255  REAL(r_std), PARAMETER :: cepdu2 = (0.1_r_std)**2         !! Squared wind shear (m^2.s^{-2})
256  REAL(r_std), PARAMETER :: ct_karman = 0.41_r_std          !! Van Karmann Constant (unitless)
257  REAL(r_std), PARAMETER :: cte_grav = 9.80665_r_std        !! Acceleration of the gravity (m.s^{-2})
258  REAL(r_std), PARAMETER :: pa_par_hpa = 100._r_std         !! Transform pascal into hectopascal (unitless)
259  REAL(r_std), PARAMETER :: RR = 8.314                      !! Ideal gas constant (J.mol^{-1}.K^{-1})
260  REAL(r_std), PARAMETER :: Sct = 1370.                     !! Solar constant (W.m^{-2})
261
262
263  INTEGER(i_std), SAVE :: testpft = 6
264!$OMP THREADPRIVATE(testpft)
265  !-
266  ! 3. Climatic constants
267  !-
268  !! Constantes of the Louis scheme
269  REAL(r_std), SAVE :: cb = 5._r_std              !! Constant of the Louis scheme (unitless);
270                                                  !! reference to Louis (1979)
271!$OMP THREADPRIVATE(cb)
272  REAL(r_std), SAVE :: cc = 5._r_std              !! Constant of the Louis scheme (unitless);
273                                                  !! reference to Louis (1979)
274!$OMP THREADPRIVATE(cc)
275  REAL(r_std), SAVE :: cd = 5._r_std              !! Constant of the Louis scheme (unitless);
276                                                  !! reference to Louis (1979)
277!$OMP THREADPRIVATE(cd)
278  REAL(r_std), SAVE :: rayt_cste = 125.           !! Constant in the computation of surface resistance (W.m^{-2})
279!$OMP THREADPRIVATE(rayt_cste)
280  REAL(r_std), SAVE :: defc_plus = 23.E-3         !! Constant in the computation of surface resistance (K.W^{-1})
281!$OMP THREADPRIVATE(defc_plus)
282  REAL(r_std), SAVE :: defc_mult = 1.5            !! Constant in the computation of surface resistance (K.W^{-1})
283!$OMP THREADPRIVATE(defc_mult)
284
285  !-
286  ! 4. Soil thermodynamics constants
287  !-
288  ! Look at constantes_soil.f90
289
290
291  !
292  ! OPTIONAL PARTS OF THE MODEL
293  !
294  LOGICAL,PARAMETER :: diag_qsat = .TRUE.         !! One of the most frequent problems is a temperature out of range
295                                                  !! we provide here a way to catch that in the calling procedure.
296                                                  !! (from Jan Polcher)(true/false)
297  LOGICAL, SAVE     :: almaoutput =.FALSE.        !! Selects the type of output for the model.(true/false)
298                                                  !! Value is read from run.def in intersurf_history
299!$OMP THREADPRIVATE(almaoutput)
300
301  !
302  ! DIVERSE
303  !
304  CHARACTER(LEN=100), SAVE :: stomate_forcing_name='NONE'  !! NV080800 Name of STOMATE forcing file (unitless)
305                                                           ! Compatibility with Nicolas Viovy driver.
306!$OMP THREADPRIVATE(stomate_forcing_name)
307  CHARACTER(LEN=100), SAVE :: stomate_Cforcing_name='NONE' !! NV080800 Name of soil forcing file (unitless)
308                                                           ! Compatibility with Nicolas Viovy driver.
309!$OMP THREADPRIVATE(stomate_Cforcing_name)
310  INTEGER(i_std), SAVE :: forcing_id                 !! Index of the forcing file (unitless)
311!$OMP THREADPRIVATE(forcing_id)
312  LOGICAL, SAVE :: allow_forcing_write=.TRUE.        !! Allow writing of stomate_forcing file.
313                                                     !! This variable will be set to false for teststomate.
314!$OMP THREADPRIVATE(allow_forcing_write)
315
316
317
318                         !------------------------!
319                         !  SECHIBA PARAMETERS    !
320                         !------------------------!
321 
322
323  !
324  ! GLOBAL PARAMETERS   
325  !
326  REAL(r_std), SAVE :: min_wind = 0.1      !! The minimum wind (m.s^{-1})
327!$OMP THREADPRIVATE(min_wind)
328  REAL(r_std), PARAMETER :: min_qc = 1.e-4 !! The minimum value for qc (qc=drag*wind) used in coupled(enerbil) and forced mode (enerbil and diffuco)
329  REAL(r_std), SAVE :: snowcri = 1.5       !! Sets the amount above which only sublimation occures (kg.m^{-2})
330!$OMP THREADPRIVATE(snowcri)
331
332
333  !
334  ! FLAGS ACTIVATING SUB-MODELS
335  !
336  LOGICAL, SAVE :: treat_expansion = .FALSE.   !! Do we treat PFT expansion across a grid point after introduction? (true/false)
337!$OMP THREADPRIVATE(treat_expansion)
338  LOGICAL, SAVE :: ok_herbivores = .FALSE.     !! flag to activate herbivores (true/false)
339!$OMP THREADPRIVATE(ok_herbivores)
340  LOGICAL, SAVE :: harvest_agri = .TRUE.       !! flag to harvest aboveground biomass from agricultural PFTs)(true/false)
341!$OMP THREADPRIVATE(harvest_agri)
342  LOGICAL, SAVE :: lpj_gap_const_mort          !! constant moratlity (true/false). Default value depend on OK_DGVM.
343!$OMP THREADPRIVATE(lpj_gap_const_mort)
344  LOGICAL, SAVE :: disable_fire = .TRUE.       !! flag that disable fire (true/false)
345!$OMP THREADPRIVATE(disable_fire)
346  LOGICAL, SAVE :: spinup_analytic = .FALSE.   !! Flag to activate analytical resolution for spinup (true/false)
347!$OMP THREADPRIVATE(spinup_analytic)
348
349  !
350  ! CONFIGURATION VEGETATION
351  !
352  LOGICAL, SAVE :: agriculture = .TRUE.    !! allow agricultural PFTs (true/false)
353!$OMP THREADPRIVATE(agriculture)
354  LOGICAL, SAVE :: impveg = .FALSE.        !! Impose vegetation ? (true/false)
355!$OMP THREADPRIVATE(impveg)
356  LOGICAL, SAVE :: impsoilt = .FALSE.      !! Impose soil ? (true/false)
357!$OMP THREADPRIVATE(impsoilt)
358  LOGICAL, SAVE :: impslope = .FALSE.      !! Impose reinf_slope ? (true/false)
359!$OMP THREADPRIVATE(impslope)
360  LOGICAL, SAVE :: do_now_stomate_lcchange = .FALSE.  !! Time to call lcchange in stomate_lpj
361!$OMP THREADPRIVATE(do_now_stomate_lcchange)
362  LOGICAL, SAVE :: do_now_stomate_woodharvest = .FALSE.  !! Time to call woodharvest in stomate_lpj
363!$OMP THREADPRIVATE(do_now_stomate_woodharvest)
364  LOGICAL, SAVE :: done_stomate_lcchange = .FALSE.    !! If true, call lcchange in stomate_lpj has just been done.
365!$OMP THREADPRIVATE(done_stomate_lcchange)
366  LOGICAL, SAVE :: read_lai = .FALSE.      !! Flag to read a map of LAI if STOMATE is not activated (true/false)
367!$OMP THREADPRIVATE(read_lai)
368  LOGICAL, SAVE :: veget_reinit = .TRUE.   !! To change LAND USE file in a run. (true/false)
369!$OMP THREADPRIVATE(veget_reinit)
370  LOGICAL, SAVE :: vegetmap_reset = .FALSE.!! Reset the vegetation map and reset carbon related variables
371!$OMP THREADPRIVATE(vegetmap_reset)
372  INTEGER(i_std) , SAVE :: veget_update    !! Update frequency in years for landuse (nb of years)
373!$OMP THREADPRIVATE(veget_update)
374  !
375  ! PARAMETERS USED BY BOTH HYDROLOGY MODELS
376  !
377  REAL(r_std), SAVE :: max_snow_age = 50._r_std !! Maximum period of snow aging (days)
378!$OMP THREADPRIVATE(max_snow_age)
379  REAL(r_std), SAVE :: snow_trans = 0.2_r_std   !! Transformation time constant for snow (m), reduced from the value 0.3 (04/07/2016)
380!$OMP THREADPRIVATE(snow_trans)
381  REAL(r_std), SAVE :: sneige                   !! Lower limit of snow amount (kg.m^{-2})
382!$OMP THREADPRIVATE(sneige)
383  REAL(r_std), SAVE :: maxmass_snow = 3000.     !! The maximum mass of snow (kg.m^{-2})
384!$OMP THREADPRIVATE(maxmass_snow)
385
386  !! Heat capacity
387  REAL(r_std), PARAMETER :: rho_water = 1000.           !! Density of water (kg/m3)
388  REAL(r_std), PARAMETER :: rho_ice = 920.              !! Density of ice (kg/m3)
389  REAL(r_std), PARAMETER :: rho_soil = 2700.            !! Density of soil particles (kg/m3), value from Peters-Lidard et al. 1998
390
391  !! Thermal conductivities
392  REAL(r_std), PARAMETER :: cond_water = 0.6            !! Thermal conductivity of liquid water (W/m/K)
393  REAL(r_std), PARAMETER :: cond_ice = 2.2              !! Thermal conductivity of ice (W/m/K)
394  REAL(r_std), PARAMETER :: cond_solid = 2.32           !! Thermal conductivity of mineral soil particles (W/m/K)
395
396  !! Time constant of long-term soil humidity (s)
397  REAL(r_std), PARAMETER :: lhf = 0.3336*1.E6           !! Latent heat of fusion (J/kg)
398
399  INTEGER(i_std), PARAMETER :: nsnow=3                  !! Number of levels in the snow for explicit snow scheme   
400  REAL(r_std), PARAMETER    :: XMD    = 28.9644E-3 
401  REAL(r_std), PARAMETER    :: XBOLTZ      = 1.380658E-23 
402  REAL(r_std), PARAMETER    :: XAVOGADRO   = 6.0221367E+23 
403  REAL(r_std), PARAMETER    :: XRD    = XAVOGADRO * XBOLTZ / XMD 
404  REAL(r_std), PARAMETER    :: XCPD   = 7.* XRD /2. 
405  REAL(r_std), PARAMETER    :: phigeoth = 0.057 ! 0. DKtest
406  REAL(r_std), PARAMETER    :: thick_min_snow = .01 
407
408  !! The maximum snow density and water holding characterisicts
409  REAL(r_std), SAVE         :: xrhosmax = 750.  ! (kg m-3)
410!$OMP THREADPRIVATE(xrhosmax)
411  REAL(r_std), SAVE         :: xwsnowholdmax1   = 0.03  ! (-)
412!$OMP THREADPRIVATE(xwsnowholdmax1)
413  REAL(r_std), SAVE         :: xwsnowholdmax2   = 0.10  ! (-)
414!$OMP THREADPRIVATE(xwsnowholdmax2)
415  REAL(r_std), SAVE         :: xsnowrhohold     = 200.0 ! (kg/m3)
416!$OMP THREADPRIVATE(xsnowrhohold)
417  REAL(r_std), SAVE         :: xrhosmin = 50. 
418!$OMP THREADPRIVATE(xrhosmin)
419  REAL(r_std), PARAMETER    :: xci = 2.106e+3 
420  REAL(r_std), PARAMETER    :: xrv = 6.0221367e+23 * 1.380658e-23 /18.0153e-3 
421
422  !! ISBA-ES Critical snow depth at which snow grid thicknesses constant
423  REAL(r_std), PARAMETER    :: xsnowcritd = 0.03  ! (m)
424
425  !! The threshold of snow depth used for preventing numerical problem in thermal calculations
426  REAL(r_std), PARAMETER    :: snowcritd_thermal = 0.01  ! (m) 
427 
428  !! ISBA-ES CROCUS (Pahaut 1976): snowfall density coefficients:
429  REAL(r_std), PARAMETER       :: snowfall_a_sn = 109.0  !! (kg/m3)
430  REAL(r_std), PARAMETER       :: snowfall_b_sn =   6.0  !! (kg/m3/K)
431  REAL(r_std), PARAMETER       :: snowfall_c_sn =  26.0  !! [kg/(m7/2 s1/2)]
432
433  REAL(r_std), PARAMETER       :: dgrain_new_max=  2.0e-4!! (m) : Maximum grain size of new snowfall
434 
435  !! Used in explicitsnow to prevent numerical problems as snow becomes vanishingly thin.
436  REAL(r_std), PARAMETER                :: psnowdzmin = .0001   ! m
437  REAL(r_std), PARAMETER                :: xsnowdmin = .000001  ! m
438
439  REAL(r_std), PARAMETER                :: ph2o = 1000.         !! Water density [kg/m3]
440 
441  ! ISBA-ES Thermal conductivity coefficients from Anderson (1976):
442  ! see Boone, Meteo-France/CNRM Note de Centre No. 70 (2002)
443  REAL(r_std), SAVE                     :: ZSNOWTHRMCOND1 = 0.02    ! [W/m/K]
444!$OMP THREADPRIVATE(ZSNOWTHRMCOND1)
445  REAL(r_std), SAVE                     :: ZSNOWTHRMCOND2 = 2.5E-6  ! [W m5/(kg2 K)]
446!$OMP THREADPRIVATE(ZSNOWTHRMCOND2)
447 
448  ! ISBA-ES Thermal conductivity: Implicit vapor diffn effects
449  ! (sig only for new snow OR high altitudes)
450  ! from Sun et al. (1999): based on data from Jordan (1991)
451  ! see Boone, Meteo-France/CNRM Note de Centre No. 70 (2002)
452  !
453  REAL(r_std), SAVE                       :: ZSNOWTHRMCOND_AVAP  = -0.06023 ! (W/m/K)
454!$OMP THREADPRIVATE(ZSNOWTHRMCOND_AVAP)
455  REAL(r_std), SAVE                       :: ZSNOWTHRMCOND_BVAP  = -2.5425  ! (W/m)
456!$OMP THREADPRIVATE(ZSNOWTHRMCOND_BVAP)
457  REAL(r_std), SAVE                       :: ZSNOWTHRMCOND_CVAP  = -289.99  ! (K)
458!$OMP THREADPRIVATE(ZSNOWTHRMCOND_CVAP)
459 
460  REAL(r_std),SAVE :: xansmax = 0.85      !! Maxmimum snow albedo
461!$OMP THREADPRIVATE(xansmax)
462  REAL(r_std),SAVE :: xansmin = 0.50      !! Miniumum snow albedo
463!$OMP THREADPRIVATE(xansmin)
464  REAL(r_std),SAVE :: xans_todry = 0.008  !! Albedo decay rate for dry snow
465!$OMP THREADPRIVATE(xans_todry)
466  REAL(r_std),SAVE :: xans_t = 0.240      !! Albedo decay rate for wet snow
467!$OMP THREADPRIVATE(xans_t)
468
469  ! ISBA-ES Thermal conductivity coefficients from Anderson (1976):
470  ! see Boone, Meteo-France/CNRM Note de Centre No. 70 (2002)
471  REAL(r_std), PARAMETER                  :: XP00 = 1.E5
472
473  ! ISBA-ES Thermal conductivity: Implicit vapor diffn effects
474  ! (sig only for new snow OR high altitudes)
475  ! from Sun et al. (1999): based on data from Jordan (1991)
476  ! see Boone, Meteo-France/CNRM Note de Centre No. 70 (2002)
477  !
478  REAL(r_std), SAVE          :: ZSNOWCMPCT_RHOD  = 150.0        !! (kg/m3)
479!$OMP THREADPRIVATE(ZSNOWCMPCT_RHOD)
480  REAL(r_std), SAVE          :: ZSNOWCMPCT_ACM   = 2.8e-6       !! (1/s
481!$OMP THREADPRIVATE(ZSNOWCMPCT_ACM)
482  REAL(r_std), SAVE          :: ZSNOWCMPCT_BCM   = 0.04         !! (1/K)
483!$OMP THREADPRIVATE(ZSNOWCMPCT_BCM)
484  REAL(r_std), SAVE          :: ZSNOWCMPCT_CCM   = 460.         !! (m3/kg)
485!$OMP THREADPRIVATE(ZSNOWCMPCT_CCM)
486  REAL(r_std), SAVE          :: ZSNOWCMPCT_V0    = 3.7e7        !! (Pa/s)
487!$OMP THREADPRIVATE(ZSNOWCMPCT_V0)
488  REAL(r_std), SAVE          :: ZSNOWCMPCT_VT    = 0.081        !! (1/K)
489!$OMP THREADPRIVATE(ZSNOWCMPCT_VT)
490  REAL(r_std), SAVE          :: ZSNOWCMPCT_VR    = 0.018        !! (m3/kg)
491!$OMP THREADPRIVATE(ZSNOWCMPCT_VR)
492
493  !
494  ! BVOC : Biogenic activity  for each age class
495  !
496  REAL(r_std), SAVE, DIMENSION(nleafages) :: iso_activity = (/0.5, 1.5, 1.5, 0.5/)     !! Biogenic activity for each
497                                                                                       !! age class : isoprene (unitless)
498!$OMP THREADPRIVATE(iso_activity)
499  REAL(r_std), SAVE, DIMENSION(nleafages) :: methanol_activity = (/1., 1., 0.5, 0.5/)  !! Biogenic activity for each
500                                                                                       !! age class : methanol (unnitless)
501!$OMP THREADPRIVATE(methanol_activity)
502
503
504 !
505 ! Parameters for irrigation scheme
506 !
507  REAL(r_std), SAVE :: irrig_dosmax = 3.              !! The maximum irrigation water injected per hour (kg.m^{-2}/hour)
508!$OMP THREADPRIVATE(irrig_dosmax)
509  REAL(r_std), SAVE :: cum_dh_thr = 0.64                !! Cumulated nroot threshoold to define root zone, and calculate water deficit for irrigation (-)
510!$OMP THREADPRIVATE(cum_dh_thr)
511  LOGICAL, SAVE :: irrigated_soiltile = .FALSE.       !! Do we introduce a new soil tile for irrigated croplands? (true/false)
512!$OMP THREADPRIVATE(irrigated_soiltile)
513  LOGICAL, SAVE :: old_irrig_scheme = .FALSE.         !! Do we run with the old irrigation scheme? (true/false)  , add to compatiblity
514!$OMP THREADPRIVATE(old_irrig_scheme)
515  INTEGER, SAVE :: irrig_st = 3                       !! Which is the soil tile with irrigation flux
516!$OMP THREADPRIVATE(irrig_st)
517  REAL(r_std), SAVE, DIMENSION(3) :: avail_reserve = (/0.9,0.9,0.9/)     !! Available water from routing reservoirs, to withdraw for irrigation
518                                                      !! IMPORTANT: As the routing model uses 3 reservoirs, dimension is set to 3
519                                                      !! IMPORTANT: Order of available water must be in this order: streamflow, fast, and slow reservoir
520!$OMP THREADPRIVATE(avail_reserve)
521  REAL(r_std), SAVE :: beta_irrig = 0.9               !! Threshold multiplier of Target SM to calculate root deficit(unitless)
522!$OMP THREADPRIVATE(beta_irrig)
523  REAL(r_std), SAVE :: lai_irrig_min = 0.1            !! Minimum LAI to trigger irrigation (kg.m^{-2}/hour)
524!$OMP THREADPRIVATE(lai_irrig_min)
525  LOGICAL, SAVE :: irrig_map_dynamic_flag = .FALSE.   !! Do we use a dynamic irrig map?
526!$OMP THREADPRIVATE(irrig_map_dynamic_flag)
527  LOGICAL, SAVE :: select_source_irrig = .FALSE.      !! Do we use the new priorization scheme, based on maps of equipped area with surface water?
528!$OMP THREADPRIVATE(select_source_irrig)
529  LOGICAL, SAVE :: Reinfiltr_IrrigField = .FALSE.     !! Do we reinfiltrate all runoff from crop soil tile?O
530!$OMP THREADPRIVATE(Reinfiltr_IrrigField)
531  REAL, SAVE :: reinf_slope_cropParam = 0.8           !! Externalized for irrigated cropland, when Reinfiltr_IrrigField=.TRUE.
532                                                      !! Max value of reinf_slope in irrig_st 
533!$OMP THREADPRIVATE(reinf_slope_cropParam)
534  REAL, SAVE :: a_stream_adduction = 0.05             !! Externalized for available volume to adduction
535!$OMP THREADPRIVATE(a_stream_adduction)
536
537
538
539  !
540  ! condveg.f90
541  !
542
543  ! 1. Scalar
544
545  ! 1.1 Flags used inside the module
546
547  LOGICAL, SAVE :: alb_bare_model = .FALSE. !! Switch for choosing values of bare soil
548                                            !! albedo (see header of subroutine)
549                                            !! (true/false)
550!$OMP THREADPRIVATE(alb_bare_model)
551  LOGICAL, SAVE :: alb_bg_modis = .TRUE.    !! Switch for choosing values of bare soil
552                                            !! albedo read from file
553                                            !! (true/false)
554!$OMP THREADPRIVATE(alb_bg_modis)
555  LOGICAL, SAVE :: impaze = .FALSE.         !! Switch for choosing surface parameters
556                                            !! (see header of subroutine). 
557                                            !! (true/false)
558!$OMP THREADPRIVATE(impaze)
559  LOGICAL, SAVE :: rough_dyn = .TRUE.       !! Chooses between two methods to calculate the
560                                            !! the roughness height : static or dynamic (varying with LAI)
561                                            !! (true/false)
562!$OMP THREADPRIVATE(rough_dyn)
563  LOGICAL, SAVE :: use_ratio_z0m_z0h = .FALSE. !! To impose a constant ratio as done in ROUGH_DYN=F
564!$OMP THREADPRIVATE(use_ratio_z0m_z0h)
565
566  LOGICAL, SAVE :: new_watstress = .FALSE.
567!$OMP THREADPRIVATE(new_watstress)
568
569  REAL(r_std), SAVE :: alpha_watstress = 1.
570!$OMP THREADPRIVATE(alpha_watstress)
571
572  ! 1.2 Others
573
574
575  REAL(r_std), SAVE :: height_displacement = 0.66        !! Factor to calculate the zero-plane displacement
576                                                         !! height from vegetation height (m)
577!$OMP THREADPRIVATE(height_displacement)
578  REAL(r_std), SAVE :: z0_bare = 0.01                    !! bare soil roughness length (m)
579!$OMP THREADPRIVATE(z0_bare)
580  REAL(r_std), SAVE :: z0_ice = 0.001                    !! ice roughness length (m)
581!$OMP THREADPRIVATE(z0_ice)
582  REAL(r_std), SAVE :: tcst_snowa = 10.0                 !! Time constant of the albedo decay of snow (days), increased from the value 5.0 (04/07/2016)
583!$OMP THREADPRIVATE(tcst_snowa)
584  REAL(r_std), SAVE :: snowcri_alb = 10.                 !! Critical value for computation of snow albedo (cm)
585!$OMP THREADPRIVATE(snowcri_alb)
586  REAL(r_std), SAVE :: fixed_snow_albedo = undef_sechiba !! To choose a fixed snow albedo value (unitless)
587!$OMP THREADPRIVATE(fixed_snow_albedo)
588  REAL(r_std), SAVE :: z0_scal = 0.15                    !! Surface roughness height imposed (m)
589!$OMP THREADPRIVATE(z0_scal)
590  REAL(r_std), SAVE :: roughheight_scal = zero           !! Effective roughness Height depending on zero-plane
591                                                         !! displacement height (m) (imposed)
592!$OMP THREADPRIVATE(roughheight_scal)
593  REAL(r_std), SAVE :: emis_scal = 1.0                   !! Surface emissivity imposed (unitless)
594!$OMP THREADPRIVATE(emis_scal)
595
596  REAL(r_std), SAVE :: c1 = 0.32                         !! Constant used in the formulation of the ratio of
597!$OMP THREADPRIVATE(c1)                                  !! friction velocity to the wind speed at the canopy top
598                                                         !! see Ershadi et al. (2015) for more info
599  REAL(r_std), SAVE :: c2 = 0.264                        !! Constant used in the formulation of the ratio of
600!$OMP THREADPRIVATE(c2)                                  !! friction velocity to the wind speed at the canopy top
601                                                         !! see Ershadi et al. (2015) for more info
602  REAL(r_std), SAVE :: c3 = 15.1                         !! Constant used in the formulation of the ratio of
603!$OMP THREADPRIVATE(c3)                                  !! friction velocity to the wind speed at the canopy top
604                                                         !! see Ershadi et al. (2015) for more info
605  REAL(r_std), SAVE :: Cdrag_foliage = 0.2               !! Drag coefficient of the foliage
606!$OMP THREADPRIVATE(Cdrag_foliage)                       !! See Ershadi et al. (2015) and Su et. al (2001) for more info
607  REAL(r_std), SAVE :: Ct = 0.01                         !! Heat transfer coefficient of the leaf
608!$OMP THREADPRIVATE(Ct)                                  !! See Ershadi et al. (2015) and Su et. al (2001) for more info
609  REAL(r_std), SAVE :: Prandtl = 0.71                    !! Prandtl number used in the calculation of Ct_star
610!$OMP THREADPRIVATE(Prandtl)                             !! See Su et. al (2001) for more info
611
612
613
614  ! 2. Arrays
615
616  REAL(r_std), SAVE, DIMENSION(2) :: alb_deadleaf = (/ .12, .35/)    !! albedo of dead leaves, VIS+NIR (unitless)
617!$OMP THREADPRIVATE(alb_deadleaf)
618  REAL(r_std), SAVE, DIMENSION(2) :: alb_ice = (/ .60, .20/)         !! albedo of ice, VIS+NIR (unitless)
619!$OMP THREADPRIVATE(alb_ice)
620  REAL(r_std), SAVE, DIMENSION(2) :: albedo_scal = (/ 0.25, 0.25 /)  !! Albedo values for visible and near-infrared
621                                                                     !! used imposed (unitless)
622!$OMP THREADPRIVATE(albedo_scal)
623  REAL(r_std) , SAVE, DIMENSION(classnb) :: vis_dry = (/0.24,&
624       &0.22, 0.20, 0.18, 0.16, 0.14, 0.12, 0.10, 0.27/)  !! Soil albedo values to soil colour classification:
625                                                          !! dry soil albedo values in visible range
626!$OMP THREADPRIVATE(vis_dry)
627  REAL(r_std), SAVE, DIMENSION(classnb) :: nir_dry = (/0.48,&
628       &0.44, 0.40, 0.36, 0.32, 0.28, 0.24, 0.20, 0.55/)  !! Soil albedo values to soil colour classification:
629                                                          !! dry soil albedo values in near-infrared range
630!$OMP THREADPRIVATE(nir_dry)
631  REAL(r_std), SAVE, DIMENSION(classnb) :: vis_wet = (/0.12,&
632       &0.11, 0.10, 0.09, 0.08, 0.07, 0.06, 0.05, 0.15/)  !! Soil albedo values to soil colour classification:
633                                                          !! wet soil albedo values in visible range
634!$OMP THREADPRIVATE(vis_wet)
635  REAL(r_std), SAVE, DIMENSION(classnb) :: nir_wet = (/0.24,&
636       &0.22, 0.20, 0.18, 0.16, 0.14, 0.12, 0.10, 0.31/)  !! Soil albedo values to soil colour classification:
637                                                          !! wet soil albedo values in near-infrared range
638!$OMP THREADPRIVATE(nir_wet)
639  REAL(r_std), SAVE, DIMENSION(classnb) :: albsoil_vis = (/ &
640       &0.18, 0.16, 0.16, 0.15, 0.12, 0.105, 0.09, 0.075, 0.25/)   !! Soil albedo values to soil colour classification:
641                                                                   !! Averaged of wet and dry soil albedo values
642                                                                   !! in visible and near-infrared range
643!$OMP THREADPRIVATE(albsoil_vis)
644  REAL(r_std), SAVE, DIMENSION(classnb) :: albsoil_nir = (/ &
645       &0.36, 0.34, 0.34, 0.33, 0.30, 0.25, 0.20, 0.15, 0.45/)  !! Soil albedo values to soil colour classification:
646                                                                !! Averaged of wet and dry soil albedo values
647                                                                !! in visible and near-infrared range
648!$OMP THREADPRIVATE(albsoil_nir)
649
650  !
651  ! diffuco.f90
652  !
653
654  ! 0. Constants
655
656  REAL(r_std), PARAMETER :: Tetens_1 = 0.622         !! Ratio between molecular weight of water vapor and molecular weight 
657                                                     !! of dry air (unitless)
658  REAL(r_std), PARAMETER :: Tetens_2 = 0.378         !!
659  REAL(r_std), PARAMETER :: ratio_H2O_to_CO2 = 1.6   !! Ratio of water vapor diffusivity to the CO2 diffusivity (unitless)
660  REAL(r_std), PARAMETER :: mol_to_m_1 = 0.0244      !!
661  REAL(r_std), PARAMETER :: RG_to_PAR = 0.5          !!
662  REAL(r_std), PARAMETER :: W_to_mol = 4.6          !! W_to_mmol * RG_to_PAR = 2.3
663
664  ! 1. Scalar
665
666  INTEGER(i_std), SAVE :: nlai = 20             !! Number of LAI levels (unitless)
667!$OMP THREADPRIVATE(nlai)
668  LOGICAL, SAVE :: ldq_cdrag_from_gcm = .FALSE. !! Set to .TRUE. if you want q_cdrag coming from GCM
669!$OMP THREADPRIVATE(ldq_cdrag_from_gcm)
670  REAL(r_std), SAVE :: laimax = 12.             !! Maximal LAI used for splitting LAI into N layers (m^2.m^{-2})
671!$OMP THREADPRIVATE(laimax)
672  LOGICAL, SAVE :: downregulation_co2 = .TRUE.             !! Set to .TRUE. if you want CO2 downregulation version used for CMIP6 6.1.0-6.1.10
673!$OMP THREADPRIVATE(downregulation_co2)
674  LOGICAL, SAVE :: downregulation_co2_new = .FALSE.        !! Set to .TRUE. if you want CO2 downregulation version revised for CMIP6 6.1.11
675!$OMP THREADPRIVATE(downregulation_co2_new)
676  REAL(r_std), SAVE :: downregulation_co2_baselevel = 380. !! CO2 base level (ppm)
677!$OMP THREADPRIVATE(downregulation_co2_baselevel)
678  REAL(r_std), SAVE :: downregulation_co2_minimum = 280.   !! CO2 value above which downregulation is taken into account
679!$OMP THREADPRIVATE(downregulation_co2_minimum)
680
681  REAL(r_std), SAVE :: gb_ref = 1./25.                     !! Leaf bulk boundary layer resistance (s m-1)
682!$OMP THREADPRIVATE(gb_ref)
683
684  ! 3. Coefficients of equations
685
686  REAL(r_std), SAVE :: lai_level_depth = 0.15  !!
687!$OMP THREADPRIVATE(lai_level_depth)
688!
689  REAL(r_std), SAVE, DIMENSION(6) :: dew_veg_poly_coeff = &            !! coefficients of the 5 degree polynomomial used
690  & (/ 0.887773, 0.205673, 0.110112, 0.014843,  0.000824,  0.000017 /) !! in the equation of coeff_dew_veg
691!$OMP THREADPRIVATE(dew_veg_poly_coeff)
692!
693  REAL(r_std), SAVE               :: Oi=210000.    !! Intercellular oxygen partial pressure (ubar)
694!$OMP THREADPRIVATE(Oi)
695  !
696  ! slowproc.f90
697  !
698
699  ! 1. Scalar
700
701  INTEGER(i_std), SAVE :: veget_year_orig = 0        !!  first year for landuse (number)
702!$OMP THREADPRIVATE(veget_year_orig)
703 
704  REAL(r_std), SAVE :: min_vegfrac = 0.001           !! Minimal fraction of mesh a vegetation type can occupy (0-1, unitless)
705!$OMP THREADPRIVATE(min_vegfrac)
706  REAL(r_std), SAVE :: frac_nobio_fixed_test_1 = 0.0 !! Value for frac_nobio for tests in 0-dim simulations (0-1, unitless)
707!$OMP THREADPRIVATE(frac_nobio_fixed_test_1)
708 
709  REAL(r_std), SAVE :: stempdiag_bid = 280.          !! only needed for an initial LAI if there is no restart file
710!$OMP THREADPRIVATE(stempdiag_bid)
711
712
713                           !-----------------------------!
714                           !  STOMATE AND LPJ PARAMETERS !
715                           !-----------------------------!
716
717
718  !
719  ! lpj_constraints.f90
720  !
721 
722  ! 1. Scalar
723
724  REAL(r_std), SAVE  :: too_long = 5.      !! longest sustainable time without
725                                           !! regeneration (vernalization) (years)
726!$OMP THREADPRIVATE(too_long)
727
728
729  !
730  ! lpj_establish.f90
731  !
732
733  ! 1. Scalar
734
735  REAL(r_std), SAVE :: estab_max_tree = 0.12   !! Maximum tree establishment rate (ind/m2/dt_stomate)
736!$OMP THREADPRIVATE(estab_max_tree)
737  REAL(r_std), SAVE :: estab_max_grass = 0.12  !! Maximum grass establishment rate (ind/m2/dt_stomate)
738!$OMP THREADPRIVATE(estab_max_grass)
739 
740  ! 3. Coefficients of equations
741
742  REAL(r_std), SAVE :: establish_scal_fact = 5.  !!
743!$OMP THREADPRIVATE(establish_scal_fact)
744  REAL(r_std), SAVE :: max_tree_coverage = 0.98  !! (0-1, unitless)
745!$OMP THREADPRIVATE(max_tree_coverage)
746  REAL(r_std), SAVE :: ind_0_estab = 0.2         !! = ind_0 * 10.
747!$OMP THREADPRIVATE(ind_0_estab)
748
749
750  !
751  ! lpj_fire.f90
752  !
753
754  ! 1. Scalar
755
756  REAL(r_std), SAVE :: tau_fire = 30.           !! Time scale for memory of the fire index (days).
757!$OMP THREADPRIVATE(tau_fire)
758  REAL(r_std), SAVE :: litter_crit = 200.       !! Critical litter quantity for fire
759                                                !! below which iginitions extinguish
760                                                !! @tex $(gC m^{-2})$ @endtex
761!$OMP THREADPRIVATE(litter_crit)
762  REAL(r_std), SAVE :: fire_resist_struct = 0.5 !!
763!$OMP THREADPRIVATE(fire_resist_struct)
764  ! 2. Arrays
765
766  REAL(r_std), SAVE, DIMENSION(nparts) :: co2frac = &    !! The fraction of the different biomass
767       & (/ .95, .95, 0., 0.3, 0., 0., .95, .95 /)       !! compartments emitted to the atmosphere
768!$OMP THREADPRIVATE(co2frac)                                                         !! when burned (unitless, 0-1) 
769
770  ! 3. Coefficients of equations
771
772  REAL(r_std), SAVE, DIMENSION(3) :: bcfrac_coeff = (/ .3,  1.3,  88.2 /)         !! (unitless)
773!$OMP THREADPRIVATE(bcfrac_coeff)
774  REAL(r_std), SAVE, DIMENSION(4) :: firefrac_coeff = (/ 0.45, 0.8, 0.6, 0.13 /)  !! (unitless)
775!$OMP THREADPRIVATE(firefrac_coeff)
776
777  !
778  ! lpj_gap.f90
779  !
780
781  ! 1. Scalar
782
783  REAL(r_std), SAVE :: ref_greff = 0.035         !! Asymptotic maximum mortality rate
784                                                 !! @tex $(year^{-1})$ @endtex
785!$OMP THREADPRIVATE(ref_greff)
786
787  !               
788  ! lpj_light.f90
789  !             
790
791  ! 1. Scalar
792 
793  LOGICAL, SAVE :: annual_increase = .TRUE. !! for diagnosis of fpc increase, compare today's fpc to last year's maximum (T) or
794                                            !! to fpc of last time step (F)? (true/false)
795!$OMP THREADPRIVATE(annual_increase)
796  REAL(r_std), SAVE :: min_cover = 0.05     !! For trees, minimum fraction of crown area occupied
797                                            !! (due to its branches etc.) (0-1, unitless)
798                                            !! This means that only a small fraction of its crown area
799                                            !! can be invaded by other trees.
800!$OMP THREADPRIVATE(min_cover)
801  !
802  ! lpj_pftinout.f90
803  !
804
805  ! 1. Scalar
806
807  REAL(r_std), SAVE :: min_avail = 0.01         !! minimum availability
808!$OMP THREADPRIVATE(min_avail)
809  REAL(r_std), SAVE :: ind_0 = 0.02             !! initial density of individuals
810!$OMP THREADPRIVATE(ind_0)
811  ! 3. Coefficients of equations
812 
813  REAL(r_std), SAVE :: RIP_time_min = 1.25      !! test whether the PFT has been eliminated lately (years)
814!$OMP THREADPRIVATE(RIP_time_min)
815  REAL(r_std), SAVE :: npp_longterm_init = 10.  !! Initialisation value for npp_longterm (gC.m^{-2}.year^{-1})
816!$OMP THREADPRIVATE(npp_longterm_init)
817  REAL(r_std), SAVE :: everywhere_init = 0.05   !!
818!$OMP THREADPRIVATE(everywhere_init)
819
820
821  !
822  ! stomate_alloc.f90
823  !
824
825  ! 0. Constants
826
827  REAL(r_std), PARAMETER :: max_possible_lai = 10. !! (m^2.m^{-2})
828  REAL(r_std), PARAMETER :: Nlim_Q10 = 10.         !!
829
830  ! 1. Scalar
831
832  LOGICAL, SAVE :: ok_minres = .TRUE.              !! [DISPENSABLE] Do we try to reach a minimum reservoir even if
833                                                   !! we are severely stressed? (true/false)
834!$OMP THREADPRIVATE(ok_minres)
835  REAL(r_std), SAVE :: reserve_time_tree = 30.     !! Maximum number of days during which
836                                                   !! carbohydrate reserve may be used for
837                                                   !! trees (days)
838!$OMP THREADPRIVATE(reserve_time_tree)
839  REAL(r_std), SAVE :: reserve_time_grass = 20.    !! Maximum number of days during which
840                                                   !! carbohydrate reserve may be used for
841                                                   !! grasses (days)
842!$OMP THREADPRIVATE(reserve_time_grass)
843
844  REAL(r_std), SAVE :: f_fruit = 0.1               !! Default fruit allocation (0-1, unitless)
845!$OMP THREADPRIVATE(f_fruit)
846  REAL(r_std), SAVE :: alloc_sap_above_grass = 1.0 !! fraction of sapwood allocation above ground
847                                                   !! for grass (0-1, unitless)
848!$OMP THREADPRIVATE(alloc_sap_above_grass)
849  REAL(r_std), SAVE :: min_LtoLSR = 0.2            !! Prescribed lower bounds for leaf
850                                                   !! allocation (0-1, unitless)
851!$OMP THREADPRIVATE(min_LtoLSR)
852  REAL(r_std), SAVE :: max_LtoLSR = 0.5            !! Prescribed upper bounds for leaf
853                                                   !! allocation (0-1, unitless)
854!$OMP THREADPRIVATE(max_LtoLSR)
855  REAL(r_std), SAVE :: z_nitrogen = 0.2            !! Curvature of the root profile (m)
856!$OMP THREADPRIVATE(z_nitrogen)
857
858  ! 3. Coefficients of equations
859
860  REAL(r_std), SAVE :: Nlim_tref = 25.             !! (C)
861!$OMP THREADPRIVATE(Nlim_tref)
862
863
864  !
865  ! stomate_data.f90
866  !
867
868  ! 1. Scalar
869
870  ! 1.1 Parameters for the pipe model
871
872  REAL(r_std), SAVE :: pipe_tune1 = 100.0        !! crown area = pipe_tune1. stem diameter**(1.6) (Reinicke's theory) (unitless)
873!$OMP THREADPRIVATE(pipe_tune1)
874  REAL(r_std), SAVE :: pipe_tune2 = 40.0         !! height=pipe_tune2 * diameter**pipe_tune3 (unitless)
875!$OMP THREADPRIVATE(pipe_tune2)
876  REAL(r_std), SAVE :: pipe_tune3 = 0.5          !! height=pipe_tune2 * diameter**pipe_tune3 (unitless)
877!$OMP THREADPRIVATE(pipe_tune3)
878  REAL(r_std), SAVE :: pipe_tune4 = 0.3          !! needed for stem diameter (unitless)
879!$OMP THREADPRIVATE(pipe_tune4)
880  REAL(r_std), SAVE :: pipe_density = 2.e5       !! Density
881!$OMP THREADPRIVATE(pipe_density)
882  REAL(r_std), SAVE :: pipe_k1 = 8.e3            !! one more SAVE
883!$OMP THREADPRIVATE(pipe_k1)
884  REAL(r_std), SAVE :: pipe_tune_exp_coeff = 1.6 !! pipe tune exponential coeff (unitless)
885!$OMP THREADPRIVATE(pipe_tune_exp_coeff)
886
887  ! 1.2 climatic parameters
888
889  REAL(r_std), SAVE :: precip_crit = 100.        !! minimum precip, in (mm/year)
890!$OMP THREADPRIVATE(precip_crit)
891  REAL(r_std), SAVE :: gdd_crit_estab = 150.     !! minimum gdd for establishment of saplings
892!$OMP THREADPRIVATE(gdd_crit_estab)
893  REAL(r_std), SAVE :: fpc_crit = 0.95           !! critical fpc, needed for light competition and establishment (0-1, unitless)
894!$OMP THREADPRIVATE(fpc_crit)
895
896  ! 1.3 sapling characteristics
897
898  REAL(r_std), SAVE :: alpha_grass = 0.5         !! alpha coefficient for grasses (unitless)
899!$OMP THREADPRIVATE(alpha_grass)
900  REAL(r_std), SAVE :: alpha_tree = 1.           !! alpha coefficient for trees (unitless)
901!$OMP THREADPRIVATE(alpha_tree)
902  REAL(r_std), SAVE :: mass_ratio_heart_sap = 3. !! mass ratio (heartwood+sapwood)/sapwood (unitless)
903!$OMP THREADPRIVATE(mass_ratio_heart_sap)
904
905  ! 1.4  time scales for phenology and other processes (in days)
906
907  REAL(r_std), SAVE :: tau_hum_month = 20.        !! (days)       
908!$OMP THREADPRIVATE(tau_hum_month)
909  REAL(r_std), SAVE :: tau_hum_week = 7.          !! (days) 
910!$OMP THREADPRIVATE(tau_hum_week)
911  REAL(r_std), SAVE :: tau_t2m_month = 20.        !! (days)     
912!$OMP THREADPRIVATE(tau_t2m_month)
913  REAL(r_std), SAVE :: tau_t2m_week = 7.          !! (days) 
914!$OMP THREADPRIVATE(tau_t2m_week)
915  REAL(r_std), SAVE :: tau_tsoil_month = 20.      !! (days)     
916!$OMP THREADPRIVATE(tau_tsoil_month)
917  REAL(r_std), SAVE :: tau_soilhum_month = 20.    !! (days)     
918!$OMP THREADPRIVATE(tau_soilhum_month)
919  REAL(r_std), SAVE :: tau_gpp_week = 7.          !! (days) 
920!$OMP THREADPRIVATE(tau_gpp_week)
921  REAL(r_std), SAVE :: tau_gdd = 40.              !! (days) 
922!$OMP THREADPRIVATE(tau_gdd)
923  REAL(r_std), SAVE :: tau_ngd = 50.              !! (days) 
924!$OMP THREADPRIVATE(tau_ngd)
925  REAL(r_std), SAVE :: coeff_tau_longterm = 3.    !! (unitless)
926!$OMP THREADPRIVATE(coeff_tau_longterm)
927  REAL(r_std), SAVE :: tau_longterm_max           !! (days) 
928!$OMP THREADPRIVATE(tau_longterm_max)
929
930  ! 3. Coefficients of equations
931
932  REAL(r_std), SAVE :: bm_sapl_carbres = 5.             !!
933!$OMP THREADPRIVATE(bm_sapl_carbres)
934  REAL(r_std), SAVE :: bm_sapl_sapabove = 0.5           !!
935!$OMP THREADPRIVATE(bm_sapl_sapabove)
936  REAL(r_std), SAVE :: bm_sapl_heartabove = 2.          !!
937!$OMP THREADPRIVATE(bm_sapl_heartabove)
938  REAL(r_std), SAVE :: bm_sapl_heartbelow = 2.          !!
939!$OMP THREADPRIVATE(bm_sapl_heartbelow)
940  REAL(r_std), SAVE :: init_sapl_mass_leaf_nat = 0.1    !!
941!$OMP THREADPRIVATE(init_sapl_mass_leaf_nat)
942  REAL(r_std), SAVE :: init_sapl_mass_leaf_agri = 1.    !!
943!$OMP THREADPRIVATE(init_sapl_mass_leaf_agri)
944  REAL(r_std), SAVE :: init_sapl_mass_carbres = 5.      !!
945!$OMP THREADPRIVATE(init_sapl_mass_carbres)
946  REAL(r_std), SAVE :: init_sapl_mass_root = 0.1        !!
947!$OMP THREADPRIVATE(init_sapl_mass_root)
948  REAL(r_std), SAVE :: init_sapl_mass_fruit = 0.3       !! 
949!$OMP THREADPRIVATE(init_sapl_mass_fruit)
950  REAL(r_std), SAVE :: cn_sapl_init = 0.5               !!
951!$OMP THREADPRIVATE(cn_sapl_init)
952  REAL(r_std), SAVE :: migrate_tree = 10.*1.E3          !!
953!$OMP THREADPRIVATE(migrate_tree)
954  REAL(r_std), SAVE :: migrate_grass = 10.*1.E3         !!
955!$OMP THREADPRIVATE(migrate_grass)
956  REAL(r_std), SAVE :: lai_initmin_tree = 0.3           !!
957!$OMP THREADPRIVATE(lai_initmin_tree)
958  REAL(r_std), SAVE :: lai_initmin_grass = 0.1          !!
959!$OMP THREADPRIVATE(lai_initmin_grass)
960  REAL(r_std), SAVE, DIMENSION(2) :: dia_coeff = (/ 4., 0.5 /)            !!
961!$OMP THREADPRIVATE(dia_coeff)
962  REAL(r_std), SAVE, DIMENSION(2) :: maxdia_coeff =(/ 100., 0.01/)        !!
963!$OMP THREADPRIVATE(maxdia_coeff)
964  REAL(r_std), SAVE, DIMENSION(4) :: bm_sapl_leaf = (/ 4., 4., 0.8, 5./)  !!
965!$OMP THREADPRIVATE(bm_sapl_leaf)
966
967
968
969  !
970  ! stomate_litter.f90
971  !
972
973  ! 0. Constants
974
975  REAL(r_std), PARAMETER :: Q10 = 10.               !!
976
977  ! 1. Scalar
978
979  REAL(r_std), SAVE :: z_decomp = 0.2               !!  Maximum depth for soil decomposer's activity (m)
980!$OMP THREADPRIVATE(z_decomp)
981
982  ! 2. Arrays
983
984  REAL(r_std), SAVE :: frac_soil_struct_aa = 0.55   !! corresponding to frac_soil(istructural,iactive,iabove)
985!$OMP THREADPRIVATE(frac_soil_struct_aa)
986  REAL(r_std), SAVE :: frac_soil_struct_ab = 0.45   !! corresponding to frac_soil(istructural,iactive,ibelow)
987!$OMP THREADPRIVATE(frac_soil_struct_ab)
988  REAL(r_std), SAVE :: frac_soil_struct_sa = 0.7    !! corresponding to frac_soil(istructural,islow,iabove)
989!$OMP THREADPRIVATE(frac_soil_struct_sa)
990  REAL(r_std), SAVE :: frac_soil_struct_sb = 0.7    !! corresponding to frac_soil(istructural,islow,ibelow)
991!$OMP THREADPRIVATE(frac_soil_struct_sb)
992  REAL(r_std), SAVE :: frac_soil_metab_aa = 0.45    !! corresponding to frac_soil(imetabolic,iactive,iabove)
993!$OMP THREADPRIVATE(frac_soil_metab_aa)
994  REAL(r_std), SAVE :: frac_soil_metab_ab = 0.45    !! corresponding to frac_soil(imetabolic,iactive,ibelow)
995!$OMP THREADPRIVATE(frac_soil_metab_ab)
996  REAL(r_std), SAVE, DIMENSION(nparts) :: CN = &    !! C/N ratio of each plant pool (0-100, unitless)
997       & (/ 40., 40., 40., 40., 40., 40., 40., 40. /) 
998!$OMP THREADPRIVATE(CN)
999  REAL(r_std), SAVE, DIMENSION(nparts) :: LC = &    !! Lignin/C ratio of different plant parts (0,22-0,35, unitless)
1000       & (/ 0.22, 0.35, 0.35, 0.35, 0.35, 0.22, 0.22, 0.22 /)
1001!$OMP THREADPRIVATE(LC)
1002
1003  ! 3. Coefficients of equations
1004
1005  REAL(r_std), SAVE :: metabolic_ref_frac = 0.85    !! used by litter and soilcarbon (0-1, unitless)
1006!$OMP THREADPRIVATE(metabolic_ref_frac)
1007  REAL(r_std), SAVE :: metabolic_LN_ratio = 0.018   !! (0-1, unitless)   
1008!$OMP THREADPRIVATE(metabolic_LN_ratio)
1009  REAL(r_std), SAVE :: tau_metabolic = 0.066        !!
1010!$OMP THREADPRIVATE(tau_metabolic)
1011  REAL(r_std), SAVE :: tau_struct = 0.245           !!
1012!$OMP THREADPRIVATE(tau_struct)
1013  REAL(r_std), SAVE :: soil_Q10 = 0.69              !!= ln 2
1014!$OMP THREADPRIVATE(soil_Q10)
1015  REAL(r_std), SAVE :: tsoil_ref = 30.              !!
1016!$OMP THREADPRIVATE(tsoil_ref)
1017  REAL(r_std), SAVE :: litter_struct_coef = 3.      !!
1018!$OMP THREADPRIVATE(litter_struct_coef)
1019  REAL(r_std), SAVE, DIMENSION(3) :: moist_coeff = (/ 1.1,  2.4,  0.29 /) !!
1020!$OMP THREADPRIVATE(moist_coeff)
1021  REAL(r_std), SAVE :: moistcont_min = 0.25  !! minimum soil wetness to limit the heterotrophic respiration
1022!$OMP THREADPRIVATE(moistcont_min)
1023
1024
1025  !
1026  ! stomate_lpj.f90
1027  !
1028
1029  ! 1. Scalar
1030
1031  REAL(r_std), SAVE :: frac_turnover_daily = 0.55  !! (0-1, unitless)
1032!$OMP THREADPRIVATE(frac_turnover_daily)
1033
1034
1035  !
1036  ! stomate_npp.f90
1037  !
1038
1039  ! 1. Scalar
1040
1041  REAL(r_std), SAVE :: tax_max = 0.8 !! Maximum fraction of allocatable biomass used
1042                                     !! for maintenance respiration (0-1, unitless)
1043!$OMP THREADPRIVATE(tax_max)
1044
1045
1046  !
1047  ! stomate_phenology.f90
1048  !
1049
1050  ! 1. Scalar
1051
1052  REAL(r_std), SAVE :: min_growthinit_time = 300.  !! minimum time since last beginning of a growing season (days)
1053!$OMP THREADPRIVATE(min_growthinit_time)
1054  REAL(r_std), SAVE :: moiavail_always_tree = 1.0  !! moisture monthly availability above which moisture tendency doesn't matter
1055                                                   !!  - for trees (0-1, unitless)
1056!$OMP THREADPRIVATE(moiavail_always_tree)
1057  REAL(r_std), SAVE :: moiavail_always_grass = 0.6 !! moisture monthly availability above which moisture tendency doesn't matter
1058                                                   !! - for grass (0-1, unitless)
1059!$OMP THREADPRIVATE(moiavail_always_grass)
1060  REAL(r_std), SAVE :: t_always                    !! monthly temp. above which temp. tendency doesn't matter
1061!$OMP THREADPRIVATE(t_always)
1062  REAL(r_std), SAVE :: t_always_add = 10.          !! monthly temp. above which temp. tendency doesn't matter (C)
1063!$OMP THREADPRIVATE(t_always_add)
1064
1065  ! 3. Coefficients of equations
1066 
1067  REAL(r_std), SAVE :: gddncd_ref = 603.           !!
1068!$OMP THREADPRIVATE(gddncd_ref)
1069  REAL(r_std), SAVE :: gddncd_curve = 0.0091       !!
1070!$OMP THREADPRIVATE(gddncd_curve)
1071  REAL(r_std), SAVE :: gddncd_offset = 64.         !!
1072!$OMP THREADPRIVATE(gddncd_offset)
1073
1074
1075  !
1076  ! stomate_prescribe.f90
1077  !
1078
1079  ! 3. Coefficients of equations
1080
1081  REAL(r_std), SAVE :: bm_sapl_rescale = 40.       !!
1082!$OMP THREADPRIVATE(bm_sapl_rescale)
1083
1084
1085  !
1086  ! stomate_resp.f90
1087  !
1088
1089  ! 3. Coefficients of equations
1090
1091  REAL(r_std), SAVE :: maint_resp_min_vmax = 0.3   !!
1092!$OMP THREADPRIVATE(maint_resp_min_vmax)
1093  REAL(r_std), SAVE :: maint_resp_coeff = 1.4      !!
1094!$OMP THREADPRIVATE(maint_resp_coeff)
1095
1096
1097  !
1098  ! stomate_soilcarbon.f90
1099  !
1100
1101  ! 2. Arrays
1102
1103  ! 2.1 frac_carb_coefficients
1104
1105  REAL(r_std), SAVE :: frac_carb_ap = 0.004  !! from active pool: depends on clay content  (0-1, unitless)
1106                                             !! corresponding to frac_carb(:,iactive,ipassive)
1107!$OMP THREADPRIVATE(frac_carb_ap)
1108  REAL(r_std), SAVE :: frac_carb_sa = 0.42   !! from slow pool (0-1, unitless)
1109                                             !! corresponding to frac_carb(:,islow,iactive)
1110!$OMP THREADPRIVATE(frac_carb_sa)
1111  REAL(r_std), SAVE :: frac_carb_sp = 0.03   !! from slow pool (0-1, unitless)
1112                                             !! corresponding to frac_carb(:,islow,ipassive)
1113!$OMP THREADPRIVATE(frac_carb_sp)
1114  REAL(r_std), SAVE :: frac_carb_pa = 0.45   !! from passive pool (0-1, unitless)
1115                                             !! corresponding to frac_carb(:,ipassive,iactive)
1116!$OMP THREADPRIVATE(frac_carb_pa)
1117  REAL(r_std), SAVE :: frac_carb_ps = 0.0    !! from passive pool (0-1, unitless)
1118                                             !! corresponding to frac_carb(:,ipassive,islow)
1119!$OMP THREADPRIVATE(frac_carb_ps)
1120
1121  ! 3. Coefficients of equations
1122
1123  REAL(r_std), SAVE :: active_to_pass_clay_frac = 0.68  !! (0-1, unitless)
1124!$OMP THREADPRIVATE(active_to_pass_clay_frac)
1125  !! residence times in carbon pools (days)
1126  REAL(r_std), SAVE :: carbon_tau_iactive = 0.149   !! residence times in active pool (days)
1127!$OMP THREADPRIVATE(carbon_tau_iactive)
1128  REAL(r_std), SAVE :: carbon_tau_islow = 7.0       !! residence times in slow pool (days)
1129!$OMP THREADPRIVATE(carbon_tau_islow)
1130  REAL(r_std), SAVE :: carbon_tau_ipassive = 300.   !! residence times in passive pool (days)
1131!$OMP THREADPRIVATE(carbon_tau_ipassive)
1132  REAL(r_std), SAVE, DIMENSION(3) :: flux_tot_coeff = (/ 1.2, 1.4, .75/)
1133!$OMP THREADPRIVATE(flux_tot_coeff)
1134
1135  !
1136  ! stomate_turnover.f90
1137  !
1138
1139  ! 3. Coefficients of equations
1140
1141  REAL(r_std), SAVE :: new_turnover_time_ref = 20. !!(days)
1142!$OMP THREADPRIVATE(new_turnover_time_ref)
1143  REAL(r_std), SAVE :: leaf_age_crit_tref = 20.    !! (C)
1144!$OMP THREADPRIVATE(leaf_age_crit_tref)
1145  REAL(r_std), SAVE, DIMENSION(3) :: leaf_age_crit_coeff = (/ 1.5, 0.75, 10./) !! (unitless)
1146!$OMP THREADPRIVATE(leaf_age_crit_coeff)
1147
1148
1149  !
1150  ! stomate_vmax.f90
1151  !
1152 
1153  ! 1. Scalar
1154
1155  REAL(r_std), SAVE :: vmax_offset = 0.3        !! minimum leaf efficiency (unitless)
1156!$OMP THREADPRIVATE(vmax_offset)
1157  REAL(r_std), SAVE :: leafage_firstmax = 0.03  !! relative leaf age at which efficiency
1158                                                !! reaches 1 (unitless)
1159!$OMP THREADPRIVATE(leafage_firstmax)
1160  REAL(r_std), SAVE :: leafage_lastmax = 0.5    !! relative leaf age at which efficiency
1161                                                !! falls below 1 (unitless)
1162!$OMP THREADPRIVATE(leafage_lastmax)
1163  REAL(r_std), SAVE :: leafage_old = 1.         !! relative leaf age at which efficiency
1164                                                !! reaches its minimum (vmax_offset)
1165                                                !! (unitless)
1166!$OMP THREADPRIVATE(leafage_old)
1167  !
1168  ! stomate_season.f90
1169  !
1170
1171  ! 1. Scalar
1172
1173  REAL(r_std), SAVE :: gppfrac_dormance = 0.2  !! report maximal GPP/GGP_max for dormance (0-1, unitless)
1174!$OMP THREADPRIVATE(gppfrac_dormance)
1175  REAL(r_std), SAVE :: tau_climatology = 20.   !! tau for "climatologic variables (years)
1176!$OMP THREADPRIVATE(tau_climatology)
1177  REAL(r_std), SAVE :: hvc1 = 0.019            !! parameters for herbivore activity (unitless)
1178!$OMP THREADPRIVATE(hvc1)
1179  REAL(r_std), SAVE :: hvc2 = 1.38             !! parameters for herbivore activity (unitless)
1180!$OMP THREADPRIVATE(hvc2)
1181  REAL(r_std), SAVE :: leaf_frac_hvc = 0.33    !! leaf fraction (0-1, unitless)
1182!$OMP THREADPRIVATE(leaf_frac_hvc)
1183  REAL(r_std), SAVE :: tlong_ref_max = 303.1   !! maximum reference long term temperature (K)
1184!$OMP THREADPRIVATE(tlong_ref_max)
1185  REAL(r_std), SAVE :: tlong_ref_min = 253.1   !! minimum reference long term temperature (K)
1186!$OMP THREADPRIVATE(tlong_ref_min)
1187
1188  ! 3. Coefficients of equations
1189
1190  REAL(r_std), SAVE :: ncd_max_year = 3.
1191!$OMP THREADPRIVATE(ncd_max_year)
1192  REAL(r_std), SAVE :: gdd_threshold = 5.
1193!$OMP THREADPRIVATE(gdd_threshold)
1194  REAL(r_std), SAVE :: green_age_ever = 2.
1195!$OMP THREADPRIVATE(green_age_ever)
1196  REAL(r_std), SAVE :: green_age_dec = 0.5
1197!$OMP THREADPRIVATE(green_age_dec)
1198
1199END MODULE constantes_var
Note: See TracBrowser for help on using the repository browser.