The amount of daily solar radiation intercepted by the leaf area of the plant is calculated using Beer’s law (Monsi and Saeki, 1953):

$H_{phosyn}=0.5*H_{day}*(1-exp(-k_{\Box}*LAI))$ 5:2.1.1

where $H_{phosyn}$ is the amount of intercepted photosynthetically active radiation on a given day (MJ m$^{-2}$), $H_{day}$ is the incident total solar (MJ m$^{-2}$), $0.5*H_{day}$ is the incident photosynthetically active radiation (MJ m$^{-2}$), $k_{\Box}$ is the light extinction coefficient, and $LAI$ is the leaf area index.

Photosynthetically active radiation is radiation with a wavelength between 400 and 700 mm (McCree, 1972). Direct solar beam radiation contains roughly 45% photosynthetically active radiation while diffuse radiation contains around 60% photosynthetically active radiation (Monteith, 1972; Ross, 1975). The fraction of photosynthetically active radiation will vary from day to day with variation in overcast conditions but studies in Europe and Israel indicate that 50% is a representative mean value (Monteith, 1972; Szeicz, 1974; Stanhill and Fuchs, 1977).

Radiation-use efficiency is the amount of dry biomass produced per unit intercepted solar radiation. The radiation-use efficiency is defined in the plant growth database and is assumed to be independent of the plant’s growth stage. The maximum increase in biomass on a given day that will result from the intercepted photosynthetically active radiation is estimated (Monteith, 1977):

$\Delta bio=RUE*H_{phosyn}$ 5:2.1.2

where $\Delta bio$ is the potential increase in total plant biomass on a given day (kg/ha), $RUE$ is the radiation-use efficiency of the plant (kg/ha⋅(MJ/m$^2$)$^{-1}$ or 10$^{-1}$ g/MJ), and $H_{phosyn}$ is the amount of intercepted photosynthetically active radiation on a given day (MJ m$^{-2}$). Equation 5:2.1.2 assumes that the photosynthetic rate of a canopy is a linear function of radiant energy.

The total biomass on a given day, $d$, is calculated as:

$bio=\sum_{i=1}^{d}\Delta bio_i$ 5:2.1.3

where $bio$ is the total plant biomass on a given day (kg ha$^{-1}$), and $\Delta bio_i$ is the increase in total plant biomass on day $i$ (kg/ha).

The amount of total plant biomass partitioned to the root system is 30-50% in seedlings and decreases to 5-20% in mature plants (Jones, 1985). SWAT+ varies the fraction of total biomass in roots from 0.40 at emergence to 0.20 at maturity. The daily root biomass fraction is calculated with the equation:

$fr_{root}=0.40-0.20*fr_{PHU}$ 5:2.1.21

where $fr_{root}$ is the fraction of total biomass partitioned to roots on a given day in the growing season, and $fr_{PHU}$ is the fraction of potential heat units accumulated for the plant on a given day in the growing season.

Calculation of root depth varies according to plant type. SWAT+ assumes perennials and trees have roots down to the maximum rooting depth defined for the soil throughout the growing season:

$z_{root}=z_{root,mx}$ 5:2.1.22

where $z_{root}$ is the depth of root development in the soil on a given day (mm), and $z_{root,mx}$ is the maximum depth for root development in the soil (mm). The simulated root depth for annuals varies linearly from 10.0 mm at the beginning of the growing season to the maximum rooting depth at $fr_{PHU}$= 0.40 using the equation:

$z_{root}=2.5*fr_{PHU}*z_{root,mx}$ if $fr_{PHU} \le 0.40$ 5:2.1.23

$z_{root}=z_{root,mx}$ if $fr_{PHU} > 0.40$ 5:2.1.24

where $z_{root}$ is the depth of root development in the soil on a given day (mm), $fr_{PHU}$ is the fraction of potential heat units accumulated for the plant on a given day in the growing season, and $z_{root,mx}$ is the maximum depth for root development in the soil (mm). The maximum rooting depth is defined by comparing the maximum potential rooting depth for the plant from the plant growth database (RDMX in crop.dat), and the maximum potential rooting depth for the soil from the soil input file (SOL_ZMX in .sol—if no value is provided for this variable the model will set it to the deepest depth specified for the soil profile). The shallower of these two depths is the value used for $z_{root,mx}$.

