GDAY model

GDAY (Generic Decomposition And Yield) is a simple ecosystem model that simulates carbon, nitrogen, and water dynamics at the stand scale (Comins and McMurtrie, 1993; Medlyn et al. 2000; Corbeels et al. 2005a,b).

The model can be run at either a daily time step, or a 30-minute time step. When the model is run at the sub-daily timescale, photosynthesis is calculated using a two-leaf (sunlit/shade) approximation (de Pury and Farquhar, 1997; Wang and Leuning, 1998), otherwise photosynthesis is calculated following Sands (1995;1996). The sub-daily approach (photosynthesis & leaf energy balance) mirrors MAESPA, without the complexity of the radiation treatment. In the standard model the water balance is represented simply, with two (fixed) soil water "buckets", which represent a top soil (e.g. 5 cm) and a larger root-zone. If you are using the sub-daily version, there is now the option to use a SPA-style representation of hydraulics. GDAY-SPA resolves multiple soil layers, soil and leaf water potential and limits gas exchange following the Emax approach (see MAESPA for more details).

GDAY uses a modified version of the CENTURY model to simulate soil carbon and nutrient dynamics (Parton et al. 1987; 1993).

Installation

To get the code your best route is probably to fork the repository, there is a nice explanation on github.

The model is coded entirely in C without any dependancies. The wrapper files for the example scripts and the script used to change parameter options, are written in python. The old python version is still online.

There is a Makefile in the src directory...

$ make clean ; make

The Makefile will need to be edited by hand to set the $ARCH flag, which sets the installation path. Currently it is hardwired to my computer.

Running the model

A simple model usage can be displayed by calling GDAY as follows:

$ gday -u

Presently, there are only two options which the user can set via the command line. All model options, including all of the model parameters are customisable via the parameter file.

To spin-up GDAY:

$ gday -s -p param_file.cfg

To run GDAY:

$ gday -p param_file.cfg

When the model is run it expects to find its "model state" (i.e. from a previous spin-up) in the parameter file. This state is automatically written the parameter file after the initial spin-up when the "print_options" flag has been set to "end", rather than "daily".

Parameter file

GDAY expects a parameter file to be supplied as an argument (-p filename) on the command line. Parameter files follow the standard .INI format, although only a simple INI parser has been coded into GDAY.

Parameter files are broken down into 6 section, namely [git], [files], [params], [control], [state] and [print]. The order of these sections shouldn't make any difference. The basic element contained in the parameter file is the key or property. Every key has a name and a value, delimited by an equals sign (=). The name appears to the left of the equals sign.

eac = 79430.0

As mentioned above, keys are grouped into sections:

[control]
alloc_model = allometric
deciduous_model = false

[files]
cfg_fname = params/NCEAS_DUKE_model_youngforest_amb.cfg
met_fname = met_data/DUKE_met_data_amb_co2.csv
out_fname = outputs/D1GDAYDUKEAMB.csv
out_param_fname = params/NCEAS_DUKE_model_simulation_amb.cfg

As all the model parameters are accessible via this file, these files can be quite long. Clearly it isn't necessary to list every parameter. The recommended approach is to use the base file and then customise whichever parameters are required via a shell script, e.g. see the python wrapper script. This file just lists the parameters which needs to be changed and calls adjust_gday_param_file.py to swap the parameters (listed as a python dictionary) with the default parameter. I have also written an equivalent version in R adjust_gday_param_file.R. I should highlight that I wouldn't necessarily trust the default values :).

Finally, the options to print different the state and flux variables on the fly is a nice hangover from the python implementation. Sadly, this functionality doesn't actually exist in the C code, instead all the state and flux variables used in the FACE intercomparisons are dumped as standard.

When I have time I will write something more extensive (ha), but information about what different variable names refer to are listed in the header file, which documents the different structures (i.e. control, state, params).

The git hash allows you to connect which version of the model code produced which version of the model output. I'd argue for maintaining this functionality, but if you don't use git or wish to ignore me, filling this line with gibberish and disabling the shell command in the Makefile should allow you to do this.

