DESCRIPTION OF THE MODEL AND APPROACHES TO MODELING N DYNAMICS

The nitrogen dynamics routines of the CERES models were designed to simulate each of the major N loss processes and the contributions to the N balance made by mineralization. The routines also describe the uptake of N by the crop and the effects of N deficiency on crop growth processes. The transformations simulated are mineralization and/or immobilization, nitrification, denitrification, and urea hydrolysis. Nitrate movement associated with water movement in both an upward and downward direction is also simulated. Since the rates of transformation of nitrogen are very much influenced by soil water status, the simulation of nitrogen dynamics requires that water balance also be simulated. Soil temperature greatly influences many of the transformation rates. Therefore, a procedure to calculate soil temperature at various depths, based on the soil temperature routine of the EPIC model (Williams et al. 1984), is also invoked in the nitrogen component of the model.

The model does not simulate losses by ammonia volatilization or ammonium exchange equilibria and fixation. Under conditions of good fertilizer practice where fertilizer is either incorporated or placed beneath the soil surface, volatile ammonia losses should be small.

INITIALIZATION

Inputs describing the amount of organic matter and the amount of mineral nitrogen present in the soil are required to initialize the model. The model requires the organic carbon concentration in each layer (OC(L)) as an input and using an assumed soil C:N ratio of 10:1 calculates the amount of organic N associated with this organic matter (HUM(L)). These initializations are performed in subroutine SOILNI. To determine the contribution of recent crop residues to the supply of nitrogen in the soil, the model also requires an estimate of the amount of crop residue (STRAW) which is present. Based on this estimate and the depth of incorporation (SDEP) of the crop residue, the fresh organic matter content of each layer (FOM(L)) is estimated. An estimate of the amount of root residue remaining from the previous crop is also required for the calculation of FOM(L). Procedures for estimating STRAW and ROOT and the default values used are described in Chapter 6. Initial partitioning of the fresh organic matter into the component pools of carbohydrate (FPOOL(L,1)), cellulose (FPOOL(L,2)) and lignin (FPOOL(L,3)) is also performed in subroutine SOILNI.

NITRATE FLUX