With annuals and perennials, the plants are able to reach full maturity within a single calendar year. With trees, a number of years are needed for a plant to transition from a seedling to a sapling to a fully-developed tree. The parameters in the plant growth database related to radiation-use efficiency represent the annual growth for a fully-developed tree. The heat units to maturity input in the management file is also used to simulate growth within a single year, defining for trees and perennials the period within a year bounded by the development of buds at the beginning of the annual growing season and the maturation of plant seeds at the end of the growing season.

To simulate the smaller amount of biomass accumulation seen in seedlings/saplings, tree growth within a single year is limited to a fixed amount determined by the age of the tree relative to the number of years for the tree species to reach full development. Parameters in the plant growth database define the total number of years for trees to reach full development as well as the biomass of a fully-developed tree. Until the trees in an HRU reach full development, the amount of biomass they can accumulate in a single year is limited to:

$bio_{annual}=1000*(\frac{yr_{cur}}{yr_{fulldev}})*bio_{fulldev}$ 5:2.1.9

where $bio_{annual}$ is the amount of biomass a tree can accumulate in a single year (kg/ha), $yr_{cur}$ is the current age of the tree (years), $yr_{fulldev}$ is the number of years for the tree species to reach full development (years), $bio_{fulldev}$ is the biomass of a fully developed tree stand for the specific tree species (metric tons/ha), and 1000 is a conversion factor.

Once the total growth in biomass in a year, $bio$, reaches the annual limit, $bio_{fulldev}$, no more growth occurs until the next year when a new annual limit is calculated. When a tree stand has reached its biomass limit in a year, the increase in plant biomass for a day,$\Delta bio_i$, is set to 0.

Plant growth is modeled by simulating leaf area development, light interception and conversion of intercepted light into biomass assuming a plant species-specific radiation-use efficiency.

Plant maturity is reached when the fraction of potential heat units accumulated, $fr_{PHU}$, is equal to 1.00. Once maturity is reached, the plant ceases to transpire and take up water and nutrients. Simulated plant biomass remains stable until the plant is harvested or killed via a management operation.

Table 5:2-1: SWAT+ input variables that pertain to optimal plant growth.

The change in canopy height and leaf area for annuals and perennials through the growing season as modeled by SWAT+ is illustrated using parameters for Alamo Switchgrass in Figures 5:2-1 and 5:2-2.

In the initial period of plant growth, canopy height and leaf area development are controlled by the optimal leaf area development curve:

$fr_{LAImx}=\frac{fr_{PHU}}{fr_{PHU}+exp(\Box_1 - \Box_2 * fr_{PHU})}$ 5:2.1.10

where $fr_{LAImx}$ is the fraction of the plant’s maximum leaf area index corresponding to a given fraction of potential heat units for the plant, $fr_{PHU}$ is the fraction of potential heat units accumulated for the plant on a given day in the growing season, and $\Box_1$and $\Box_2$ are shape coefficients. The fraction of potential heat units accumulated by a given date is calculated:

$fr_{PHU}=\frac{\sum_{i=1}^d HU_i}{PHU}$ 5:2.1.11

where $fr_{PHU}$ is the fraction of potential heat units accumulated for the plant on day d in the growing season, $HU$ is the heat units accumulated on day $i$ (heat units), and $PHU$ is the total potential heat units for the plant (heat units).

The shape coefficients are calculated by solving equation 5:2.1.10 using two known points ($fr_{LAI,1}$,$fr_{PHU,1}$) and ($fr_{LAI,2}$,$fr_{PHU,2}$):

$\Box_1=1n[\frac{fr_{PHU,1}}{fr_{LAI,1}}-fr_{PHU,1}]+\Box_2*fr_{PHU,1}$ 5:2.1.12

$\Box_2=\frac{(1n[\frac{fr_{PHU,1}}{fr_{LAI,1}}-fr_{PHU,1}]-1n[\frac{fr_{PHU,2}}{fr_{LAI,2}}-fr_{PHU,2}])}{fr_{PHU,2}-fr_{PHU,1}}$ 5:2.1.13

where $\Box_1$ is the first shape coefficient, $\Box_2$ is the second shape coefficient, $fr_{PHU,1}$ is the fraction of the growing season (i.e. fraction of total potential heat units) corresponding to the 1st point on the optimal leaf area development curve, $fr_{LAI,1}$ is the fraction of the maximum plant leaf area index (i.e. fraction of $LAI_{mx}$) corresponding to the 1st point on the optimal leaf area development curve, $fr_{PHU,2}$ is the fraction of the growing season corresponding to the 2nd point on the optimal leaf area development curve, and $fr_{LAI,2}$ is the fraction of the maximum plant leaf area index corresponding to the 2nd point on the optimal leaf area development curve.