Potential gotchas

The parameter alpha_j which represents the quantum yield of electron transport (mol mol-1) is the intrinsic quantum yield (i.e. per unit APAR). For the two-leaf version of the model, alpha_j should be divided by (1.0 - omega), where omega is the leaf scattering coefficient of PAR (leaf reflectance and transmittance combined). Currently we are assuming omega = 0.15 (radiation.c), this is currently hardwired.
The deciduous phenology scheme does not currently work with the two-leaf version of the model (can be fixed).
Wind speed must be > 0.0, some flux files have timesteps where this isn't the case and it will crash the model. Similarly, with the Medlyn gs model, VPD must be greater than 0.05 kPa.

Meteorological driving file

30-minute file:

Variable	Description	Units
year
doy	day of year	[1-365/6]
hod	hour of day	[0-47]
rain	rainfall	mm 30 min^-1
par	photosynthetically active radiation	umol m^-2 s^-1
tair	air temperature	deg C
tsoil	soil temperature	deg C
vpd	vapour pressure deficit	kPa
co2	CO₂ concentration	ppm
ndep	nitrogen deposition	t ha^-1 30 min^-1
nfix	biological nitrogen fixation	t ha^-1 30 min^-1
wind	wind speed	m s^-1
press	atmospheric pressure	kPa

Day file:

Variable	Description	Units
year
doy	day of year	[0-365/6]
tair	(daylight) air temperature	deg C
rain	rainfall	mm day^-1
tsoil	soil temperature	deg C
tam	morning air temperature	deg C
tpm	afternoon air temperature	deg C
tmin	minimum (day) air temperature	deg C
tmax	minimum (day) air temperature	deg C
tday	day average air temperature (24 hrs)	deg C
vpd_am	morning vapour pressure deficit	kPa
vpd_pm	afternoon vapour pressure deficit	kPa
co2	CO₂ concentration	ppm
ndep	nitrogen deposition	t ha^-1 day^-1
nfix	biological nitrogen fixation	t ha^-1 day^-1
wind	wind speed	m s^-1
press	atmospheric pressure	kPa
wind_am	morning wind speed	m s^-1
wind_pm	afternoon wind speed	m s^-1
par_am	morning photosynthetically active radiation	MJ m^-2 d^-1
par_pm	afternoon photosynthetically active radiation	MJ m^-2 d^-1

Nitrogen inputs

Nitrogen (N) entering the system via biological N fixation (BNF; tonnes ha^-1 yr^-1) and N deposition (tonnes ha^-1 yr^-1) are prescribed and passed via the met file. If information isn't available from the experiment GDAY is being applied to, BNF can be calculated as a function of evapotranspiration (ET) based on Cleveland et al. 1999.

Following Smith et al. (2014), Biogeosciences and Wieder et al. (2015), ERL, we also suggest the conservation BNF equation (Fig. 1). For estimates of ET you can either use values by PFT based on Table 1 in Cleveland or use the sum of canopy evaporation and transpiration. Wieder et al (pg 3) argued that using total ET leads to a high bias in BNF estimates in arid regions.

BNF (kg N ha-1 yr-1) is then calculated as a function of ET:

BNF = 0.102 * (ET * mm_2_cm) + 0.524

Hydraulics

From SPA we borrow the multi-layer soil scheme, which considers infiltration and drainage between layers. We also implement the soil-to-leaf hydraulics from SPA, which includes weighting soil water potential. We limit gas exchange following the Emax approach (Duursma et al. 2008). This approach therefore assumes isohydric behaviour, i.e. the plant maintains a leaf water potential above a critical minimum value.

We do not currently implement the thermal calculations which would allow you to estimate soil temperature.

Model spinup

There are now two options to spin the model up: (i) brute force, i.e. continusly recycling the met data and re-running GDAY (Flag - spinup_method = brute); and (ii) an implementation of the Semi-Analytical Solution (SAS) for steady-state approximation (Flag - spinup_method = sas). The SAS approach is about 60% faster than the slow spinup so I'd suggest using it.