Leaching of nitrates is probably the most common and best understood N loss process. Nitrates leaching from soil often become a source of contamination of groundwater and has recently generated interest in leaching from an environmental standpoint. There have been many approaches to modelling leaching based on numerical techniques which require solution in a manner inappropriate for use in a management level model such as CERES. In the CERES model, leaching is simulated using a simple approach based on the cascading system for drainage described in the previous chapter. Nitrate N may move between layers of the soil profile in the CERES models, but the movement of ammonium is not considered. Nitrate flux calculations are performed in subroutine NFLUX. Nitrate movement in the soil profile is highly dependent upon water movement. Therefore, the volume of water present in each layer (SW(L) * DLAYR(L)) and the water draining from each layer ((FLUX(L)) in the profile is used to calculate the nitrate lost from each layer (NOUT) as follows:

NOUT = SNO3(L) * FLUX(L)/(SW(L) * DLAYR(L) + FLUX(L))

A fraction of the mass of nitrate (SNO3(L)) present in each layer thus moves with each drainage event. A simple cascading approach is used where the nitrate lost from one layer is added to the layer below. When the concentration of nitrate in a layer falls to 0.25 g NO3 per g of soil no further leaching from that layer is deemed to occur. The method used may be termed a "reservoir mixing model" and is similar to the approach used by Burns (1974), but water movement is controlled by the SWCON variable in the drainage routine. The implicit assumption is that all the nitrate present in a layer is uniformly and instantaneously in solution in all of the water in the layer. Thus no attempt is made to separate nitrate in solution between the retained water and the mobile water. Differences in the relative volumes of retained water and mobile water between clays and sands occur as a function of the relative magnitudes of LL(L), DUL(L), and SAT(L). The rate of nitrate flux is also sensitive to changes in SWCON since this variable determines the rate of drainage. Nitrate is more readily displaced from sands since the volume of water which can move ((SAT(L) - DUL(L)) * DLAYR(L)) is large in comparison to the retained water (DUL(L) * DLAYR(L)). Most of the difference in the simulated leaching rate between soils of different texture is explained by this difference in proportion of water which is mobile. Some difference is also attributable to the rate at which the soil profile can drain (SWCON). The upward flow of water in the top four soil layers will also cause some redistribution of nitrate. A second loop, commencing in the deepest layer of evaporative water loss (MU), is used to calculate this redistribution. Nitrate moving from a layer (NUP) is calculated as a function of upward movement (FLOW(L)) in a manner identical to leaching:

NUP = SNO3(L) * FLOW(L)/(SW(L) * DLAYR(L) + FLOW(L)) * 0.5

No upward loss from the top layer occurs by this process. Since there will occasionally be instances when this slowly moving water can move in a downward direction (negative values of FLOW(L)) a third loop is set up with calculations commencing in the top layer and running to lower layers. This is achieved by first reinitializing the array FLUX to 0 and reversing the sign (to make it positive) at the FLOW array and copying it to FLUX. When this has been done the normal leaching calculations used in the first loop can be used again. These instances would occur when a small rainfall wets the top layer of a very dry soil profile. There may have been insufficient water for drainage to occur but a moisture potential between the top layer and the second layer initiates this flow. The resultant movement of nitrate will be very small.

SOIL NITROGEN TRANSFORMATIONS

The CERES model simulates the decay of organic matter and the subsequent mineralization and/or immobilization of N, the nitrification of ammonium and denitrification in subroutine NTRANS. Fertilizer addition and transformations (assumed to be instantaneous) are also performed in this subroutine.

Fertilizer Additions

Fertilizer N is partitioned in the model between nitrate and ammonium pools according to the nature of the fertilizer used. Fertilizer products are specified by a numeric code IFTYPE. In addition to the numeric code for fertilizer type, inputs required to describe the fertilizer are: the date of application (FDAY), the amount of N applied (AFERT) and the depth of placement (DFERT). For any placement depth the assumption is made that the fertilizer is uniformly incorporated into the layer. Layer thicknesses are supplied as input and are usually based on natural horizonation in the profile. These must correspond with those used to describe the soil water inputs. Surface fertilizer applications are treated as being uniformly incorporated into the top layer. Up to 10 split applications can be accommodated by the model.

Mineralization and Immobilization

Mineralization refers to the net release of mineral nitrogen with the decay of organic matter and immobilization refers to the transformation of mineral nitrogen to the organic state. Both processes are microbial in origin. Immobilization occurs when soil microorganisms assimilate inorganic N compounds and utilize them in the synthesis of the organic constituents of their cells. A balance exists between the two processes. When crop residues with a high C:N ratio are added to soil, the balance can shift resulting in net immobilization for a period of time. After some of the soil carbon has been consumed by respiration, net mineralization may resume. N mineralized from the soil organic pool can often constitute a large part of the nitrogen available to the crop.

The perceived application for the CERES models in studies examining crop growth and fertilizer management requires that a mineralization model be simple, require few inputs, and work on a diversity of soils. Simulation studies examining the affects of crop residues also requires that the model be capable of simulating the fate of residues of different compositions. Other studies examining the potential role of nitrification inhibitors require a model wherein the processes of ammonification and nitrification are separated. The approach used in the CERES-WHEAT model is based on a modified version of the mineralization and immobilization component of the PAPRAN model (Seligman and van Keulen, 1981). This model is an attempt at maintaining some of the functionality of the microbiological level models but doing so at a very simplified level. The model's modifications have been to simulate nitrification separately and to partition the simulated fresh organic matter pools differently. Modifications were also made to temperature and water indices to fit the CERES water balance and soil temperature routines. Unless otherwise indicated, the coefficients used for the mineralization/immobilization functions described below were drawn from the PAPRAN model.

The mineralization and immobilization routine simulates the decay of two types of organic matter: Fresh organic matter (FOM) which comprises crop residues or green manure and a stable organic or humic pool (HUM). Three pools comprise the FOM pool in each layer (L), vis:

FPOOL(L,1) = carbohydrate

FPOOL(L,2) = cellulose

FPOOL(L,3) = lignin.

In PAPRAN, FOM is simulated as one pool and the decay rate constant is selected according to the proportion of the initial amount of FOM remaining. The CERES model separates FOM into three pools giving a better estimate of soluble carbon which is used in the denitrification routine. These three pools are initialized as a fraction of the FOM(L) pool in subroutine SOILNI. Initially, the FOM(L) contains 20% carbohydrate, 70% cellulose and 10% lignin. The model requires as input data, the amount of straw added, its C:N ratio and its depth of incorporation (if any) and an estimate of the amount of root residue from the previous crop. Based upon these data, initial values of FOM and the N contained within it (FON) for each layer are calculated in subroutine SOILNI. The soil organic carbon in each layer (OC(L)) is also required by the mineralization routine. This is used to calculate HUM(L), and together with a simplifying assumption of a bulk soil C:N ratio of 10, is used to estimate the N associated with this fraction (NHUM(L)). Each of the three FOM pools (FPOOL (L,1 to 3)) has a different decay rate (RDECR

(1 to 3)). Under nonlimiting conditions the decay constants as reported by Seligman and van Keulen (1981) are 0.80, 0.05, and 0.0095 for carbohydrate, cellulose, and lignin, respectively. A decay constant at 0.20 for the carbohydrate fraction has since been found to be more appropriate. The decay constant for carbohydrate implies that under nonlimiting conditions 20% of the pool will decay in one day. Nonlimiting conditions very seldom occur in soils since one or all of soil temperature, soil moisture, or residue composition will limit the decay process. To quantify these limits three zero to unity dimensionless factors are calculated. A water factor (MF) is first determined from the volumetric soil water content (SW(L)) relative to the lower limit (LL), and drained upper limit (DUL). In accordance with the soil water balance model, provision is made for the water content of the uppermost layer to be lower than the lower limit. The variable SWEF determines the lowest possible value the uppermost soil layer water content have. When the soil is drier than DUL, MF is calculated as:

AD = LL(L)

IF (L.EQ.1) AD = LL(1) *SWEF

MF = (SW(L)-AD)/(DUL(L)-AD)

where:

AD = lowest moisture content for a layer (volume fraction)

When the soil is wetter than DUL, MF is calculated as:

MF = 1.0-(SW(L)-DUL(L))/(SAT(L)-DUL(L))*0.5

The functions follow the observations reported by Myers et al. (1982)

and Linn and Doran (1984) on moisture effects on ammonification. Under very wet conditions (100% of water filled porosity) ammonification proceeds at approximately half of the rate of ammonification at field capacity (Linn and Doran, 1984). The comparative effects of soil moisture on the simulated rates of ammonification, nitrification, and denitrification can be seen in Fig. 5.1. A temperature factor (TF) is calculated directly from soil temperature (ST(L)):

TF = (ST(L)-5.0)/30.0

This approximates the soil temperature effects on ammonification reported by others (Stanford et al., 1973; Myers, 1975). If the soil temperature (ST(L)) is less than 5^o C then TF is set to zero and no decay occurs. The C:N ratio (CNR) imposes the third limit on decay rate. In this case C:N ratio is calculated as the C contained in FOM divided by the N "available" for the decay process. This N available for decay is the sum of the N contained in the FOM, which is FON, and the extractable mineral N present in the layer (TOTN). Thus,

CNR = (0.4*FOM(l))/(FON(L)+TOTN)

From CNR an index (CNRF) is calculated which has a critical C:N ratio of 25.

CNRF = EXP(-0.693*(CNR-25)/25.0)

Thus, in low N containing residues (e.g., freshly incorporated wheat straw) with a high C:N ratio, the N available for the decay process will greatly limit the decay rate (Fig. 3.5). For each of the FOM pools a decay rate (GRCOM) appropriate for that pool (JP) can be calculated. G1 = TF*MF*CNRF*RDECR(JP)

GRCOM = G1*FPOOL1(L,JP)

The gross mineralization of N associated with this decay (GRNOM) is then calculated according to the proportion of the pool which is decaying.

GRNOM = G1 * FPOOL(L,JP)/FOM(L) * FON(L)

GRCOM and GRNOM are summed for each of three pools in each layer. The procedure used for calculating the N released from the humus (RHMIN) also utilizes TF and MF. In this case CNRF is not used and the potential decay rate constant (DMINR) is very small (8.3E-05). A further index (DMOD) was added to the RHMIN calculations to adjust the mineralization rate for certain atypical soils. On soils with chemically protected organic matter, a less than unity value of DMOD is required so that mineralization is not overestimated. On freshly cultivated virgin soils, a slightly greater than unity value has been found necessary to account for the sudden increase in mineralization activity. In all other circumstances a value of 1.0 is used for DMOD. Satisfactory alternatives for estimating DMOD are currently being sought. The procedure for calculating RHMIN, then is the product of the various indices and the N contained within the humus (NHUM(L)).

RHMIN=NHUM(L)*DMINR*TF*MF*DMOD

After calculating the gross mineralization rate, HUM(L) and NHUM(L) are updated.

HUM(L)=HUM(L)-RHMIN*10.0+0.2*GRNOM/0.04

NHUM(L)=NHUM(L)-RHMIN+0.2*GRNOM

These calculations also allow for the transfer of 20% of the gross amount of N released by mineralization of FON(L) (0.2*GRNOM) to be incorporated into NHUM(L). This accounts for N incorporated into microbial biomass and has a concentration of 4% (0.04) determined as 0.1 g N/g C (soil C:N ratio of 10) multiplied by 0.4 g C/g OM (40% of OM is C). As organic matter decomposes some N is required by the decay process and may be incorporated into microbial biomass. The N which is immobilized in this way (RNAC) is calculated as the minimum of the soil extractable mineral N (TOTN) and the demand for N by the decaying FOM(L).

RNAC=AMIN1(TOTN,GRCOM*(0.02-FON(L)/FOM(L))

where 0.02 is the N requirement for microbial decay of a unit of FOM(I). The value of 0.02 is the product of the fraction of C in the FOM(L) (40%), the biological efficiency of C turnover by the microbes (40%) and the N:C ratio of the microbes (0.125). FOM(L) and FON(L) are then updated.

FOM(L)=FOM(L)-GRCOM

FON(L)=FON(L)+RNAC-GRNOM

The balance between RNAC and GRNOM determines whether net mineralization or immobilization occurs. The net N released from all organic sources (NNOM) is:

NNOM=0.8*GRNOM+RHMIN-RNAC.

Note that only 80% of GRNOM enters this pool since the remaining 20% was incorporated into NHUM(L). NNOM can then be used to update the ammonium pool (SNH4(L)).

SNH4(L)=SNH4(L)+NNOM

If net immobilization occurs (NNOM negative) ammonium is first immobilized and if there is not a sufficient amount to retain this pool with a concentration of 0.5 ppm, withdrawals are made from the nitrate pool.

Nitrification

Nitrification refers to the process of oxidation of ammonium to nitrate. It is a biological process and occurs under aerobic conditions. The main factors which limit nitrification are: Substrate NH4+, oxygen, soil pH, and temperature. The approach used in the CERES models has been to calculate a potential nitrification rate and a series of zero to unity environmental indices to reduce this rate. This potential nitrification rate is a Michaelis-Menten kinetic function dependent only on ammonium concentration and is thus independent of soil type. A further index, termed a "nitrification capacity" index, is introduced which was designed to introduce a lag effect on nitrification if conditions in the immediate past (last 2 days) have been unfavorable for nitrification. Actual nitrification capacity is calculated by reducing the potential rate by the most limiting of the environmental indices and the capacity index. The capacity index is an arbitrary term introduced to accommodate an apparent lag in nitrification observed in some data sets. The functions reported below were found to be appropriate across the range of data sets tested. The nitrification routine in subroutine NTRANS calculates the nitrification of ammonium in each layer. First, an ammonium concentration factor (SANC) is calculated.

SANC=1.0-EXP(-0.01363*SNH4(L))

This is a zero to unity index which has approximately zero values when there is less than 1 ppm of ammonium present and has a value of 0.75 at 100 ppm. The temperature factor calculated above for mineralization (TF) and a soil water factor for nitrification (WFD) (Fig. 3.4) are used together with SANC to determine an environmental limit on nitrification capacity (ELNC).

ELNC=AMIN1(TF,WFD,SANC)

To accommodate lags which occur in nitrifier populations ELNC and the previous day's relative microbial nitrification potential in the layer (CNI(L)) are used to calculate the interim variable RP2 which represents the relative nitrification potential for the day.

RP2=CNI(L)*EXP(2.302*ELNC)

RP2 is constrained between 0.01 and 1.0. Today's value of the nitrification potential (CNI(L)) is then set equal to RP2. Since EXP(2.302*ELNC) varies from 1.0 to 10.0 when ELNC varies from 0.0 to 1.0, relative nitrification potential can increase rapidly, up to tenfold per day. An interim variable A is then determined from these indices and an index for pH effect on nitrification. This pH index is calculated in subroutine SOILNI and represents the conclusions drawn by Schmidt (1982) on the pH effect in nitrification.

A=AMIN1(RP2,WFD,TF,PHN(L))

This interim variable A is used together with the ammonium concentration (NH(L)) in a Michaelis-Menten function described by McLaren (1970) to estimate the rate of nitrification. The function has been modified to estimate the proportion of the pool of ammonium (SNH4(L)) which is nitrified on a day.

B=(A*40.0*NH4(L)/NH4(L)+90.0))*SNH4(L)

A maximum of 80% at the ammonium pool is allowed to nitrify in one day. A check is made to ensure some ammonium is retained in the layer and thus the daily rate of nitrification (RNTRF) is

RNTRF=AMIN1(B,SNH4(L))

Following this calculation, soil nitrate and ammonium pools can be updated.

SNH4(L)=SNH4(L)-RNTRF

SNO3(L)=SNO3(L)+RNTRF

Finally, the soil temperature, moisture and NH4 after nitrification are used to update (CNI(L)), which is used in the subsequent day's calculations.

SARNC=1.0-EXP(-0.1363*SNH4(L))

XW=AMAX1(WF,WFY(L))

XT=AMAX1(TF,TFY(L))

CNI(L)=CNI(L)*AMIN1(XW,XT,SARNC)

SARNC is a zero to unity factor for ammonium availability. WFD and WFY(L) are today's and yesterday's soil water factors, respectively, and TF and TFY(L) are today's and yesterday's soil temperature factors, respectively. The least limiting of the current day's and the previous day's water and temperature factors are used in the calculation of the new value of CNI(L). This prevents a single day of low soil temperature or water from severely reducing CNI(L). It is important to note that the relative nitrification potential CNI(L) is calculated twice each day. Since (EXP(2.302*ELNC)) varies from 1.0 to 10.0, CNI(L) increases prior to the calculation of the nitrification rate. After the nitrification calculations when the level of ammonium has declined, CNI(L) is reduced. The relative magnitudes of (EXP(2.302*ELNC) and AMIN1(XW,XT,SARNC)) determine whether relative nitrification potential increases or decreases over the short term.

Denitrification

Denitrification is the dissimilatory reduction of nitrate (or nitrite) to gaseous products including N0, N20, and N2 (Knowles, 1981).

Denitrification is a microbial process which occurs under anaerobic conditions and is influenced by organic carbon content, soil aeration, temperature and soil pH. The approach adopted in the CERES models has been to adapt the functions described by Rolston et al. (1980) to fit within the framework of the model and to match inputs derived from the water balance and mineralization components of CERES. The basic function used by these authors was also used by Davidson et al. (1978a) and was the subject of field testing under a variety of conditions in California. Predicted rates of denitrification compared favorably with direct measures of gaseous losses in the field experiments. Denitrification calculations are only performed when the soil water content (SW) exceeds the drained upper limit (DUL). A zero to unity index (FW) (see Fig. 5.1) for soil water in the range from DUL to saturation (SAT) is calculated.

FW = 1.0 - (SAT(L)-SW(L))/(SAT(L)-DUL(L))

Linn and Doran (1983) used percentage of water filled porosity as an index of soil water availability effects on soil N transformations. In their studies, denitrification commenced with a water-filled porosity of 60% and increased linearly up to 100% water filled porosity. This approximates the linear increase in FW as SW increases from DUL to SAT. A factor for soil temperature (FT) is also calculated.

FT=0.1*EXP(0.046*ST(L))

Rolston et al. (1980) using the data of Burford and Bremner (1976) and Reddy et al. (1971) to estimate the water-extractable C in soil organic matter (CW) as:

CW=24.5 + 0.0031*SOILC

In the CERES model, SOILC is calculated as 58% of the stable humic fraction. To this is added the carbon contained in the carbohydrate fraction organic matter pool (40% of FPOOL(L,1)). Appropriate unit conversions are made using FAC(L) and the total water extractable carbon (CW) estimated.

CW = FAC(L)*(SOILC*0.0031+0.4*FPOOL(L,1))+24.5

Denitrification rate (DNRATE) is then calculated from the nitrate concentration and converted to a kg N/ha basis for the mass balance calculations.

DNRATE = 6.0*1.0E-05*CW*NO3(L)*FW*FT*DLAYR(L)

Following the calculation of DNRATE the nitrate pool in the layer is updated with appropriate checks to ensure that a minimum concentration of nitrate is retained in the layer.

SNO3(L)=SNO3(L)-DNRATE

SOIL TEMPERATURE

The soil temperature in each layer is used in the functions describing most of the major soil N transformations. The soil temperature model used in CERES is based on that used in the EPIC model (Williams et

al., 1984). This method is based upon some simple empiricisms and requires only two additional inputs to those soil parameters required by the water balance and N transformation routines. These inputs are: TAV, the annual average ambient temperature and AMP the annual amplitude in mean monthly temperature. The method used to calculate the soil temperature at various depths in the profile requires the determination of a damping depth (the depth at which no diurnal variation in temperature is experienced). At depths more shallow than this, diurnal change in temperature occurs with the greatest fluctuation happening near the surface. The location of this damping depth (DD) is dependent upon parameters which influence the flux of heat in the soil, notably the bulk density and the moisture content. DD is updated daily to allow for changes in soil moisture content. Soil surface temperatures are modelled as a function of the ambient temperature, the solar radiation, and the albedo. The 5-day moving average surface temperature is used to compute the temperatures in each layer as follows:

TMA(1) = (1.0-ALBEDO) * (TEMPM + (TEMPMX - TEMPM) *

SQRT(SOLRAD * 0.03)) + ALBEDO * TMA(1)

where:

TMA(1) = Daily surface temperature

ALBEDO = The albedo of the soil surface and is an input variable for bare soils. As the crop canopy develops ALBEDO becomes a function of the leaf area. These calculations of albedo are performed in the water balance routine as they are a fundamental component of the evaporation model.

SOLRAD = Solar radiation in MJ/square metre.

TEMPMX,TEMPM = Daily maximum and mean temperature C, respectively.

The long-term average daily ambient temperature (TA) for the current day of the year can be estimated from TAV and AMP.

TA = TAV + AMP * COS(ALX)/2.0

ALX is a variable (in units of radians) to relate the current day of the year (XI) to the time of the hottest day of the year (HDAY). In the northern hemisphere this is assumed to be day 200 and in the southern hemisphere day 20.

ALX = (XI - HDAY) * 0.0174

The coefficient 0.0174 is 1/365 days multiplied by 2 radians. Deviations in the actual dates of the hottest day of the year in lower latitudes are of little importance since the volumes of AMP will be small and hence TA will approximate TAV. The departure (DT) of the moving average temperature from TA is used in the calculation of the soil temperature in each layer (ST(L)) as follows:

ST(L) = TAV + (AMP/2.0 * COS(ALX + ZD) + DT * EXP(ZD)

Where ZD = depth of layer L/current day's damping depth.

PLANT CRITICAL N CONCENTRATIONS AND N DEFICIT FACTORS

Plant growth is greatly affected by the supply of N. Typically the supply of N to plants at the beginning of the season is often relatively high and becomes lower as the plant reaches maturity. The concentration of N in plant tissues also changes as the plant ages. During early growth, N concentrations are usually high due to synthesis of large amounts of organic N compounds required by the biochemical processes constituting photosynthesis and growth. As the plant ages, less of this new material is required and export from old tissues to new tissues occurs lowering the whole plant N concentration. At any point in time there exists a critical N concentration in the plant tissue below which growth will be reduced. These concentrations are determined as a function of crop ontogenetic age and are used within the model as part of the procedure to simulate the effects of N deficiency. The model's critical concentration functions are based upon the often used Zadoks' growth scale (Zadoks et al. 1974). Zadoks' growth scale is a decimal index of crop development generalized for all cereals. The intervals between growth scale index values are based on crop morphological observations and are not related to a thermal time concept. To incorporate the Zadoks' scale, a scheme to provide a conversion between the integer growth stages recognized by the model (ISTAGE) and a functional form of the Zadoks' scale had to be devised. XSTAGE is a fractional growth stage which is used to determine an approximate value for the corresponding Zadoks' stage (ZSTAGE). The conversions were performed using several functions which are tabulated below (Table 1). The functions are located in subroutine NFACTO.

Table 1. Functions Used for Converting From Fractional Growth Stage (XSTAGE) to Zadoks'Growth Stage (ZSTAGE)

Morphological Stage XSTAGE Range Function

Emergence to terminal spikelet 0.0-2.0 ZSTAGE = XSTAGE

Terminal spikelet to booting 2.0-3.0 ZSTAGE = 2.0 + 2.0*(XSTAGE-2.0)

Booting to ear emergence 3.0-4.0 ZSTAGE = 4.0 + 1.7*(XSTAGE-3.0)

Ear emergence to anthesis 4.0-4.4 ZSTAGE = 5.7 + 0.8*(XSTAGE-4.0)

Anthesis to maturity 4.4-6.0 ZSTAGE = 6.02 + 1.86*(XSTAGE-4.4)

To develop appropriate relationships for critical N concentrations in wheat, published data from field experiments that met the following criteria were assembled:

1. Experiments had a series of N rates with sufficient range to define optimal or near-optimal growth patterns.

2. Experiments were considered to have been conducted under conditions where the potential effects of other interacting factors (e.g., heat stress, moisture stress, frost, supply of other nutrients, etc.) were minimized.

3. Plant tops N concentration was reported at several times during the growing season.

4. The growth stage or phenological age of the crop was reported for the times of plant sampling.

In some cases, critical concentrations were defined by the authors and where appropriate were adopted. In two studies only one N rate was used but was described as being an optimal rate by the authors. Data were drawn from the following sources (Table 2) representing a diversity of wheat genotypes and wheat-growing environments.

Table 2. Data Sources Used for Determination of Critical N Concentration Relationships

Author Spring or Winter Wheat Location

Engel and Zubriski (1982) Spring North Dakota

Campbell et al. (1977a) Spring Canada

Wagger et al. (1981) Winter Kansas

Leitch and Vaidanathan (1983) Winter U.K.

Wagger (1983) Winter Kansas

Waldren and Flowerday (1979) Winter Nebraska (?)

Page et al. (1977) Winter U.K.

Alessi et al. (1979) Spring North Dakota

Mugwira and Bishnoi (1980) Winter Alabama

Boatwright and Haas (1961) Spring North Dakota

Gasser and Thorburn (1972) Spring U.K.

Bhargava and Motiramani (1967) Spring Australia

Walia et al. (1980) Spring India

McNeal et al. (1968) Spring Montana

Spratt and Gasser (1970) Spring U.K.

From these data, relationships defining critical N concentration as a function of Zadoks' growth stage were determined. The critical N concentration was defined as the N concentration in the plant tissues at optimal or near optimal growth (as defined by biomass, yield or leaf area from the response data). The relationship thus determined is defined as the concentration above which no further increases in crop growth occur and below which some effect on a growth process will occur. Winter wheats and spring wheats were found to have different relationships (Fig. 3.7). The differences between winter and sprin wheats may be an artifact created by the different growing conditions of the experiments cited above. It has been difficult to characterize critical concentrations particularly for the period of rapid growth in the spring when phenological age, N uptake and biomass are all increasing rapidly. These relationships for the tops critical N percentage TCNP) appear in subroutine NFACTO as a function of Zadoks' growth stage (ZSTAGE).

For winter wheats:

TCNP = -5.0112 - 6.3507 * ZSTAGE + 14.9578 * SQRT(ZSTAGE) +

0.2238 * (ZSTAGE * ZSTAGE)

For spring wheats:

TCNP = 7.4532 - 1.7908 * ZSTAGE + 0.6093 * SQRT(ZSTAGE) +

0.0934 * ZSTAGE * ZSTAGE

Root critical N concentration (RCNP) relationships were derived from the greenhouse data of Peterson et al. (1983) and Day et al. (1985).

RCNP = 2.10 - 0.14 * SQRT(ZSTAGE)

The minimum concentration of N in plant tissues as a function of plant age is seldom reported. To formulate an appropriate relationship for use in the model, some of the minimum concentrations reported in the above studies were used as well as those reported from an extensive survey of N concentration in wheat crops spanning several years and locations in South Australia by Schultz and French (1976). In the model the tops minimum concentration (TMNC) is calculated as a function of model growth stage (XSTAGE):

TMNC = 2.97 - 0.455 * XSTAGE

Root critical minimum N concentration (RMNC) is used during the grain filling calculations (in subroutine GROSUB) and is assumed to be a constant 75% of the critical concentration.

RMNC = 0.75 * RCNP

The coupling of these functions to the phenology routines thus enables critical concentrations to be determined for any variety growing in any environment. The critical and minimum concentrations are used to define a nitrogen factor (NFAC) which ranges from zero to slightly above unity. NFAC is the primary mechanism used within the model to determine the effect of N on plant growth. It is an index of deficiency relating the actual concentration (TANC) to these critical concentrations. NFAC has a value of zero when TANC is at its minimum value of TMNC and increases to 1.0 as concentration increases toward the critical concentration. NFAC is calculated as:

NFAC = 1.0 - (TCNP - TANC)/(TCNP - TMNC)

Since all plant growth processes are not equally affected by N stress, a series of indices based on NFAC are used. For photosynthetic rate (NDEF1) the index is calculated as:

NDEF1 = 0.10 + 2.0 *NFAC (NDEF<1.0)

For leaf expansion growth (NDEF2) a more sensitive factor is used:

NDEF2 = NFAC

For tillering (NDEF3) the index is calculated as:

NDEF3=NFAC*NFAC

For the calculation of these indices NFAC has a maximum value of 1.0. This implies that when TANC exceeds TCNP no extra growth occurs. A fourth factor used to modify the rate of grain N accumulation (NDEF4) is also calculated from NFAC, and can range from 0.0 to 1.5.

NDEF4 = NFAC * NFAC

These relations are depicted in Fig. 3.8. In the growth subroutine, GROSUB, the law of the minimum is used extensively to modify rates of plant growth. For each of the major functions (e.g., photosynthetic rate, leaf expansion rate, tiller number determination) the minimum of several zero to unit stress indices is used to modify a potential rate for the process.

N UPTAKE

The approach used in the CERES models has been to separately calculate the components of demand and supply and then use the lesser of these two to determine the actual rate of uptake. Demand can be considered as having two components. First there is a "deficiency demand." This is the amount of N required to restore the actual N concentration in the plant (TANC for tops) to the critical concentration (TCNP for tops). Critical concentrations for shoots and roots are defined in section 3.7. This deficiency demand can be quantified as the product of the existing biomass and the concentration difference as below:

TNDEM = TOPWT * (TCNP - TANC)

Similarly for roots the discrepancy in concentration (difference between RCNP and RANC) is multiplied by the root biomass (RTWT) to calculate the root N demand.

RNDEM = RTWT * (RCNP-RANC)

If luxury consumption of N has occurred such that TANC is greater than TCNP then these demand components have negative values. If total N demand is negative then no uptake is performed on that day. The second component of N demand is the demand for N by the new growth. Here the assumption is made that the plant would attempt to maintain a critical N concentration in the newly formed tissues. To calculate the new growth demand, a potential amount of new growth is first estimated in the GROSUB subroutine. New growth is estimated from potential photosynthesis (PCARB) and is partitioned into a potential root growth (PGRORT) and a potential tops growth (PDWI). Partitioning between potential shoot and root growth occurs as a function of phenological age:

PGRORT = PCARB * (60 - XSTAGE * 8)/100

PDWI = PCARB - PGRORT

These potential growth increments provide a mechanism for the tops actual N concentration (TANC) to exceed TCNP. This occurs when some stress prevails and the actual growth increment is less than the potential. New growth demand for tops (DNG) is calculated as

DNG = PDWI * TCNP

and the new growth demand for roots is calculated as

PGRORT * RCNP.

During the early stages of plant growth the new growth component of N

demand will be a large proportion of the total demand. As the crop biomass increases the deficiency demand becomes the larger component. During grain filling, the N required by the grain is removed from the vegetative and root pools to form a grain N pool. The resultant lowering of concentration in these pools may lead to increased demand. The total plant N demand (NDEM) is the sum of all of these demand components. Calculations of soil supply of N are on a per hectare basis which necessitates recalculation of the per plant demand into a per hectare demand (ANDEM).

ANDEM = NDEM * PLANTS * 10.0

To calculate the potential supply of N to the crop, zero to unity availability factors for each of nitrate (FNO3) and ammonium (FNH4) are calculated from the soil concentrations of the respective ions:

FNO3 = 1.0 - EXP(-0.0275 * NO3(L))

FNH4 = 1.0 - EXP(-0.025 * NH4(L))

The coefficients used in these two functions, obtained by trial and error, were found to be appropriate over a range of data sets. The greater mobility of nitrate ions in soil is reflected by the larger coefficient (0.0275) in these equations. A zero to unity soil water factor (SMDFR) which reduces potential uptake is calculated as a function of the relative availability of soil water:

SMDFR = (SW(L) - LL(L)/ESW(L)

To account for increased anaerobiosis and declining root function at moisture contents above the drained upper limit, SMDFR is reduced as

saturation is approached.

IF(SW(L) . GT . DUL(L))SMDFR = 1.0 - (SW(L) - DUL(L))/(SAT(L) - DUL(L))

The maximum potential N uptake from a layer may be calculated as a function of the maximum uptake per unit length of root and the total amount of root present in the layer. The first of these is a temporary variable (RFAC) which integrates the effects of root length density (RLV(L)), the soil water factor described above, and the depth of the layer:

RFAC = RLV(L) * SMDFR * SMDRF * DLAYR(L) * 100.0

The second of these equations incorporates the ion concentration effect (FNO3) and the maximum uptake per unit length of root (0.009 kg N/ha cm root) to yield a potential uptake of nitrate from the layer (RNO3U(L)).

RNO3U(L) = RFAC * FNO3 * 0.009

(RNO3U(L)) is thus the potential uptake of nitrate from layer L in kg N/ha constrained by the availability of water, the root length density and the concentration of nitrate. Initial estimates for the maximum uptake per unit length of root coefficient were obtained from the maize root data of Warncke and Barber (1974). This estimate was the subject of continuing modification during early model development. The value reported here appears to be appropriate across a broad range of data sets. The effect of each of these parameters on determining potential uptake can be seen in Fig. 3.9. A similar function is employed to calculate the potential uptake of ammonium (RNH4U(L)).

RNH4U(L) = RFAC * FNH4 * 0.009

Potential N uptake from the whole profile (TRNU) is the sum of RNO3U(L) and RNH4U(L) from all soil layers where roots occur. Thus TRNU represents an integrated value which is sensitive to (a) rooting density, (b) the concentration of the two ionic species, and (c) their ease of extraction as a function of the soil water status of the different layers. This method of determining potential uptake enables the common condition, where N is concentrated in the upper layers of the profile, where most of the roots are present and where a nutritional drought due to shortage of water in these upper layers may occur, to be simulated. This can occur when the crops demand for water is satisfied from soil water located deeper in the profile but where there may be little N present. If the potential N supply from the whole profile (TRNU) is greater than the crop N demand (ANDEM) an N uptake factor (NUF) is calculated and used to reduce the N uptake from each layer to the level of demand.

NUF = ANDEM/TRNU

This could occur when plants are young and have a high N supply. If the demand is greater than the supply then NUF has a value of 1.0. When NUF is less than 1.0, uptake from each layer is reduced as follows:

UNO3 = RNO3U(L) * NUF

UNH4 = RNH4U(L) * NUF

Following these calculations the soil mineral N pools can be updated for the actual uptake which has occurred.

SNO3(L) = SNO3(L) - UNO3

SNH4(L) = SNH4(L) - UNH4

Under conditions of luxury N uptake (TANC > TCNP) exudation of organic N compounds can occur. Rovira (1969) found changes in the shoot environment which cause more rapid growth can increase exudation. Bowen (1969) reported that N deficiency can cause exudation to decrease. In the CERES-N model this exuded N is added to the fresh organic N pool (FON(L)) and can be mineralized and subsequently made available to the plant again. The amount of N which can be lost from the plant in this manner is calculated as 5% of the N contained in the roots/day. These losses are distributed to the FON(L)) pool according to the differing root length densities present in each layer as a proportion of the total root length.

IF(TANC . GT . TCNP)RNLOSS = RANC * RTWT * 0.05 * PLANTS * RLV(L)/TRLV

Following uptake, concentrations of N in each of the shoots and roots are updated. To do this TRNU is converted from kg N/ha to a g N/plant basis.

TRNU = TRNU/(PLANTS * 10.0)

The proportion of the total plant demand (NDEM) arising from shoots (TNDEM) and roots (RNDEM) and the total root N loss (TRNLOS) are used to calculate the changes in N content of the shoots (DTOPSN) and roots (DROOTN).

DTOPSN = TNDEM/NDEM * TRNU - PTF * TRNLOS/(PLANTS * 10.0)

DROOTN = RNDEM/NDEM * TRNU - (1.0 - PTF) * TRNLOS/(PLANTS * 10.0)

TRNLOS is distributed over shoots and roots according to the plant top fraction (PTF) and must also be converted from a unit area basis to a per plant basis. Shoot and root N pools (TOPSN and ROOTN, respectively) can then be updated and new concentrations calculated:

TOPSN = TOPSN + DTOPSN

ROOTN = ROOTN + DROOTN

TANC = TOPSN/TOPWT

RANC = ROOTN/(RTWT - 0.01 * RTWT)

When updating the root concentration allowance is made for the losses in root biomass occurring due to root exudation.

N REDISTRIBUTION DURING GRAIN GROWTH AND GRAIN N DETERMINATION

In many wheat-growing areas when the crop reaches the grain-filling stage soil supplies of N are very low. In these cases the nitrogen requirement of the developing grains is largely satisfied by remobilization of protein from vegetative organs. When nitrogen supply is increased, the proportion of grain N arising from remobilization declines, and the proportion from uptake increases (Vos 1981). Many studies (e.g., Benzian et al., 1983, Terman et al., 1969) have found negative correlations between grain yield and grain protein concentration. Temperature and soil moisture also affect the grain nitrogen content. When constructing the N grain-filling routines, procedures were adopted to closely mimick those predicting grain mass (or carbon) accumulation. In this procedure the rate of grain filling (RGFILL) (mg/day) is determined by temperature and thermal time (DTT).

To define similar functions for the rate of grain N accumulation (RGNFIL) (in micrograms per kernel per degree C day), the controlled environment studies of Sofield et al. (1977), Vos (1981) and Bhullar and Jenner (1985) were used. These studies examined various cultivars over a range of temperature conditions and other treatments. The relationship which best described these studies and mimicked the grain mass accumulation functions was:

RGNFIL = 4.8297 - 3.2488 * DTT + 0.2503 *(TEMPMX - TEMPMN) +

4.3067 * TEMPM

and when the mean temperature is less than 10

RGNFIL = 0.483 * TEMPM

Where TEMPMX, TEMPMN, TEMPM are the maximum, minimum, and mean temperatures (C), respectively. A whole plant grain N sink (NSINK) can then be determined in similar manner to GROGRN.

NSINK = RGNFIL * GPP * 1.E-6 (g N/plant)

Since N stress will affect the rate at which plant tissues can mobilize N and supply it to the grain, an N stress factor NDEF4 from subroutine NFACTO is also introduced.

NSINK = NSINK * NDEF4

If N is present in the plant vegetative tissues (TANC greater than TCNP) the size of the sink is increased. If there is no grain N demand (NSINK = 0) on a day then no grain N accumulation occurs. Two pools of N within the plant are available for translocation, a shoot pool (NPOOL1) and a root pool (NPOOL2). These pools are determined from the N concentration (VANC or RANC) relative to the critical concentration (VMNC or RMNC) and the biomass of the pool (RTWT or TOPWT).

NPOOL1 = TOPWT * (VANC-VMNC)

and

NPOOL2 = RTWT * (RANC-RMNC)

Not all of the N contained within these pools can be immediately mobilized. The fraction of these pools which is labile will depend on the N status of the plant. this fraction (XNF) is calculated by considering the N stress index NDEF2 used for vegetative growth and senescence.

XNF = 0.15 + 0.2 * NDEF2

The labile fraction will range between 15% and 35% of each of the pools depending on the plant N status. The labile poos can be calculated as:

For tops:

TNLAB = XNF * NPOOL1

and roots:

RNLAB = XNF * NPOOL2

The total N available for translocation (NPOOL) is the sum of these two labile pools. When NPOOL is not sufficient to supply the grain N demand (NSINK), NSINK is reduced to NPOOL. If NSINK is greater than that which can be supplied by the tops (TNLAB), then TNLAB is removed from TOPSN and the remaining NSINK which must come from the root pool (RNOUT) is calculated. If (NSINK.GT.TNLAB) Then

TOPSN = TOPSN - TNLAB

RNOUT = NSINK - TNALB

TNLAB = 0

ROOTN = ROOTN - RNOUT

When NSINK is less than TNLAB it can be totally satisfied from the shoot pool and the root pool need not be modified.

TOPSN = TOPSN - NSINK

Following the removal of N from shoot and root pools the simulated tissue concentrations (VANC and RANC) are updated. The total amount of N contained in the grain can then be accumulated.

GRAINN = GRAINN + NSINK

The grain nitrogen concentration will vary daily but is only calculated at the end of the simulatin run (in subroutine PHENOL) as:

GNP = GRAINN/GRNWT

These procedures together with the remainder of the growth routine and the N deficiency indices can provide several pathways by which N stress during grain filling can affect grain yield and grain protein content. First, as N is removed from the vegetative tissues NFAC will become lower. This will in turn lower NDEF4 and lower the sink size for N thus providing for the capability of reduced grain N concentration. Lowering NFAC will also lower NDEF1 which will cause the rate of crop photosynthesis to fall thus lowering the assimilate available for grain filling. A declining NFAC will also speed the rate of senescence which will reduce the leaf area available for photosynthesis. Different temperature regimes during grain filling will also affect the final grain N concentration since the function for RGNFIL is more sensitive to temperature than RGFILL. Soil water stress during grain filling can also increase the grain N concentration since SWDF1 will reduce photosynthesis, lowering assimilate availability and thus not diluting grain N as much as would occur in an unstressed crop.