The canopy height on a given day is calculated:

$h_c=h_{c,mx}*\sqrt{fr_{LAImx}}$ 5:2.1.14

where $h_c$ is the canopy height for a given day (m), $h_{c,mx}$ is the plant’s maximum canopy height (m), and $fr_{LAI,mx}$ is the fraction of the plant’s maximum leaf area index corresponding to a given fraction of potential heat units for the plant. As can be seen from Figure 5:2-1, once the maximum canopy height is reached, $h_c$ will remain constant until the plant is killed.

For tree stands, the canopy height varies from year to year rather than day to day:

$h_c=h_{c,mx}*(\frac{yr_{cur}}{yr_{fulldev}})$ 5:2.1.15

where $h_c$ is the canopy height for a given day (m), $h_{c,mx}$ is the plant’s maximum canopy height (m), $yr_{cur}$ is the age of the tree (years), and $yr_{fulldev}$ is the number of years for the tree species to reach full development (years).

The amount of canopy cover is expressed as the leaf area index. For annuals and perennials, the leaf area added on day $i$ is calculated:

$\Delta LAI_i=(fr_{LAImx,i}-fr_{LAImx,i-1})*LAI_{mx}*(1-exp(5*(LAI_{i-1}-LAI_{mx})))$

5:2.1.16

while for trees, the leaf area added on day $i$ is calculated:

$\Delta LAI_i=(fr_{LAImx,i}-fr_{LAImx,i-1})*(\frac{yr_{cur}}{yr_{fulldev}})*LAI_{mx}*(1-exp(5*(LAI_{i-1}-(\frac{yr_{cur}}{yr_{fulldev}})*LAI_{mx})))$

5:2.1.17

The total leaf area index is calculated:

$LAI_i=LAI_{i-1}+\Delta LAI_{i}$ 5:2.1.18

where $\Delta LAI_i$ is the leaf area added on day $i$, $LAI_i$ and $LAI_{i-1}$ are the leaf area indices for day $i$ and $i-1$ respectively, $fr_{LAImx,i}$ and $fr_{LAImx,i-1}$ are the fraction of the plant’s maximum leaf area index calculated with equation 5:2.1.10 for day $i$ and $i-1$,$LAI_{mx}$ is the maximum leaf area index for the plant, $yr_{cur}$ is the age of the tree (years), and $yr_{fulldev}$ is the number of years for the tree species to reach full development (years).

Leaf area index is defined as the area of green leaf per unit area of land (Watson, 1947). As shown in Figure 5:2-2, once the maximum leaf area index is reached, $LAI$ will remain constant until leaf senescence begins to exceed leaf growth. Once leaf senescence becomes the dominant growth process, the leaf area index for annuals and perrenials is calculated:

$LAI=LAI_{mx}*\frac{(1-fr_{PHU})}{(1-fr_{PHU,sen})}$ $fr_{PHU}>fr_{PHU,sen}$ 5:2.1.19

while for trees, the calculation is

$LAI=(\frac{yr_{cur}}{yr_{fulldev}})*LAI_{mx}*\frac{(1-fr_{PHU})}{(1-fr_{PHU,sen})}$ $fr_{PHU}>fr_{PHU,sen}$ 5:2.1.20

where $LAI$ is the leaf area index for a given day, $LAI_{mx}$ is the maximum leaf area index, $fr_{PHU}$ is the fraction of potential heat units accumulated for the plant on a given day in the growing season, $fr_{PHU,sen}$ is the fraction of growing season ($PHU$) at which senescence becomes the dominant growth process, $yr_{cur}$ is the number of years of development the tree has accrued (years), and $yr_{fulldev}$ is the number of years for the tree species to reach full development (years).

Radiation-use efficiency is sensitive to variations in atmospheric concentrations and equations have been incorporated into SWAT+ to modify the default radiation-use efficiency values in the plant database for climate change studies. The relationship used to adjust the radiation-use efficiency for effects of elevated is (Stockle et al., 1992):

5:2.1.4

where is the radiation-use efficiency of the plant (kg/ha⋅(MJ/m) or 10 g/MJ), is the concentration of carbon dioxide in the atmosphere (ppmv), and and are shape coefficients.