The SAS method is after:

Xia, J., Y. Luo, Y.-P. Wang, and O. Hararuk. 2013. Traceable components of terrestrial carbon storage capacity in biogeochemical models. Global Change Biology 19:2104-2116

Example run

The example directory has two python scripts which provide an example of how one might set about running the model. example.py simulates the DUKE FACE experiment and run_experiment.py is just nice a wrapper script around this which produces a plot at the end comparing the data to the observations.

cd example/
python run_experiment.py

This should pop a plot of NPP, LAI and transpiration onto your screen. This example tends to break from time to time when I change various options, so please let me know if it isn't working!

NB to use this wrapper script you will need to have an installation of the Pandas and Matplotlib libraries installed. If you are a python user this is fairly standard.

Key References

Comins, H. N. and McMurtrie, R. E. (1993) Long-Term Response of Nutrient-Limited Forests to CO2 Enrichment; Equilibrium Behavior of Plant-Soil Models. Ecological Applications, 3, 666-681.
Medlyn, B. E., McMurtrie, R. E., Dewar, R. C. and Jeffreys, M. P. (2000), Soil processes dominate the long-term response of forest net primary productivity to increased temperature and atmospheric CO2 concentration, Canadian Journal of Forest Research, 30, 873–888.
Corbeels M, McMurtrie RE, Pepper DA, O’Connell AM (2005a) A process-based model for nitrogen cycling in forest planta- tions Part I. Structure, calibration and analysis of decomposi- tion model. Ecological Modelling, 187, 426–448.
Corbeels M, McMurtrie RE, Pepper DA, O’Connell AM (2005b) A process-based model for nitrogen cycling in forest planta- tions Part II. Simulating growth and nitrogen mineralisation of Eucalyptus globulus plantations in south-western Australia. Ecological Modelling, 187, 449–474.
Sands PJ (1995) Modelling canopy production. II. From single-leaf photosynthetic parameters to daily canopy photosynthesis. Australian Journal of Plant Physiology, 22, 603-614.
Sands PJ (1996) Modelling canopy production. III. Canopy light-utilisation efficiency and its sensitivity to physiological environmental variables. Australian Journal of Plant Physiology, 23, 103-114.
de Pury, D.G.G., Farquhar, G.D. (1997) Simple scaling of photosynthesis from leaves to canopies without the errors of big-leaf models. Plant, Cell and Environment, 20, 537-557.
Wang, Y-P. and Leuning R. (1998) A two-leaf model for canopy conductance, photosynthesis and portioning of available energy I: Model description and comparison with a multi-layered model. Agricultural and Forest Meteorology, 91, 89–111.
Parton, W.J., Schimel, D.S., Cole, C.V., and Ojima, D.S. (1987) Analysis of factors controlling soil organic matter levels in Great Plains grasslands. Soil Sc. Soc. Am. J. 51: 1173–1179.
Parton, W.J., Scurlock, J.M.O., Ojima, D.S., Gilmanov, T.G., Scholes, R.J., Schimel, D.S., Kirchner, T., Menaut, J.-C., Seastedt, T., Garcia Moya, E. Kamnalrut, A., and Kinyamario, J.I. (1993) Observations and modeling of biomass and soil or- ganic matter dynamics for the grassland biome worldwide. Global Biogeochem. Cycles, 7: 785–809.
Duursma, RA and Kolari, P and Perämäki, M and Nikinmaa, E and Hari, P and Delzon, S and Loustau, D and Ilvesniemi, H and Pumpanen, J and Mäkelä, A. (2008) Predicting the decline in daily maximum transpiration rate of two pine stands during drought based on constant minimum leaf water potential and plant hydraulic conductance. Tree physiology, 28, 265-276.
Duusma, R. A. and Medlyn, B. E. (2012) MAESPA: a model to study interactions between water limitation, environmental drivers and vegetation function at tree and stand levels, with an example application to [CO2] x drought interactions. Geoscientific Model Development, 5, 919-940.

Jinyan-Yang/GDAY