The shape coefficients are calculated by solving equation 5:2.1.4 using two known points (, ) and (, ):

5:2.1.5

5:2.1.6

where is the first shape coefficient, is the second shape coefficient, is the ambient atmospheric concentration (ppmv), is the radiation-use efficiency of the plant at ambient atmospheric concentration (kg/ha⋅(MJ/m) or 10 g/MJ), is an elevated atmospheric concentration (ppmv), is the radiation-use efficiency of the plant at the elevated atmospheric concentration, , (kg/ha⋅(MJ/m) or 10 g/MJ). Equation 5:2.1.4 was developed when the ambient atmospheric concentration was 330 ppmv and is valid for carbon dioxide concentrations in the range 330-660 ppmv. Even though the ambient atmospheric concentration of carbon dioxide is now higher than 330 ppmv, this value is still used in the calculation. If the concentration used in the simulation is less than 330 ppmv, the model defines RUE = .

Stockle and Kiniry (1990) have shown that a plant’s radiation-use efficiency is affected by vapor pressure deficit. For a plant, a threshold vapor pressure deficit is defined at which the plant’s radiation-use efficiency begins to drop in response to the vapor pressure deficit. The adjusted radiation-use efficiency is calculated:

if 5:2.1.7

if 5:2.1.8

where is the radiation-use efficiency adjusted for vapor pressure deficit (kg/ha⋅(MJ/m) or 10 g/MJ), is the radiation-use efficiency for the plant at a vapor pressure deficit of 1 kPa (kg/ha⋅(MJ/m) or 10 g/MJ), is the rate of decline in radiation-use efficiency per unit increase in vapor pressure deficit (kg/ha⋅(MJ/m)⋅kPa or (10 g/MJ)⋅kPa), is the vapor pressure deficit (kPa), and is the threshold vapor pressure deficit above which a plant will exhibit reduced radiation-use efficiency (kPa). The radiation-use efficiency value reported for the plant in the plant growth database, , or adjusted for elevated carbon dioxide levels (equation 5:2.1.4) is the value used for . The threshold vapor pressure deficit for reduced radiation-use efficiency is assumed to be 1.0 kPa for all plants ().

The radiation-use efficiency is never allowed to fall below 27% of . This minimum value was based on field observations (Kiniry, personal communication, 2001).

Once the potential water uptake has been modified for soil water conditions, the actual amount of water uptake from the soil layer is calculated:

5:2.2.7

where is the actual water uptake for layer (mm HO), is the amount of water in the soil layer on a given day (mm HO), and is the water content of layer at wilting point (mm HO). The total water uptake for the day is calculated:

5:2.2.8

where is the total plant water uptake for the day (mm HO), is the actual water uptake for layer (mm HO), and n is the number of layers in the soil profile. The total plant water uptake for the day calculated with equation 5:2.2.8 is also the actual amount of transpiration that occurs on the day.

5:2.2.9

where is the actual amount of transpiration on a given day (mm HO) and is the total plant water uptake for the day (mm HO).

Table 5:2-2: SWAT+ input variables that pertain to plant water uptake.

If nitrate levels in the root zone are insufficient to meet the demand of a legume, SWAT+ allows the plant to obtain additional nitrogen through nitrogen fixation. Nitrogen fixation is calculated as a function of soil water, soil nitrate content and growth stage of the plant.

5:2.3.9

where is the amount of nitrogen added to the plant biomass by fixation (kg N/ha), is the plant nitrogen demand not met by uptake from the soil (kg N/ha), is the growth stage factor (0.0-1.0), is the soil water factor (0.0-1.0), and is the soil nitrate factor (0.0-1.0). The maximum amount of nitrogen that can be fixed by the plant on a given day is .

Growth stage exerts the greatest impact on the ability of the plant to fix nitrogen. The growth stage factor is calculated:

when 5:2.3.10

when 5:2.3.11

when 5:2.3.12

when 5:2.3.13

when 5:2.3.14

where is the growth stage factor and is the fraction of potential heat units accumulated for the plant on a given day in the growing season. The growth stage factor is designed to reflect the buildup and decline of nitrogen fixing bacteria in the plant roots during the growing season.

The soil nitrate factor inhibits nitrogen fixation as the presence of nitrate in the soil goes up. The soil nitrate factor is calculated:

when 5:2.3.15

when 5:2.3.16

when 5:2.3.17

where is the soil nitrate factor and is the nitrate content of the soil profile (kg NO-N/ha).

The soil water factor inhibits nitrogen fixation as the soil dries out. The soil water factor is calculated:

5:2.3.18

Plant nitrogen uptake is controlled by the plant nitrogen equation. The plant nitrogen equation calculates the fraction of nitrogen in the plant biomass as a function of growth stage given optimal growing conditions.

5:2.3.1

where is the fraction of nitrogen in the plant biomass on a given day, is the normal fraction of nitrogen in the plant biomass at emergence, is the normal fraction of nitrogen in the plant biomass at maturity, is the fraction of potential heat units accumulated for the plant on a given day in the growing season, and and are shape coefficients.

The shape coefficients are calculated by solving equation 5:2.3.1 using two known points (, ) and (, ):

5:2.3.2

5:2.3.3

where is the first shape coefficient, is the second shape coefficient, is the normal fraction of nitrogen in the plant biomass at emergence, is the normal fraction of nitrogen in the plant biomass at 50% maturity, is the normal fraction of nitrogen in the plant biomass at maturity, is the normal fraction of nitrogen in the plant biomass near maturity, is the fraction of potential heat units accumulated for the plant at 50% maturity (=0.5), and is the fraction of potential heat units accumulated for the plant at maturity (=1.0). The normal fraction of nitrogen in the plant biomass near maturity () is used in equation 5:2.3.3 to ensure that the denominator term does not equal 1. The model assumes

To determine the mass of nitrogen that should be stored in the plant biomass on a given day, the nitrogen fraction is multiplied by the total plant biomass:

5:2.3.4

where is the optimal mass of nitrogen stored in plant material for the current growth stage (kg N/ha), is the optimal fraction of nitrogen in the plant biomass for the current growth stage, and is the total plant biomass on a given day (kg ha).

Originally, SWAT+ calculated the plant nitrogen demand for a given day by taking the difference between the nitrogen content of the plant biomass expected for the plant’s growth stage and the actual nitrogen content . This method was found to calculate an excessive nitrogen demand immediately after a cutting (i.e. harvest operation). The equation used to calculate plant nitrogen demand is now

5:2.3.5

where is the potential nitrogen uptake (kg N/ha), is the optimal mass of nitrogen stored in plant material for the current growth stage (kg N/ha), is the actual mass of nitrogen stored in plant material (kg N/ha), is the normal fraction of nitrogen in the plant biomass at maturity, and is the potential increase in total plant biomass on a given day (kg/ha).

The depth distribution of nitrogen uptake is calculated with the function:

5:2.3.6

where is the potential nitrogen uptake from the soil surface to depth (kg N/ha), is the potential nitrogen uptake (kg N/ha), is the nitrogen uptake distribution parameter, is the depth from the soil surface (mm), and is the depth of root development in the soil (mm). Note that equation 5:2.3.6 is similar in form to the depth distribution for water uptake described by equation 5:2.2.1. The potential nitrogen uptake for a soil layer is calculated by solving equation 5:2.3.6 for the depth at the upper and lower boundaries of the soil layer and taking the difference.

5:2.3.7

The actual amount if nitrogen removed from a soil layer is calculated:

For each day of simulation, potential plant growth, i.e. plant growth under ideal growing conditions, is calculated. Ideal growing conditions consist of adequate water and nutrient supply and a favorable climate. Differences in growth between plant species are defined by the parameters contained in the plant growth database.

SWAT+ monitors plant uptake of nitrogen and phosphorus.

If upper layers in the soil profile do not contain enough water to meet the potential water uptake calculated with equation 5:2.2.2, users may allow lower layers to compensate. The equation used to calculate the adjusted potential water uptake is:

$w'_{up,ly}=w_{up,ly}+w_{demand}*epco$ 5:2.2.3

where $w'_{up,ly}$is the adjusted potential water uptake for layer $ly$ (mm H$_2$O), $w_{up,ly}$ is the potential water uptake for layer $ly$ calculated with equation 5:2.2.2 (mm H$_2$O), $w_{demand}$ is the water uptake demand not met by overlying soil layers (mm H$_2$O), and $epco$ is the plant uptake compensation factor. The plant uptake compensation factor can range from 0.01 to 1.00 and is set by the user. As $epco$ approaches 1.0, the model allows more of the water uptake demand to be met by lower layers in the soil. As $epco$ approaches 0.0, the model allows less variation from the depth distribution described by equation 5:2.2.1 to take place.

As the water content of the soil decreases, the water in the soil is held more and more tightly by the soil particles and it becomes increasingly difficult for the plant to extract water from the soil. To reflect the decrease in the efficiency of the plant in extracting water from dryer soils, the potential water uptake is modified using the following equations:

$w''_{up,ly}=w'_{up,ly}*exp[5*(\frac{SW_{ly}}{(.25*AWC_{ly})}-1)]$ when $SW_{ly}< (.25*AWC_{ly})$ 5:2.2.4

$w''_{up,ly}=w'_{up,ly}$ when $SW_{ly} \ge (.25*AWC_{ly})$ 5:2.2.5

where $w''_{up,ly}$is the potential water uptake adjusted for initial soil water content(mm H$_2$O), $w'_{up,ly}$ is the adjusted potential water uptake for layer $ly$ (mm H$_2$O), $SW_{ly}$ is the amount of water in the soil layer on a given day (mm H$_2$O), and $AWC_{ly}$ is the available water capacity for layer $ly$ (mm H$_2$O). The available water capacity is calculated:

$AWC_{ly}=FC_{ly}-WP_{ly}$ 5:2.2.6

where $AWC_{ly}$ is the available water capacity for layer $ly$ (mm H$_2$O), $FC_{ly}$ is the water content of layer $ly$ at field capacity (mm H$_2$O), and $WP_{ly}$ is the water content of layer $ly$ at wilting point (mm H$_2$O).

Plant phosphorus uptake is controlled by the plant phosphorus equation. The plant phosphorus equation calculates the fraction of phosphorus in the plant biomass as a function of growth stage given optimal growing conditions.

$fr_P=(fr_{P,1}-fr_{P,3})*[1-\frac{fr_{PHU}}{fr_{PHU}+exp(p_1-p_2*fr_{PHU})}]+fr_{P,3}$ 5:2.3.19

where $fr_P$ is the fraction of phosphorus in the plant biomass on a given day, $fr_{P,1}$ is the normal fraction of phosphorus in the plant biomass at emergence, $fr_{P,3}$ is the normal fraction of phosphorus in the plant biomass at maturity, $fr_{PHU}$ is the fraction of potential heat units accumulated for the plant on a given day in the growing season, and $p_1$ and $p_2$ are shape coefficients.

The shape coefficients are calculated by solving equation 5:2.3.19 using two known points ($fr_{P,2},fr_{PHU,50\%}$) and ($fr_{P,3},fr_{PHU,100\%}$):

$p_1=1n[\frac{fr_{PHU,50\%}}{(1-\frac{(fr_{P,2}-fr_{P,3})}{fr_{P,1}-fr_{P,3})})}-fr_{PHU,50\%}]+p_2*fr_{PHU,50\%}$ 5:2.3.20

$p_2=\frac{(1n[\frac{fr_{PHU,50\%}}{(1-\frac{(fr_{P,2}-fr_{P,3})}{(fr_{P,1}-fr_{P,3})})}-fr_{PHU,50\%}]-1n[\frac{fr_{PHU,100\%}}{(1-\frac{(fr_{P,\sim3}-fr_{P,3})}{(fr_{P,1}-fr_{P,3})})}-fr_{PHU,100\%}])}{fr_{PHU,100\%}-fr_{PHU,50\%}}$ 5:2.3.21

where $p_1$ is the first shape coefficient, $p_2$ is the second shape coefficient, $fr_{P,1}$ is the normal fraction of phosphorus in the plant biomass at emergence, $fr_{P,2}$ is the normal fraction of phosphorus in the plant biomass at 50% maturity, $fr_{P,3}$ is the normal fraction of phosphorus in the plant biomass at maturity, $fr_{P,\sim 3}$ is the normal fraction of phosphorus in the plant biomass near maturity, $fr_{PHU,50\%}$ is the fraction of potential heat units accumulated for the plant at 50% maturity ($fr_{PHU,50\%}$=0.5), and $fr_{PHU,100\%}$ is the fraction of potential heat units accumulated for the plant at maturity ($fr_{PHU,100\%}$=1.0). The normal fraction of phosphorus in the plant biomass near maturity ($fr_{N,\sim 3}$) is used in equation 5:2.3.21 to ensure that the denominator term $(1-\frac{(fr_{P,\sim3}-fr_{P,3})}{(fr_{P,1}-fr_{P,3})})$does not equal 1. The model assumes $(fr_{P,\sim 3}-fr_{P,3})=0.00001$

To determine the mass of phosphorus that should be stored in the plant biomass for the growth stage, the phosphorus fraction is multiplied by the total plant biomass:

$bio_{P,opt}=fr_P*bio$ 5:2.3.22

where $bio_{P,opt}$ is the optimal mass of phosphorus stored in plant material for the current growth stage (kg P/ha), $fr_P$ is the optimal fraction of phophorus in the plant biomass for the current growth stage, and $bio$ is the total plant biomass on a given day (kg ha$^{-1}$).

Originally, SWAT+ calculated the plant nutrient demand for a given day by taking the difference between the nutrient content of the plant biomass expected for the plant’s growth stage and the actual nutrient content. This method was found to calculate an excessive nutrient demand immediately after a cutting (i.e. harvest operation). The plant phosphorus demand for a given day is calculated:

$P_{up}=1.5*Min \begin{cases} bio_{P,opt}-bio_P \\ 4*fr_{P,3}* \Delta bio \end {cases}$ 5:2.3.23

where $P_{up}$ is the potential phosphorus uptake (kg P/ha), $bio_{P,opt}$ is the optimal mass of phosphorus stored in plant material for the current growth stage (kg P/ha), and $bio_P$ is the actual mass of phosphorus stored in plant material (kg P/ha), $fr_{P,3}$ is the normal fraction of phosphorus in the plant biomass at maturity, and $\Delta bio$ is the potential increase in total plant biomass on a given day (kg/ha). The difference between the phosphorus content of the plant biomass expected for the plant’s growth stage and the actual phosphorus content is multiplied by 1.5 to simulate luxury phosphorus uptake.

The depth distribution of phosphorus uptake is calculated with the function:

$P_{up,z}=\frac{P_{up}}{[1-exp(-\beta_p)]}*[1-exp(-\beta_p*\frac{z}{z_{root}})]$ 5:2.3.24

where $P_{up,z}$ is the potential phosphorus uptake from the soil surface to depth $z$ (kg P/ha), $P_{up}$ is the potential phosphorus uptake (kg P/ha), $\beta _P$ is the phosphorus uptake distribution parameter,$z$ is the depth from the soil surface (mm), and $z_{root}$ is the depth of root development in the soil (mm). The potential phosphorus uptake for a soil layer is calculated by solving equation 5:2.3.24 for the depth at the upper and lower boundaries of the soil layer and taking the difference.

$P_{up,ly}=P_{up,zl}-P_{up,zu}$ 5:2.3.25

where $P_{up,ly}$ is the potential phosphorus uptake for layer $ly$ (kg P/ha), $P_{up,zl}$ is the potential phosphorus uptake from the soil surface to the lower boundary of the soil layer (kg P/ha), and $P_{up,zu}$ is the potential phosphorus uptake from the soil surface to the upper boundary of the soil layer (kg P/ha).

Root density is greatest near the surface, and phosphorus uptake in the upper portion of the soil will be greater than in the lower portion. The depth distribution of phosphorus uptake is controlled by $\beta_p$, the phosphorus uptake distribution parameter, a variable users are allowed to adjust. The illustration of nitrogen uptake as a function of depth for four different uptake distribution parameter values in Figure 5:2-4 is valid for phosphorus uptake as well.

Phosphorus removed from the soil by plants is taken from the solution phosphorus pool. The importance of the phosphorus uptake distribution parameter lies in its control over the maximum amount of solution $P$ removed from the upper layers. Because the top 10 mm of the soil profile interacts with surface runoff, the phosphorus uptake distribution parameter will influence the amount of labile phosphorus available for transport in surface runoff. The model allows lower layers in the root zone to fully compensate for lack of solution P in the upper layers, so there should not be significant changes in phosphorus stress with variation in the value used for $\beta _p$.

The actual amount if phosphorus removed from a soil layer is calculated:

$P_{actualup,ly}=min\lfloor P_{up,ly}+P_{demand},P_{solution,ly}\rfloor$ 5:2.3.26

where $P_{actualup,ly}$ is the actual phosphorus uptake for layer $ly$ (kg P/ha), $P_{up,ly}$ is the potential phosphorus uptake for layer $ly$ (kg P/ha), $P_{demand}$ is the phosphorus uptake demand not met by overlying soil layers (kg P/ha), and $P_{solution,ly}$ is the phosphorus content of the soil solution in layer $ly$ (kg P/ha).